-
-
6.45 KB, 下载次数: 1
TRAJEC.xyz
program analysis Hbond between ON and NN
implicit none
real(4)::N(16,3),O(24,3),H(63,3)
real(4)::S(8,3),C(32,3)
character::u
real(4)::NX(16),NY(16),NZ(16)
real(4)::OX(24),OY(24),OZ(24)
real(4)::HX(63),HY(63),HZ(63)
real(4)::RNO(384),RNH(1008)
integer::i,j,k,nn,m
c open(17,file='NUMH.txt')
open(11,file='TRAJEC.xyz')
do i=1,10
read(11,*)u
read(11,*)u,u
open(14,file='N.txt')
do k=1,16
read(11,*)u,N(k,:)
write(14,*)'N',N(k,:)
enddo
close(14)
do k=1,8
read(11,*)u,S(k,:)
enddo
open(15,file='O.txt')
do k=1,24
read(11,*)u,O(k,:)
write(15,*)'O',O(k,:)
enddo
close(15)
do k=1,32
read(11,*)u,C(k,:)
enddo
open(16,file='H.txt')
nn=0
do k=1,63
nn=1+nn
read(11,*)u,H(k,:)
write(16,*)'H',nn,H(k,:)
enddo
close(16)
call distin(NX,NY,NZ,OX,OY,OZ,HX,HY,HZ,nn,RNO,RNH,u,m)
write(17,*)m
enddo
close(17)
close(11)
end
subroutine distin(NX,NY,NZ,OX,OY,OZ,HX,HY,HZ,nn,RNO,RNH,u,m)
implicit none
real(4)::NX(16),NY(16),NZ(16)
real(4)::OX(24),OY(24),OZ(24)
real(4)::HX(63),HY(63),HZ(63)
real(4)::nn(63),RNO(384),RNH(63)
real(4)::THETA(63),COSYTHETA(63)
real(4)::VNOX(63),VNOY(63),VNOZ(63)
real(4)::VNHX(63),VNHY(63),VNHZ(63)
character::u
integer::i,j,k,m
open(14,file='N.txt')
open(15,file='O.txt')
open(16,file='H.txt')
open(17,file='NUMH.txt')
do i=1,16
read(14,*)u,NX(i),NY(i),NZ(i)
enddo
close(14)
do i=1,24
read(15,*)u,OX(i),OY(i),OZ(i)
enddo
close(15)
do i=1,63
read(16,*)u,nn(i),HX(i),HY(i),HZ(i)
enddo
close(16)
do i=1,16
do j=1,24
RNO(i)=sqrt((NX(i)-OX(j))**2+(NY(i)-OY(j))**2+
X (NZ(i)-OZ(j))**2)
if(RNO(i).lt.2.6) then
do k=1,63
m=0.d0+k
RNH(k)=sqrt((NX(i)-HX(k))**2+(NY(i)-HY(k))**2+
X (NZ(i)-HZ(k))**2)
if(RNH(k).lt.1.4) then
VNOX(k)=OX(j)-NX(i)
VNOY(k)=OY(j)-NY(i)
VNOZ(k)=OZ(j)-NZ(i)
VNHX(k)=HX(k)-NX(i)
VNHY(k)=HY(k)-NY(i)
VNHZ(k)=HZ(k)-NZ(i)
COSYTHETA(k)=(VNOX(k)*VNHX(k)+VNOY(k)*VNHY(k)+VNOZ(k)*VNHZ(k))/
X RNO(i)/RNH(k)
THETA(k)=ACOS(COSYTHETA(k))
if(THETA(k).lt.0.5) then
m=m
write(17,*)m
endif
endif
enddo
endif
enddo
enddo
return
end
fcode 发表于 2017-10-29 10:14
你在主程序里 open(17,file='NUMH.txt') ，在子程序里就不需要Open了。
一般也不需要在子程序里write(17
否 ...
program analysis Hbond
implicit none
real(4)::NX(16),NY(16),NZ(16)
real(4)::OX(24),OY(24),OZ(24)
real(4)::HX(63),HY(63),HZ(63)
real(4)::S(8,3),C(32,3)
real(4)::RNO(384),RNH(63)
real(4)::THETA(63),COSYTHETA(63)
real(4)::VNOX(63),VNOY(63),VNOZ(63)
real(4)::VNHX(63),VNHY(63),VNHZ(63)
character::u
integer::i,j,k,m,l
open(11,file='TRAJEC.xyz')
do l=1,10
read(11,*)u,u,u
do j=1,16
read(11,*)u,NX(j),NY(j),NZ(j)
enddo
do j=1,8
read(11,*)u,S(j,:)
enddo
do j=1,24
read(11,*)u,OX(j),OY(j),OZ(j)
enddo
do j=1,32
read(11,*)u,C(j,:)
enddo
do j=1,63
read(11,*)u,HX(j),HY(j),HZ(j)
enddo
do i=1,16
do j=1,24
open(12,file='NUMH.txt')
RNO(i)=sqrt((NX(i)-OX(j))**2+(NY(i)-OY(j))**2+
X (NZ(i)-OZ(j))**2)
if(RNO(i).lt.2.6) then
do k=1,63
m=0.d0+k
RNH(k)=sqrt((NX(i)-HX(k))**2+(NY(i)-HY(k))**2+
X (NZ(i)-HZ(k))**2)
if(RNH(k).lt.1.4) then
VNOX(k)=OX(j)-NX(i)
VNOY(k)=OY(j)-NY(i)
VNOZ(k)=OZ(j)-NZ(i)
VNHX(k)=HX(k)-NX(i)
VNHY(k)=HY(k)-NY(i)
VNHZ(k)=HZ(k)-NZ(i)
COSYTHETA(k)=(VNOX(k)*VNHX(k)+VNOY(k)*VNHY(k)+VNOZ(k)*VNHZ(k))/
X RNO(i)/RNH(k)
THETA(k)=ACOS(COSYTHETA(k))
if(THETA(k).lt.0.5) then
write(12,*)l
write(12,*)k
endif
endif
enddo
endif
enddo
enddo
enddo
end
fcode 发表于 2017-10-29 14:16
open(12,file='NUMH.txt')
放到循环外面去，你要慎重考虑open语句的位置。
program analysis Hbond
implicit none
real(4)::NX(16),NY(16),NZ(16)
real(4)::OX(24),OY(24),OZ(24)
real(4)::HX(63),HY(63),HZ(63)
real(4)::S(8,3),C(32,3)
real(4)::RNO(384),RNH(63)
real(4)::THETA(63),COSYTHETA(63)
real(4)::VNOX(63),VNOY(63),VNOZ(63)
real(4)::VNHX(63),VNHY(63),VNHZ(63)
character::u
integer::i,j,k,m,l
open(11,file='TRAJEC.xyz')
open(12,file='NUMH.txt')
c do l=1,10
read(11,*)u,u,u
do j=1,16
read(11,*)u,NX(j),NY(j),NZ(j)
enddo
do j=1,8
read(11,*)u,S(j,:)
enddo
do j=1,24
read(11,*)u,OX(j),OY(j),OZ(j)
enddo
do j=1,32
read(11,*)u,C(j,:)
enddo
do j=1,63
read(11,*)u,HX(j),HY(j),HZ(j)
enddo
do i=1,16
do j=1,24
RNO(i)=sqrt((NX(i)-OX(j))**2+(NY(i)-OY(j))**2+
X (NZ(i)-OZ(j))**2)
if(RNO(i).lt.2.6) then
do k=1,63
m=0.d0+k
RNH(k)=sqrt((NX(i)-HX(k))**2+(NY(i)-HY(k))**2+
X (NZ(i)-HZ(k))**2)
if(RNH(k).lt.1.4) then
VNOX(k)=OX(j)-NX(i)
VNOY(k)=OY(j)-NY(i)
VNOZ(k)=OZ(j)-NZ(i)
VNHX(k)=HX(k)-NX(i)
VNHY(k)=HY(k)-NY(i)
VNHZ(k)=HZ(k)-NZ(i)
COSYTHETA(k)=(VNOX(k)*VNHX(k)+VNOY(k)*VNHY(k)+VNOZ(k)*VNHZ(k))/
X RNO(i)/RNH(k)
THETA(k)=ACOS(COSYTHETA(k))
if(THETA(k).lt.0.5) then
write(12,*)l
write(12,*)k
endif
endif
enddo
endif
enddo
enddo
c enddo
end
6.45 KB, 下载次数: 1
TRAJEC.xyz
| 欢迎光临 Fortran Coder (http://bbs.fcode.cn/) | Powered by Discuz! X3.2 |