就是部落格

以下程式

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為防呆，避免指定步數超過最大步數。