 |
 |
 |
 |
 |
 |
||
 |
 |
 |
 |
 |
| 5 | 1/1 | 返回列表 |
| 查看: 1051 | 回復: 7 | |||
| 當前只顯示滿足指定條件的回帖，點擊這里查看本話題的所有回帖 | |||
[交流]
【原創(chuàng)】使用MS生成mcm-41，并且插胺基的程序，與大家討論下 已有2人參與
|
|||
|
integer natom,natom0,nho,namino ************************************************************ * Number of atoms in the original structure is 9992. * But the parameter natom should include the number of atoms added subsequently. * So here the value of natom is set to 15000 ************************************************************ parameter (natom=15000,natom0=9992,nho=2319,namino=435) character a(natom)*5,fx(natom)*4,fft(natom)*5,atomname(natom)*2 integer occupation(natom),kjishu,ron,nOtemp,kstop,natom_add,Tatom integer Ohydroxy_SN(nho),Temp_SN double precision xc(natom),yc(natom),zc(natom) double precision rox,roy,roz,Ohydroxy(nho,4),Otemp(nho,4) double precision OTa(4),OTb(4),temp(3),list3(nho*3,12) integer templist3(nho*3,3),Nra double precision xtemp,ytemp,ztemp,xfinal,yfinal,zfinal * integer NSiT,NSi_S,NSi_C,NCT,NC_S,NC_C,NNT,NN_S,NN_C * character NSi_CC*1,NC_CC*1,NN_CC*1 real distance,search_step,dis1,dis2,dis3,lbond,charge(natom) ****************************************************************** * define global variables ****************************************************************** common charge common xc,yc,zc ****************************************************************** * Read the input file ****************************************************************** open(10,file='MCM41-final.car',status='old') do 20 i=1,natom0 read(10,*)a(i),xc(i),yc(i),zc(i),fx(i),occupation(i),fft(i), & atomname(i),charge(i) 20 continue close(10) ******************************************************************* * write the initial file to check whether the initial structure is read correctly. ******************************************************************* open(30,file='output.car',access='append') do 40 i=1,natom0 write(30,888)a(i),xc(i),yc(i),zc(i),fx(i),occupation(i),fft(i), & atomname(i),charge(i) 40 continue close(30) natom_add=0 Tatom=natom0+natom_add * NSiT=0 * NSi_S=0 * NSi_C=0 * NCT=0 * NC_S=0 * NC_C=0 * NNT=0 * NN_S=0 * NN_C=0 do 45 itt=1,namino ********************************************************* * add Si to the chosen Oxygen atom ********************************************************* call HO_list(Ohydroxy,Ohydroxy_SN,kjishu) Nra=int(RAN2(IDUM)*kjishu) Temp_SN=Ohydroxy_SN(Nra) xtemp=Ohydroxy(Nra,2) ytemp=Ohydroxy(Nra,3) ztemp=Ohydroxy(Nra,4) lbond=1.90 call addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) natom_add=natom_add+1 Tatom=natom0+natom_add * NSiT=NSiT+1 * NSi_C=INT((NSiT+60)/99) * NSi_S=NSiT+60-99*NSi_C * * if(NSi_C.eq.0)NSi_CC='T' * if(NSi_C.eq.1)NSi_CC='U' * if(NSi_C.eq.2)NSi_CC='V' * if(NSi_C.eq.3)NSi_CC='W' * if(NSi_C.eq.4)NSi_CC='X' * if(NSi_C.eq.5)NSi_CC='Y' * if(NSi_C.eq.6)NSi_CC='Z' * a(Tatom)='Si'//'NSi_S'//'NSi_CC' a(Tatom)='Si' xc(Tatom)=xfinal yc(Tatom)=yfinal zc(Tatom)=zfinal fx(Tatom)='XXXX' occupation(Tatom)=1 fft(Tatom)='Si3' atomname(Tatom)='Si' charge(Tatom)=5.000 open(140,file='amino.car',access='append') write(140,888)a(Tatom),xc(Tatom), & yc(Tatom),zc(Tatom), & fx(Tatom),occupation(Tatom), & fft(Tatom),atomname(Tatom), & charge(Tatom) close(140) * do 2060 i=1,nho-2 * * do 2065 ix=-2,kjishu+2 * if(Ohydroxy(i,1).gt.templist(ix)-0.1.and. * & Ohydroxy(i,1).lt.templist(ix)+0.1)then * goto 2060 * endif *2065 continue * * do 2070 j=i+1,nho-1 * * do 2075 ix=-2,kjishu+2 * if(Ohydroxy(j,1).gt.templist(ix)-0.1.and. * & Ohydroxy(j,1).lt.templist(ix)+0.1)then * goto 2070 * endif *2075 continue * * do 2080 k=j+1,nho * * do 2085 ix=-2,kjishu+2 * if(Ohydroxy(k,1).gt.templist(ix)-0.1.and. * & Ohydroxy(k,1).lt.templist(ix)+0.1)then * goto 2080 * endif *2085 continue * * dis1=sqrt((Ohydroxy(i,2)-Ohydroxy(j,2))**2+ * & (Ohydroxy(i,3)-Ohydroxy(j,3))**2+ * & (Ohydroxy(i,4)-Ohydroxy(j,4))**2) * * dis2=sqrt((Ohydroxy(i,2)-Ohydroxy(k,2))**2+ * & (Ohydroxy(i,3)-Ohydroxy(k,3))**2+ * & (Ohydroxy(i,4)-Ohydroxy(k,4))**2) * dis3=sqrt((Ohydroxy(j,2)-Ohydroxy(k,2))**2+ * & (Ohydroxy(j,3)-Ohydroxy(k,3))**2+ * & (Ohydroxy(j,4)-Ohydroxy(k,4))**2) * if (dis1.gt.1.8.and.dis1.lt.2.6.and. * & dis2.gt.1.8.and.dis1.lt.2.6.and. * & dis3.gt.1.8.and.dis1.lt.2.6)then * kjishu=kjishu+1 * do 2090 imm=1,4 * list3(kjishu,imm)=Ohydroxy(i,imm) *2090 continue * templist3(kjishu,1)=int(Ohydroxy(i,1)) * charge(Ohydroxy_SN(i))=2.0 * do 2100 imm=1,4 * list3(kjishu,imm+4)=Ohydroxy(j,imm) *2100 continue * templist3(kjishu,2)=int(Ohydroxy(j,1)) * charge(Ohydroxy_SN(j))=2.0 * do 2110 imm=1,4 * list3(kjishu,imm+8)=Ohydroxy(k,imm) *2110 continue * templist3(kjishu,3)=int(Ohydroxy(k,1)) * charge(Ohydroxy_SN(k))=2.0 |
回復此樓» 猜你喜歡
求調劑
已經有5人回復
-
材料工程專碩 348分求調劑
已經有5人回復
-
一志愿南大，0703化學，分數(shù)336，求調劑
已經有3人回復
-
322求調劑
已經有3人回復
-
0703化學調劑
已經有4人回復
-
化學調劑
已經有4人回復
-
297求調劑
已經有3人回復
-
生物學調劑
已經有3人回復
-
工科0856求調劑
已經有3人回復
-
求調劑
已經有3人回復
» 本主題相關價值貼推薦，對您同樣有幫助:
- 國家自然基金pdf格式的生成
已經有12人回復
- 經典網(wǎng)格生成書籍【轉載】
已經有20人回復
- 求助：討論不同壓力下蒸汽消耗的關系
已經有9人回復
- ms5.5 在win8 總是出問題
已經有13人回復
- 格式試劑怎么做GC或GC-MS？
已經有10人回復
- MCM-41表面接枝
已經有7人回復
- 單氰胺在酸性條件下（濃鹽酸）會不會生成氫氰酸
已經有20人回復
- 不知道怎么回答審稿人的問題，求助大家,(關于分子篩MCM-41的）
已經有10人回復
- 【求助成功】MS模擬MCM-41過程中forcite filed計算的問題
已經有12人回復
|
elseif(xtemp.gt.xcenter+5)then ia=-20 ib=0 ic=1 minlbond=lbond-0.1 maxlbond=lbond+0.1 search_step=0.1 do 1001 i=ia,ib,ic do 1011 j=ja,jb,jc do 1021 k=-20,20 targetx=xtemp+search_step*i targety=ytemp+search_step*j targetz=ztemp+search_step*k dis1=sqrt((targetx-xtemp)**2+ & (targety-ytemp)**2+ & (targetz-ztemp)**2) if(dis1.gt.minlbond.and.dis1.lt.maxlbond)then do 1031 im=1,Tatom if(im.ne.Temp_SN)then dis2=sqrt((targetx-xc(im))**2+ & (targety-yc(im))**2+ & (targetz-zc(im))**2) if(dis2.lt.3.0)then goto 1021 endif endif * if(dis1.gt.minlbond.and.dis1.lt.maxlbond)then * do 1031 im=1,Tatom * dis2=sqrt((targetx-xc(im))**2+ * & (targety-yc(im))**2+ * & (targetz-zc(im))**2) * if(im.ne.Temp_SN.and.dis2.lt.3.0)then * goto 1021 * endif 1031 continue goto 1041 endif 1021 continue 1011 continue 1001 continue 1041 xfinal=targetx yfinal=targety zfinal=targetz else ia=20 ib=-20 ic=-1 minlbond=lbond-0.1 maxlbond=lbond+0.1 search_step=0.1 do 1002 j=ja,jb,jc do 1012 i=ia,ib,ic do 1022 k=-20,20 targetx=xtemp+search_step*i targety=ytemp+search_step*j targetz=ztemp+search_step*k dis1=sqrt((targetx-xtemp)**2+ & (targety-ytemp)**2+ & (targetz-ztemp)**2) if(dis1.gt.minlbond.and.dis1.lt.maxlbond)then do 1032 im=1,Tatom if(im.ne.Temp_SN)then dis2=sqrt((targetx-xc(im))**2+ & (targety-yc(im))**2+ & (targetz-zc(im))**2) if(dis2.lt.3.0)then goto 1022 endif endif * if(dis1.gt.minlbond.and.dis1.lt.maxlbond)then * do 1032 im=1,Tatom * dis2=sqrt((targetx-xc(im))**2+ * & (targety-yc(im))**2+ * & (targetz-zc(im))**2) * if(im.ne.Temp_SN.and.dis2.lt.3.0)then * goto 1022 * endif 1032 continue goto 1042 endif 1022 continue 1012 continue 1002 continue 1042 xfinal=targetx yfinal=targety zfinal=targetz endif endsubroutine *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& * Random Number Generator *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& FUNCTION RAN2(IDUM) INTEGER idum,IM1,IM2,IMM1,IA1,IA2,IQ1,IQ2,IR1,IR2,NTAB,NDIV REAL ran2,AM,EPS,RNMX PARAMETER (IM1=2147483563,IM2=2147483399,AM=1./IM1,IMM1=IM1-1, *IA1=40014,IA2=40692,IQ1=53668,IQ2=52774,IR1=12211,IR2=3791, *NTAB=32,NDIV=1+IMM1/NTAB,EPS=1.2e-7,RNMX=1.-EPS) INTEGER idum2,j,k,iv(NTAB),iy SAVE iv,iy,idum2 DATA idum2/123456789/, iv/NTAB*0/, iy/0/ if (idum.le.0) then idum=max(-idum,1) idum2=idum do 11 j=NTAB+8,1,-1 k=idum/IQ1 idum=IA1*(idum-k*IQ1)-k*IR1 if (idum.lt.0) idum=idum+IM1 if (j.le.NTAB) iv(j)=idum 11 continue iy=iv(1) endif k=idum/IQ1 idum=IA1*(idum-k*IQ1)-k*IR1 if (idum.lt.0) idum=idum+IM1 k=idum2/IQ2 idum2=IA2*(idum2-k*IQ2)-k*IR2 if (idum2.lt.0) idum2=idum2+IM2 j=1+iy/NDIV iy=iv(j)-idum2 iv(j)=idum if(iy.lt.1)iy=iy+IMM1 ran2=min(AM*iy,RNMX) return END *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& |
4樓2010-05-31 22:20:06
|
* goto 2061 * endif *2080 continue *2070 continue *2060 continue * open(2120,file='list3.car',access='append') * do 2130 i=1,kjishu * write(2120,777)int(list3(i,1)),list3(i,2),list3(i,3),list3(i,4), * & templist3(i,1) * write(2120,777)int(list3(i,5)),list3(i,6),list3(i,7),list3(i,8), * & templist3(i,2) * write(2120,777)int(list3(i,9)),list3(i,10),list3(i,11),list3(i,12) * & ,templist3(i,3) * write(2120,*)'' *2130 continue *777 Format(I4,3X,F12.9,2X,F13.9,3X,F12.9,3X,I4) * close(2120) *************************************************** * renew hydroxy oxygen list *************************************************** * do 2135 i=1,nho * Ohydroxy_SN(i)=0 * do 2136 j=1,4 * Ohydroxy(i,j)=0 *2136 continue *2135 continue * kjishu=0 * do 2140 i=1,natom0 * if (charge(i).eq.1.0) then * kjishu=kjishu+1 * Ohydroxy_SN(kjishu)=i * Ohydroxy(kjishu,1)=kjishu * Ohydroxy(kjishu,2)=xc(i) * Ohydroxy(kjishu,3)=yc(i) * Ohydroxy(kjishu,4)=zc(i) * endif *2140 continue * open(2150,file='hydrogen oxygen-1.car',access='append') * do 2160 i=1,nho * write(2150,*)(Ohydroxy(i,j),j=1,4),Ohydroxy_SN(i) *2160 continue * close(2150) * stop ************************************ * choose a initial oxygen atom randomly ************************************ * ron=int(ran2(idum)*nho) * rox=Ohydroxy(nho,2) * roy=Ohydroxy(nho,3) * roz=Ohydroxy(nho,4) * nOtemp=0 * do 60 i=1,nho * distance=sqrt((rox-Ohydroxy(i,2))**2+(roy-Ohydroxy(i,3))**2+ * & (roz-Ohydroxy(i,4))**2) * if(distance.gt.1.0.and.distance.lt.2.6)then * nOtemp=nOtemp+1 * Otemp(nOtemp,1)=nOtemp * Otemp(nOtemp,2)=Ohydroxy(i,2) * Otemp(nOtemp,3)=Ohydroxy(i,3) * Otemp(nOtemp,4)=Ohydroxy(i,4) * endif *60 continue * kstop=0 * do 70 i=1,nOtemp-1 * do 80 j=i+1,nOtemp * distance=sqrt((Ohydroxy(i,2)-Ohydroxy(j,2))**2+ * & (Ohydroxy(i,3)-Ohydroxy(j,3))**2+ * & (Ohydroxy(i,4)-Ohydroxy(j,4))**2) * if(distance.gt.1.0.and.distance.lt.2.6)then ************************************************* * Mark of interrupt service routine 2 ************************************************* * kstop=1 * OTa(1)=i * OTa(2)=Otemp(i,2) * OTa(3)=Otemp(i,3) * OTa(4)=Otemp(i,4) * OTb(1)=j * OTb(2)=Otemp(j,2) * OTb(3)=Otemp(j,3) * OTb(4)=Otemp(j,4) * goto 90 * endif *80 continue *70 continue ******************************************** * warning 1 ******************************************** * if(kstop.eq.0)then * write(*,*)'Warning',kstop,': no tri-grafting any more' * stop * endif ******************************************** * warning 1 ******************************************** * kstop=0 *90 search_step=0.01 * do 100 i=-300,300 * do 110 j=-300,300 * do 120 k=-300,300 * temp(1)=rox+search_step*i * temp(2)=roy+search_step*j * temp(3)=roz+search_step*k * dis1=sqrt((temp(1)-rox)**2+ * & (temp(2)-roy)**2+ * & (temp(3)-roz)**2) * dis2=sqrt((temp(1)-OTa(2))**2+ * & (temp(2)-OTa(3))**2+ * & (temp(3)-OTa(4))**2) * dis3=sqrt((temp(1)-OTb(2))**2+ * & (temp(2)-OTb(3))**2+ * & (temp(3)-OTb(4))**2) * if (dis1.gt.1.8.and.dis1.lt.2.0.and. * & dis2.gt.1.8.and.dis1.lt.2.0.and. * & dis3.gt.1.8.and.dis1.lt.2.0)then * do 130 im=1,natom+natom_add * distance=sqrt((temp(1)-xc(im))**2+ * & (temp(2)-yc(im))**2+ * & (temp(3)-zc(im))**2) * if(distance.le.3.0)then * goto 120 * endif *130 continue ************************************************* * Mark of interrupt service routine 2 ************************************************* * kstop=2 * natom_add=natom_add+1 * Tatom=natom0+natom_add * a(Tatom)='Si' * xc(Tatom)=temp(1) * yc(Tatom)=temp(2) * zc(Tatom)=temp(3) * fx(Tatom)='XXXX' * occupation(Tatom)=1 * fft(Tatom)='Si3' * atomname(Tatom)='Si' * charge(Tatom)=3.000 * open(140,file='amino.car',access='append') * write(140,888)a(Tatom),xc(Tatom), * & yc(Tatom),zc(Tatom), * & fx(Tatom),occupation(Tatom), * & fft(Tatom),atomname(Tatom), * & charge(Tatom) * close(140) ****************************************************************** * add C1 atom on the chosen Si atom ****************************************************************** lbond=1.93 xtemp=xc(Tatom) ytemp=yc(Tatom) ztemp=zc(Tatom) Temp_SN=Tatom call addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) natom_add=natom_add+1 Tatom=natom0+natom_add * NCT=NCT+1 * NC_C=INT(NCT/99) * NC_S=NCT-99*NC_C * if(NC_C.eq.0)NC_CC='' * if(NC_C.eq.1)NC_CC='A' * if(NC_C.eq.2)NC_CC='B' * if(NC_C.eq.3)NC_CC='C' * if(NC_C.eq.4)NC_CC='D' * if(NC_C.eq.5)NC_CC='E' * if(NC_C.eq.6)NC_CC='F' * a(Tatom)='C'//'NC_S'//'NC_CC' a(Tatom)='C' xc(Tatom)=xfinal yc(Tatom)=yfinal zc(Tatom)=zfinal fx(Tatom)='XXXX' occupation(Tatom)=1 fft(Tatom)='C_3' atomname(Tatom)='C' charge(Tatom)=6.000 open(150,file='amino.car',access='append') write(150,888)a(Tatom),xc(Tatom), & yc(Tatom),zc(Tatom), & fx(Tatom),occupation(Tatom), & fft(Tatom),atomname(Tatom), & charge(Tatom) close(150) ****************************************************************** * add C2 atom on the chosen C1 atom ****************************************************************** lbond=1.50 xtemp=xc(Tatom) ytemp=yc(Tatom) ztemp=zc(Tatom) Temp_SN=Tatom call addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) natom_add=natom_add+1 Tatom=natom0+natom_add * NCT=NCT+1 * NC_C=INT(NCT/99) * NC_S=NCT-99*NC_C * if(NC_C.eq.0)NC_CC='' * if(NC_C.eq.1)NC_CC='A' * if(NC_C.eq.2)NC_CC='B' * if(NC_C.eq.3)NC_CC='C' * if(NC_C.eq.4)NC_CC='D' * if(NC_C.eq.5)NC_CC='E' * if(NC_C.eq.6)NC_CC='F' * a(Tatom)='C'//'NC_S'//'NC_CC' a(Tatom)='C' xc(Tatom)=xfinal yc(Tatom)=yfinal zc(Tatom)=zfinal fx(Tatom)='XXXX' occupation(Tatom)=1 fft(Tatom)='C_3' atomname(Tatom)='C' charge(Tatom)=7.000 open(160,file='amino.car',access='append') write(160,888)a(Tatom),xc(Tatom), & yc(Tatom),zc(Tatom), & fx(Tatom),occupation(Tatom), & fft(Tatom),atomname(Tatom), & charge(Tatom) close(160) ****************************************************************** * add C3 atom on the chosen C2 atom ****************************************************************** lbond=1.50 xtemp=xc(Tatom) ytemp=yc(Tatom) ztemp=zc(Tatom) Temp_SN=Tatom call addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) natom_add=natom_add+1 Tatom=natom0+natom_add * NCT=NCT+1 * NC_C=INT(NCT/99) * NC_S=NCT-99*NC_C * if(NC_C.eq.0)NC_CC='' * if(NC_C.eq.1)NC_CC='A' * if(NC_C.eq.2)NC_CC='B' * if(NC_C.eq.3)NC_CC='C' * if(NC_C.eq.4)NC_CC='D' * if(NC_C.eq.5)NC_CC='E' * if(NC_C.eq.6)NC_CC='F' * a(Tatom)='C'//'NC_S'//'NC_CC' a(Tatom)='C' xc(Tatom)=xfinal yc(Tatom)=yfinal zc(Tatom)=zfinal fx(Tatom)='XXXX' occupation(Tatom)=1 fft(Tatom)='C_3' atomname(Tatom)='C' charge(Tatom)=8.000 open(170,file='amino.car',access='append') write(170,888)a(Tatom),xc(Tatom), & yc(Tatom),zc(Tatom), & fx(Tatom),occupation(Tatom), & fft(Tatom),atomname(Tatom), & charge(Tatom) close(170) |
2樓2010-05-31 22:19:24
|
****************************************************************** * add N atom on the chosen C3 atom ****************************************************************** lbond=1.53 xtemp=xc(Tatom) ytemp=yc(Tatom) ztemp=zc(Tatom) Temp_SN=Tatom call addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) natom_add=natom_add+1 Tatom=natom0+natom_add * NNT=NNT+1 * NN_C=INT(NNT/99) * NN_S=NNT-99*NN_C * if(NN_C.eq.0)NN_CC='' * if(NN_C.eq.1)NN_CC='A' * if(NN_C.eq.2)NN_CC='B' * if(NN_C.eq.3)NN_CC='C' * if(NN_C.eq.4)NN_CC='D' * if(NN_C.eq.5)NN_CC='E' * if(NN_C.eq.6)NN_CC='F' * a(Tatom)='N'//'NN_S'//'NN_CC' a(Tatom)='N' xc(Tatom)=xfinal yc(Tatom)=yfinal zc(Tatom)=zfinal fx(Tatom)='XXXX' occupation(Tatom)=1 fft(Tatom)='N_3' atomname(Tatom)='N' charge(Tatom)=9.000 open(180,file='amino.car',access='append') write(180,888)a(Tatom),xc(Tatom), & yc(Tatom),zc(Tatom), & fx(Tatom),occupation(Tatom), & fft(Tatom),atomname(Tatom), & charge(Tatom) close(180) 45 continue **************************************************** * rewrite the atom serial. **************************************************** open(46,file='Satom.car',status='old') do 47 i=natom0+1,Tatom read(46,1888)a(i) 47 continue close(46) 1888 format(A5) * endif * *120 continue *110 continue *100 continue *************************************************** * write the output file *************************************************** open(190,file='output-final.car',access='append') do 200 i=1,Tatom write(190,888)a(i),xc(i),yc(i),zc(i),fx(i),occupation(i),fft(i), & atomname(i),charge(i) 200 continue close(190) 888 format(A5,3X,F12.9,2X,F13.9,3X,F12.9,1X,A4,1X,I1,6X,A5,3X,A2,F7.3) **************************************************** **************************************************** end *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& * To make a list of hydroxy oxygen. *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& subroutine HO_list(Ohydroxy,Ohydroxy_SN,kjishu) integer nho,natom,natom0 parameter (natom=15000,natom0=9992,nho=2319) integer kjishu,Ohydroxy_SN(nho) double precision Ohydroxy(nho,4) ******************************************************************* * declare the global varialbe ******************************************************************* real charge(natom) double precision xc(natom),yc(natom),zc(natom) common charge common xc,yc,zc kjishu=0 do 50 i=1,natom0 if (charge(i).eq.1.0) then kjishu=kjishu+1 Ohydroxy_SN(kjishu)=i Ohydroxy(kjishu,1)=kjishu Ohydroxy(kjishu,2)=xc(i) Ohydroxy(kjishu,3)=yc(i) Ohydroxy(kjishu,4)=zc(i) endif 50 continue ****************************************************************** * write out the hydroxy oxygen list ****************************************************************** open(51,file='HO_list.car',access='append') do 52 i=1,kjishu write(51,999)int(Ohydroxy(i,1)),Ohydroxy(i,2), & Ohydroxy(i,3),Ohydroxy(i,4),Ohydroxy_SN(i) 52 continue write(51,*)'***********************************************' write(51,*)'' close(51) 999 Format(I4,3X,F12.9,2X,F13.9,3X,F12.9,3X,I4) endsubroutine *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& * add new atom *&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& subroutine addatom(xtemp,ytemp,ztemp,lbond,Tatom,Temp_SN, & xfinal,yfinal,zfinal) integer natom parameter(natom=15000) real search_step,lbond,minlbond,maxlbond double precision xtemp,ytemp,ztemp,targetx,targety,targetz double precision dis1,dis2,xc(natom),yc(natom),zc(natom) double precision xfinal,yfinal,zfinal double precision center(14,2),xcenter,ycenter,lc,lct integer ia,ib,ic,ja,jb,jc integer Tatom,Temp_SN ***************************************************************** * difine global variable ***************************************************************** common xc,yc,zc ***************************************************************** ***************************************************************** open(3000,file='center.txt',status='old') do 3010 i=1,14 read(3000,*)(center(i,j),j=1,2) 3010 continue close(3000) * open(3011,file='center-1.txt',access='append') * do 3012 i=1,14 * write(3011,*)(center(i,j),j=1,2) *3012 continue * close(3011) lc=1000 xcenter=0 ycenter=0 do 3020 i=1,14 lct=sqrt((xtemp-center(i,1))**2+(ytemp-center(i,2))**2) if (lct.lt.lc)then lc=lct xcenter=center(i,1) ycenter=center(i,2) endif 3020 continue if(ytemp.le.ycenter)then ja=20 jb=0 jc=-1 jd=j else ja=-20 jb=0 jc=1 jd=j endif if(xtemp.le.xcenter-5)then ia=20 ib=0 ic=-1 minlbond=lbond-0.1 maxlbond=lbond+0.1 search_step=0.1 do 1000 i=ia,ib,ic do 1010 j=ja,jb,jc do 1020 k=-20,20 targetx=xtemp+search_step*i targety=ytemp+search_step*j targetz=ztemp+search_step*k dis1=sqrt((targetx-xtemp)**2+ & (targety-ytemp)**2+ & (targetz-ztemp)**2) if(dis1.gt.minlbond.and.dis1.lt.maxlbond)then do 1030 im=1,Tatom if(im.ne.Temp_SN)then dis2=sqrt((targetx-xc(im))**2+ & (targety-yc(im))**2+ & (targetz-zc(im))**2) if(dis2.lt.3.0)then goto 1020 endif endif * if(im.ne.Temp_SN.and.dis2.lt.3.0)then * goto 1020 * endif 1030 continue goto 1040 endif 1020 continue 1010 continue 1000 continue 1040 xfinal=targetx yfinal=targety zfinal=targetz |
3樓2010-05-31 22:19:49
5樓2010-05-31 22:22:23
| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 0703化學調劑 +4 | 妮妮ninicgb 2026-03-21 | 4/200 |
|
|---|---|---|---|---|
|
[考研] 317求調劑 +9 | 申子申申 2026-03-19 | 15/750 |
|
|
|
[考研] 296求調劑 +4 | www_q 2026-03-20 | 4/200 |
|
|
|
[考研] 求調劑 +3 | .m.. 2026-03-21 | 4/200 |
|
|
|
[考研] 0805材料320求調劑 +3 | 深海物語 2026-03-20 | 3/150 |
|
|
|
[考研] 301求調劑 +10 | yy要上岸呀 2026-03-17 | 10/500 |
|
|
|
[考研] 化學求調劑 +4 | 臨澤境llllll 2026-03-17 | 5/250 |
|
|
|
[考研] 考研調劑求學校推薦 +3 | 伯樂29 2026-03-18 | 5/250 |
|
|
|
[考研] 求調劑一志愿南京航空航天大學289分 +3 | @taotao 2026-03-19 | 3/150 |
|
|
|
[考研] 本人考085602 化學工程 專碩 +19 | 不知道叫什么！ 2026-03-15 | 21/1050 |
|
|
|
[考研] 260求調劑 +3 | 朱芷琳 2026-03-20 | 3/150 |
|
|
|
[考研] 求調劑 +3 | eation27 2026-03-20 | 3/150 |
|
|
|
[考研] 廣西大學家禽遺傳育種課題組2026年碩士招生（接收計算機專業(yè)調劑） +3 | 123阿標 2026-03-17 | 3/150 |
|
|
|
[考研] 0856調劑，是學校就去 +8 | sllhht 2026-03-19 | 9/450 |
|
|
|
[考研] 一志愿中國海洋大學，生物學，301分，求調劑 +5 | 1孫悟空 2026-03-17 | 6/300 |
|
|
|
[考研]
 320求調劑0856
+3
|
不想起名字112 2026-03-19 | 3/150 |
|
|
|
[考研] 材料工程專碩調劑 +5 | 204818@lcx 2026-03-17 | 6/300 |
|
|
|
[考研] 0854可跨調劑，一作一項核心論文五項專利，省、國級證書40+數(shù)一英一287 +8 | 小李0854 2026-03-16 | 8/400 |
|
|
|
[考研] 277調劑 +5 | 自由煎餅果子 2026-03-16 | 6/300 |
|
|
|
[考研] 288求調劑 +4 | 奇點0314 2026-03-14 | 4/200 |
|