以下程式

program geom_fo_step_structure
! This program cut the coordinate of each atom at some one step in a
trajectory.
character(16) :: filename,Ofilename
logical :: file_exist
integer :: Natom,Nstep,Nrow,step,i,j,k,l,m,n
character(51) :: C(30000)
character(3) :: atoms(50)
print *, "Please enter a [filename] for loading your data"
print *, " Warning : the length of the filenam must be less than 16"
read(*,*) filename
inquire (file=filename,exist=file_exist)
if (file_exist .eq. .true.) then
print *, "Enter the [number of atoms] in your job"
print *, " Warning : the number of atoms should be less than
50"
read(*,*) Natom
print *, "Enter the [maximum step No.] in your trajectory job"
print *, " Warning : the product of atoms No. and steps No.
should be less than 30000"
read(*,*) Nstep
Nrow = Natom*Nstep
print *, "Enter the [atom symbols] in the label order"
read(*,*) (atoms(i), i = 1, Natom)
998 print *, "Enter the [number of the step] what you want to
load"
read(*,*) step
if (step .gt. Nstep) then
print *, "The step No. what you want to load is larger
than the maximum step in your job"
go to 998
end if
open (unit=9,status='scratch')
write(9,"(i5,a4)") step,"step"
rewind(9)
read(9,*) Ofilename
open (unit=7,file=filename,status='old',action='read')
open (unit=8,file=Ofilename,status='new')
do j = 1, Nrow
read(7,2,advance='yes',end=3) C(j)
2 format(a51)
3 end do
k = ((step-1)*Natom) + 1
l = step*Natom
do m = k, l
if (mod(m,Natom) .ne. 0) then
n = mod(m,Natom)
else
n = Natom
end if
write(8,4,advance='yes') atoms(n),C(m)
4 format(a3,a51)
end do
print *, "The coordinates what you want were written in :
",Ofilename
else
print *, filename," doesn't exist"
go to 999
end if
999 stop
end

--
最後一個剪輯程式,功能是剪輯前一個程式製造出來的座標檔案,對其剪輯出所指定的步

數,並在所有座標前端加上原子符號。自動命檔名為"n"step。

行號998為防呆,避免指定步數超過最大步數。

mongqiu 發表在 痞客邦 PIXNET 留言(1) 人氣()


留言列表 (1)

發表留言
  • 曾可愛
  • 我看的懂才有鬼XD
【 X 關閉 】

【PIXNET 痞客邦】國外旅遊調查
您是我們挑選到的讀者!

填完問卷將有機會獲得心動好禮哦(注意:關閉此視窗將不再出現)

立即填寫取消