$debug
program main
dimension t(120)
character*12 fin,fou1,fou2
double precision z1,z2,z3,t1,t,pi
common r0,n,rht(10),h(10)
c data t/.0005,.0010,.0015,.0020,.0024,.0035,.0049,.0069,.0094,
c $ .0124,.0159,.0203,.0263,.0338,.0427,.0541,.0685,.0869,.1102,.1395
c $ ,.1763,.2225,.2806,.3536,.4455,.5612,.7073,.8911,1.1225,1.4141
c $ ,1.7811,2.2431,2.8242,3.5559,4.4772,5.6371,7.0969,7.9815/
PI=3.1415926535897d0
write(*,*) 'input the file name of model'
read(*,100) fin
100 format(a)
open(9,file=fin,STATUS='unknown')
2 READ(9,*) r0,N,t0,nt
DO 55 I=1,N
READ(9,*) RHT(I),dummy,H(I)
55 write(*,*) rht(i),h(i)
IF(RHT(1).LE.0.0) GO TO 2
CLOSE(9)
write(*,*) 'input the name of each outfile'
read(*,100) fou1
OPEN(11,FILE=fou1,STATUS='UNKNOWN')
c OPEN(12,FILE=fou2,STATUS='UNKNOWN')
c
dt=0.1
t(1)=t0
tl0=alog10(t0)
do i=2,nt
t(i)=10**(tl0+(i-1)*dt)
enddo
ICT=1
DO 10 I=1,nt
t1=t(i)
c X=DSQRT(PI*2.D7*T1*RHT(1))/DBLE(H(1))
CALL RES(T1,Z1,Z2,Z3,ICT)
ICT=2
WRITE(11,40) t1,z1,z2,z3
10 CONTINUE
40 FORMAT(4(1x,e12.5))
CLOSE(11)
END
C***********************************************
SUBROUTINE RES(T1,Z1,Z2,Z3,ICT)
DOUBLE PRECISION T1,Z1,Z2,Z3,RX,M,P,Q,PI,EX1,HZ,HZ0,S1,V1,V2,amu
DOUBLE PRECISION M1(144),H1(144),K(12),AA1(90),AA2(54)
COMMON R0,N,RHT1(10),H(10)
data AA1/
$1.1667d-14,1.800256d-14,2.930619d-14,4.522038d-14,7.36138d-14,
$1.135885d-13,1.849096d-13,2.853213d-13,4.644718d-13,7.166947d-13,
$1.1667d-12,1.800256d-12,2.930619d-12,4.522038d-12,7.36138d-12,
$1.135885d-11,1.849096d-11,2.853213d-11,4.644718d-11,7.166947d-11,
$1.1667d-10,1.800256d-10,2.930619d-10,4.522038d-10,7.36138d-10,
$1.135885d-9,1.849096d-9,2.853213d-9,4.644718d-9,7.166947d-9,
$1.166700d-8,1.800256d-8,2.930619d-8,4.522038d-8,7.361381d-8,
$1.135885d-7,1.849095d-7,2.853213d-7,4.644716d-7,7.166942d-7,
&1.166699d-6,1.800252d-6,2.930610d-6,4.522015d-6,7.361325d-6,
&1.135870d-5,1.849059d-5,2.853122d-5,4.644490d-5,7.166375d-5,
&1.166557d-4,1.799894d-4,2.929711d-4,4.519758d-4,7.355654d-4,
&1.134446d-3,1.845483d-3,2.844143d-3,4.621947d-3,7.109806d-3,
&1.152369d-2,1.764345d-2,2.840762d-2,4.297706d-2,6.803326d-2,
&9.978459d-2,1.510705d-1,2.035406d-1,2.712354d-1,2.760739d-1,
&2.166920d-1,-7.837237d-2,-3.406756d-1,-3.606937d-1,5.130245d-1,
&-5.947247d-2,-1.951171d-1,1.992356d-1,-1.385216d-1,8.793209d-2,
&-5.506971d-2,3.456378d-2,-2.175272d-2,1.371003d-2,-8.646564d-3,
&5.454628d-3,-3.441389d-3,2.171307d-3,-1.369986d-3,8.643990d-4/
data AA2/
&-5.453979d-4,3.441225d-4,-2.171266d-4,1.369976d-4,-8.643964d-5,
&5.453972d-5,-3.441224d-5,2.171265d-5,-1.369976d-5,8.643963d-6,
&-5.453972d-6,3.441224d-6,-2.171265d-6,1.369976d-6,-8.643964d-7,
&5.453972d-7,-3.441224d-7,2.171265d-7,-1.369976d-7,8.643963d-8,
&-5.453972d-8,3.441224d-8,-2.171265d-8,1.369976d-8,-8.643963d-9,
program main
dimension t(120)
character*12 fin,fou1,fou2
double precision z1,z2,z3,t1,t,pi
common r0,n,rht(10),h(10)
c data t/.0005,.0010,.0015,.0020,.0024,.0035,.0049,.0069,.0094,
c $ .0124,.0159,.0203,.0263,.0338,.0427,.0541,.0685,.0869,.1102,.1395
c $ ,.1763,.2225,.2806,.3536,.4455,.5612,.7073,.8911,1.1225,1.4141
c $ ,1.7811,2.2431,2.8242,3.5559,4.4772,5.6371,7.0969,7.9815/
PI=3.1415926535897d0
write(*,*) 'input the file name of model'
read(*,100) fin
100 format(a)
open(9,file=fin,STATUS='unknown')
2 READ(9,*) r0,N,t0,nt
DO 55 I=1,N
READ(9,*) RHT(I),dummy,H(I)
55 write(*,*) rht(i),h(i)
IF(RHT(1).LE.0.0) GO TO 2
CLOSE(9)
write(*,*) 'input the name of each outfile'
read(*,100) fou1
OPEN(11,FILE=fou1,STATUS='UNKNOWN')
c OPEN(12,FILE=fou2,STATUS='UNKNOWN')
c
dt=0.1
t(1)=t0
tl0=alog10(t0)
do i=2,nt
t(i)=10**(tl0+(i-1)*dt)
enddo
ICT=1
DO 10 I=1,nt
t1=t(i)
c X=DSQRT(PI*2.D7*T1*RHT(1))/DBLE(H(1))
CALL RES(T1,Z1,Z2,Z3,ICT)
ICT=2
WRITE(11,40) t1,z1,z2,z3
10 CONTINUE
40 FORMAT(4(1x,e12.5))
CLOSE(11)
END
C***********************************************
SUBROUTINE RES(T1,Z1,Z2,Z3,ICT)
DOUBLE PRECISION T1,Z1,Z2,Z3,RX,M,P,Q,PI,EX1,HZ,HZ0,S1,V1,V2,amu
DOUBLE PRECISION M1(144),H1(144),K(12),AA1(90),AA2(54)
COMMON R0,N,RHT1(10),H(10)
data AA1/
$1.1667d-14,1.800256d-14,2.930619d-14,4.522038d-14,7.36138d-14,
$1.135885d-13,1.849096d-13,2.853213d-13,4.644718d-13,7.166947d-13,
$1.1667d-12,1.800256d-12,2.930619d-12,4.522038d-12,7.36138d-12,
$1.135885d-11,1.849096d-11,2.853213d-11,4.644718d-11,7.166947d-11,
$1.1667d-10,1.800256d-10,2.930619d-10,4.522038d-10,7.36138d-10,
$1.135885d-9,1.849096d-9,2.853213d-9,4.644718d-9,7.166947d-9,
$1.166700d-8,1.800256d-8,2.930619d-8,4.522038d-8,7.361381d-8,
$1.135885d-7,1.849095d-7,2.853213d-7,4.644716d-7,7.166942d-7,
&1.166699d-6,1.800252d-6,2.930610d-6,4.522015d-6,7.361325d-6,
&1.135870d-5,1.849059d-5,2.853122d-5,4.644490d-5,7.166375d-5,
&1.166557d-4,1.799894d-4,2.929711d-4,4.519758d-4,7.355654d-4,
&1.134446d-3,1.845483d-3,2.844143d-3,4.621947d-3,7.109806d-3,
&1.152369d-2,1.764345d-2,2.840762d-2,4.297706d-2,6.803326d-2,
&9.978459d-2,1.510705d-1,2.035406d-1,2.712354d-1,2.760739d-1,
&2.166920d-1,-7.837237d-2,-3.406756d-1,-3.606937d-1,5.130245d-1,
&-5.947247d-2,-1.951171d-1,1.992356d-1,-1.385216d-1,8.793209d-2,
&-5.506971d-2,3.456378d-2,-2.175272d-2,1.371003d-2,-8.646564d-3,
&5.454628d-3,-3.441389d-3,2.171307d-3,-1.369986d-3,8.643990d-4/
data AA2/
&-5.453979d-4,3.441225d-4,-2.171266d-4,1.369976d-4,-8.643964d-5,
&5.453972d-5,-3.441224d-5,2.171265d-5,-1.369976d-5,8.643963d-6,
&-5.453972d-6,3.441224d-6,-2.171265d-6,1.369976d-6,-8.643964d-7,
&5.453972d-7,-3.441224d-7,2.171265d-7,-1.369976d-7,8.643963d-8,
&-5.453972d-8,3.441224d-8,-2.171265d-8,1.369976d-8,-8.643963d-9,
