這是一個剪輯由G03計算BOMD後fchk檔案剪出座標檔案的小程式。
由於G03的fchk都是用5行(column)在存數據，可是座標卻是X, Y, Z讀取。
因此必須做轉置動作轉存成3行的數據，例如
C11 C12 C13 C14 C15
C21 C22 C23 C24 C25
˙
˙
˙
轉置成
C11 C12 C13
C14 C15 C21
C22 C23 C24
˙
˙
˙
以下程式
--
program collecting_geom_structure
! This program collectes the coordinate of each atom at each step of a
trajectory.
character(16) :: filename,Ofilename
logical :: file_exist
integer :: Nstep,Natom,Ndata,i,j,k,l,m
real :: Nrow
character(16) ::
C1(10000),C2(10000),C3(10000),C4(10000),C5(10000),C(50000),Cr(3)
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 *, "Please enter a filename for output file"
read(*,*) Ofilename
print *, "Enter the number of atoms in your job"
read(*,*) Natom
print *, "Enter the maximum step No. in your trajectory job"
print *, "Warning : the product of No. atoms and No. steps
should be less than 16666"
read(*,*) Nstep
Ndata = Natom*3*Nstep
Nrow = Ndata/5
if (mod(Ndata,5) .eq. 0) then
Nrow = int(Nrow)
else
Nrow = int(Nrow) + 1
end if
open (unit=7,file=filename,status='old',action='read')
open (unit=8,file=Ofilename,status='new')
do i = 1, Nrow
read(7,2,advance='yes',end=4) C1(i), C2(i), C3(i),
C4(i), C5(i)
j = (i-1)*5
k = j + 1
C(k) = C1(i)
C(k+1) = C2(i)
C(k+2) = C3(i)
C(k+3) = C4(i)
C(k+4) = C5(i)
end do
do l = 1, Ndata
m = mod(l,3)
Cr(m) = C(l)
if (m .eq. 0) then
write(8,3,advance='yes') Cr(1),Cr(2),Cr(0)
end if
end do
2 format (5a16)
3 format (3(1x,a16))
else
print *, filename," doesn't exist"
go to 999
end if
4 print *, "The coordinate of each atom at each step on the time
evolution were written in ", Ofilename, ". Bye."
999 stop
end
--
上色部份為轉置功能的主要程式語法，我的程式設計概念是讀取和寫入分開。
紅色部份是利用fortran的讀取順序，由左而右由上而下將所有數據依序存成變數陣列。
變數陣列為C(Ndata)
讀入電腦記憶體後，黃色部份則是將陣列依照3的倍數依序寫入檔案。
因每三筆數據須換行，因此利用陣列的引數(index) = 1~Ndata做回圈依序除以3，可得有
序結果1,4,7,....餘數為1，2,5,8,....餘數為2，3,6,9,....餘數為0
因此將變數陣列C(Ndata)依序改成Cr(1,2,0)再寫出則為結果。
