Line data Source code
1 1 : program main
2 : c
3 : c convert model given in standard ASCII GONG format into
4 : c input model for adiabatic oscillation package.
5 : c
6 : c Note: it is assumed that model corresponds to a complete
7 : c model. Hence it is extended to zero radius, if that point is
8 : c not included.
9 : c
10 : c If second derivatives are not set in glob(11) and glob(12),
11 : c they are estimated from the central values of the other
12 : c model variables
13 : c
14 : c Original version: 25/7/97
15 : c
16 : implicit double precision(a-h, o-z)
17 : parameter (nnmax=5000, iconmx=30, ivarmx=50, iaa=6)
18 : character*280 fin, fout,head
19 : dimension glob(iconmx), var(ivarmx,nnmax), var1(ivarmx,nnmax),
20 : * data(8), aa(iaa,nnmax), ireset(16),q(nnmax),x(nnmax)
21 : data ireset /3,4,5,6,8,9,10,11,12,13,14,16,17,18,19,20/
22 : c
23 1 : cgrav=6.67232d-8
24 1 : pi4=16.d0*atan(1.d0)
25 : c
26 1 : write(6,*) 'Enter input GONG file'
27 1 : read(5,'(a)') fin
28 1 : write(6,*) 'Enter output amdl file'
29 1 : read(5,'(a)') fout
30 : c
31 1 : open(2,file=fin,status='old')
32 1 : open(3,file=fout,status='unknown',form='unformatted')
33 : c
34 1 : write(6,110)
35 5 : do 15 i=1,4
36 4 : read(2,'(a)') head
37 5 : 15 write(6,'(a)') head
38 : c
39 1 : read(2,120) nn, iconst, ivar, ivers
40 16 : read(2,130) (glob(i),i=1,iconst)
41 : c
42 1468 : do 20 n=1,nn
43 60148 : 20 read(2,130) (var(i,n),i=1,ivar)
44 : c
45 1 : if(var(1,1).gt.var(1,nn)) then
46 1 : nn1=nn+1
47 41 : do 30 i=1,ivar
48 58720 : do 25 n=1,nn
49 58720 : 25 var1(i,n)=var(i,nn1-n)
50 58721 : do 30 n=1,nn
51 58720 : 30 var(i,n)=var1(i,n)
52 : end if
53 : c
54 1 : if(var(1,1).gt.1.d6) then
55 41 : do 32 i=1,ivar
56 58721 : do 32 n=1,nn
57 58720 : 32 var1(i,n+1)=var(i,n)
58 : c
59 41 : do 34 i=1,ivar
60 41 : 34 var1(i,1)=0
61 : c
62 17 : do 36 ir=1,16
63 16 : i=ireset(ir)
64 17 : 36 var1(i,1)=var1(i,2)
65 : c
66 1 : nn=nn+1
67 41 : do 38 i=1,ivar
68 58761 : do 38 n=1,nn
69 58760 : 38 var(i,n)=var1(i,n)
70 : end if
71 : c
72 1469 : do 40 n=1,nn
73 1468 : q(n)=exp(var(2,n))
74 1469 : 40 x(n)=var(1,n)/glob(2)
75 : c
76 1 : x(1)=0
77 1 : q(1)=0
78 : c
79 1468 : do 45 n=2,nn
80 1467 : aa(1,n)=x(n)
81 1467 : aa(2,n)=q(n)/x(n)**3
82 : aa(3,n)=cgrav*glob(1)*q(n)*var(5,n)/
83 1467 : * (var(10,n)*var(4,n)*var(1,n))
84 1467 : aa(4,n)=var(10,n)
85 1467 : aa(5,n)=var(15,n)
86 1468 : 45 aa(6,n)=pi4*var(5,n)*var(1,n)**3/(glob(1)*q(n))
87 : c
88 1 : aa(1,1)=0
89 1 : aa(2,1)=pi4/3.d0*var(5,1)*glob(2)**3/glob(1)
90 1 : aa(3,1)=0
91 1 : aa(4,1)=var(10,1)
92 1 : aa(5,1)=0
93 1 : aa(6,1)=3.d0
94 : c
95 1 : if(aa(5,nn).le.10) then
96 0 : nn=nn-1
97 0 : write(6,*) 'Chop off outermost point'
98 : end if
99 : c
100 1 : data(1)=glob(1)
101 1 : data(2)=glob(2)
102 1 : data(3)=var(4,1)
103 1 : data(4)=var(5,1)
104 : c
105 1 : if(glob(11).lt.0.and.glob(11).gt.-10000) then
106 0 : data(5)=-glob(11)/var(10,1)
107 0 : data(6)=-glob(12)
108 : else
109 : data(5)=pi4/3.d0*cgrav*(var(5,1)*glob(2))**2/
110 1 : * (var(4,1)*var(10,1))
111 1 : d2amax=0.d0
112 308 : do 50 n=2,nn
113 308 : d2amax=max(d2amax,aa(5,n)/x(n)**2)
114 308 : if(x(n).ge.0.05d0) go to 55
115 0 : 50 continue
116 1 : 55 data(6)=d2amax+data(5)
117 1 : write(6,140) data(5), data(6)
118 : end if
119 : c
120 1 : data(7)=-1.d0
121 1 : data(8)=0.d0
122 : c
123 1 : nmodel=1
124 10285 : write(3) nmodel,nn,(data(i),i=1,8),((aa(i,n),i=1,6),n=1,nn)
125 1 : stop
126 : 110 format(/' GONG model header:'/)
127 : 120 format(4i10)
128 : 130 format(5e16.9)
129 : 140 format(' data(5), data(6) reset to ',1p2e13.5)
130 1 : end
|
