[Fortran] 纯文本查看 复制代码
PROGRAM test
IMPLICIT NONE
REAL, ALLOCATABLE :: a(:,:),b(:,:),hhs(:)
INTEGER :: clock_rate, clock_max, clock_1, clock_2
INTEGER :: i,ios,z, j, k, tn, num,numm, step,m,n
REAL :: summ,up,below ,cos2
real :: cos22,r1,q1,w1
character (160):: fn1, fn2, filename
Real *8 time
Integer *4 time0, time1, dtime
Real :: summc , avg , sc
fn1='x'
fn2='y'
num=100
numm=num*num
tn=1000
ALLOCATE(a(1,num),b(1,num),ss(num),rr(num),hh(num),hhs(numm))
Print*, 'Running...'
Call system_clock(time0)
OPEN(10, FILE='output.dat',status='unknown' )
Open(101,File='IDlist.dat', Status='old') ! contains 1000 filenames, eg: 1,2,3,..1000
Do i=1,tn
Read (101, *, Iostat=ios) filename
If (ios/=0) Exit
Open(200,File="D:\x\"//trim(adjustl(fn1))//trim(adjustl(filename))//'.dat',status='old')
READ(200,*,IOSTAT=ios) a
If (ios/=0) Exit
Open(300,File="D:\y\"//trim(adjustl(fn2))//trim(adjustl(filename))//'.dat',status='old')
READ(300,*,IOSTAT=ios) b
If (ios/=0) Exit
do z = 1, numm
do j = 1, num
do k = 1, num
r1 = 20
q1 = b(1,j)-a(1,k)
cos2 = (r1*q1+r1*r1+q1*q1)**2/3
w1=sqrt(q1*q1+r1*r1)
enddo
enddo
if (w1.GT.0.and.w1.LE.20) then
hhs(z)=cos22/w1/w1/w1
endif
summ=sum(hhs)
enddo
WRITE(10,'(1X,F12.6)') summ
enddo
DEALLOCATE(a,b,ss,rr,hh)
close(300)
close(200)
close(101)
close(10)
Call system_clock(time1, dtime)
time = 1D0*(time1-time0)/dtime
Write (*, '(a7,f16.7)') 'Time = ', time
END PROGRAM
[Fortran] 纯文本查看 复制代码
PROGRAM test
IMPLICIT NONE
REAL, ALLOCATABLE :: a(:),b(:),hhs(:)
INTEGER :: i, z, j, k, tn, num, numm
REAL :: summ, cos2, r1, q1, w1
character(160):: filename
Real(8) time
Integer(4) time0, time1, dtime
num=100
numm=num*num
tn=1000
ALLOCATE(a(num),b(num),hhs(numm))
Print*, 'Running...'
Call system_clock(time0)
OPEN(10, FILE='output.dat',status='unknown' )
Do i=1,tn
write(filename,'(i6)') i
Open(200,File="D:\x\"//trim(adjustl(filename))//'.dat',status='old')
READ(200,*) a
close(200)
Open(300,File="D:\y\"//trim(adjustl(filename))//'.dat',status='old')
READ(300,*) b
close(300)
z=0
summ=0.0
r1 = 20
do j = 1, num
do k = 1, num
z=z+1
q1 = b(j)-a(k)
cos2 = (r1*q1+r1*r1+q1*q1)**2/3
w1=sqrt(q1*q1+r1*r1)
if (w1.GT.0.and.w1.LE.20) then
hhs(z)=cos2/(w1*w1*w1)
endif
enddo
enddo
summ=sum(hhs)
WRITE(10,'(1X,F12.6)') summ
enddo
DEALLOCATE(a,b)
close(10)
Call system_clock(time1, dtime)
time = 1D0*(time1-time0)/dtime
Write (*, '(a7,f16.7)') 'Time = ', time
END PROGRAM