aboutsummaryrefslogtreecommitdiff
path: root/Chabund.f
blob: dcc8d2df050be2f76ae9439bac5addbe80ca6074 (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253

      subroutine chabund
c*****************************************************************************
c  this routine changes the input abundances according to the users request
c  and the re-runs the synthesis with the new abundances OVERWRITING the 
c  original output files, but using the same model, linelist, etc
c*****************************************************************************

      implicit real*8 (a-h,o-z)
      include 'Atmos.com'
      include 'Factor.com'
      include 'Mol.com'
      include 'Pstuff.com'
      include 'Multistar.com'
      real*8 xnum
      real*4 shortnum
      character choice2


c*****set the local temporary parameters to what they were for the last synth
      newnumiso = numiso
      newnumisosyn = numisosyn
      newnumpecatom = numpecatom
      newnumatomsyn = numatomsyn
      if (newnumatomsyn == 0) newnumatomsyn = 1
      do i=3,95
         newpec(i) = pec(i)
         do j=1,newnumatomsyn
            newpecabund(i,j) = pecabund(i,j)
         enddo
      enddo
      do i=1,newnumiso
         if(newnumiso >= 1) then
            newisotope(i)=isotope(i)
            do j=1,newnumisosyn
               newisoabund(i,j)=isoabund(i,j)
            enddo
         endif
      enddo
      if (newnumatomsyn <= 0) newnumatomsyn = 1


c*****SPECIAL CASE
c     for the "binary" driver, present a limited set of options; only
c     the abundances of the elements already chosen for variation in
c     the parameter file can be altered
11    if (control == 'binary ') then
         istat = ivcleof(4,1)
         if (syncount == 1) then
            array = 'Abundance Alterations for the PRIMARY STAR:   '
         else
            array = 'Abundance Alterations for the SECONDARY STAR: '
         endif
         istat = ivwrite(5,5,array,46)
         array = 'Atom   Abundance (relative to initial model)'
         istat = ivwrite(6,1,array,49)
         line = 6
         do i=3,95
            if(newpec(i) == 1) then
               line = line + 1
               write (array,1001) i,(newpecabund(i,j),j=1,newnumatomsyn)
               istat = ivwrite(line,1,array,57)
            endif
         enddo
         array = 'Options: c = change abundances '
         istat = ivwrite(line+1,1,array,41)
         array = '         x = exit no changes    q = rerun synth'
         istat = ivwrite(line+2,1,array,47)
         array = 'What is your choice? '
         nchars = 21
         call getasci (nchars,line+4)
         choice = chinfo(1:1)
         if (choice/='c' .and. choice/='x' .and.
     .       choice/='q') go to 11
      endif
         
        
c*****for the "synth" driver, present the user with the current abundance 
c     alterations, and the options; find out what is desired
1     if (control == 'synth'  .or.
     .    control == 'isotop' .or.
     .    control == 'synpop') then
         line = 4
         istat = ivcleof(line,1)
         write (array,1006)
         line = line + 1
         istat = ivwrite(6,1,array,72)
         line = line + 1
         do i=3,95
            if(newpec(i) == 1) then
               line = line + 1
               write (array,1001) i,(newpecabund(i,j),j=1,newnumatomsyn)
               istat = ivwrite(line,1,array,57)
            endif
         enddo
         if (newnumiso > 0) then
            do i=1,newnumiso
               line = line + 1
               write (array,1002) i, newisotope(i), 
     .                            (newisoabund(i,j),j=1,newnumisosyn)
               istat = ivwrite(line,1,array,68)
            enddo
         endif
         write (array,1007)
         istat = ivwrite(line+1,1,array,61)
         write (array,1008)
         istat = ivwrite(line+2,1,array,72)
         array = 'What is your choice? '
         nchars = 21
         call getasci (nchars,line+3)
         choice = chinfo(1:1)
         if (choice/='c' .and. choice/='i' .and.
     .       choice/='n' .and. choice/='x' .and.
     .       choice/='q') go to 1
      endif

         
c*****for "synth", change elemental abundances
      if     (choice == 'c') then
20       istat = ivcleof(line+1,1)
         array = 'Which element to change? '
         nchars = 25
         call getnum (nchars,line+1,xnum,shortnum)
         if (xnum<2.0 .or. xnum>95.) go to 20
         j = nint(xnum)
         if (newpec(j) == 1) then
25          array = 'n = new abundances, or z = zero offsets? '
            nchars = 41
            call getasci (nchars,line+2)
            choice2 = chinfo(1:1)
            if     (choice2 == 'z') then
               newpec(j) = 0
               do i=1,5
                  newpecabund(j,i) = 0.0
               enddo
               newnumpecatom = newnumpecatom - 1
            elseif (choice2 == 'n') then
               write (array,1004)
               istat = ivwrite(line+3,1,array,40)
               read (*,*) (newpecabund(j,i),i=1,newnumatomsyn)
            else
               go to 25
            endif
         else
            newpec(j) = 1
            write (array,1004)
            istat = ivwrite(line+2,1,array,40)
            read (*,*) (newpecabund(j,i),i=1,newnumatomsyn) 
            newnumpecatom = newnumpecatom + 1
         endif


c*****for "synth", change isotopic factors
      elseif (choice == 'i') then
         istat = ivcleof(line+1,1)
         array = 'Options: c = change an isotopic factor'
         istat = ivwrite(line+1,1,array,49)
         array = '         n = enter a new isotope'
         istat = ivwrite(line+2,1,array,48)
30       array = 'What is your choice? '
         nchars = 22
         call getasci (nchars,line+3)
         choice2 = chinfo(1:1)
         if (choice2 == 'c') then
35          istat = ivcleof(line+1,1)
            array =  'Which isotope number from the list? '
            nchars = 36
            call getnum (nchars,line+1,xnum,shortnum)
            j = nint(xnum)
            if (j<1 .or. j>newnumiso) go to 35
            istat = ivcleof(line+2,1)
            array = 'What are the new division factors? '
            istat = ivwrite(line+2,1,array,35)
            read (*,*) (newisoabund(j,i),i=1,newnumisosyn)
         elseif (choice2 == 'n') then
            newnumiso = newnumiso + 1
            istat = ivcleof(line+1,1)
            array = 'What is the new isotope designation? '
            nchars = 37
            call getnum (nchars,line+4,xnum,shortnum)
            newisotope(newnumiso) = xnum
            array = 'What are its division factors? '
            istat = ivwrite(line+2,1,array,31)
            read (*,*) (newisoabund(newnumiso,i),i=1,newnumisosyn)
         else
            go to 30
         endif
            

c*****for "synth", change the number of syntheses
      elseif (choice == 'n') then
55       array = 'How many synths? '
         nchars = 17
         call getnum (nchars,line+5,xnum,shortnum)
         if (xnum > 5.) go to 55
         newnumatomsyn = nint(xnum)
         go to 1


c*****for "synth", exit the routine without changing anything
      elseif (choice == 'x') then
         return


c*****for "synth", make the proposed alterations permanent; then
c     return to the calling routine, which allegedly will
c     redo the syntheses.
      elseif (choice == 'q') then
         numiso = newnumiso
         numisosyn = newnumisosyn
         numpecatom = newnumpecatom
         numatomsyn = newnumatomsyn
         do i=3,95
            pec(i) = newpec(i)
            do j=1,numatomsyn
               pecabund(i,j) = newpecabund(i,j)
            enddo
         enddo
         do i=1,numiso
            if(numiso >= 1) then
               isotope(i) = newisotope(i)
               do j=1,numisosyn
                  isoabund(i,j)=newisoabund(i,j)
               enddo
            endif
         enddo
         return
      endif


c*****loop back and print out the main menu again
      if (control == 'synth  ') then
         go to 1
      else
         go to 11
      endif


c*****format statements
1001  format (i3, 4x, 5(f8.3, 2x))
1002  format (i2, 2x, f10.5, 5(3x, f8.3))
1004  format ('Enter the new offsets on the line below:')
1006     format ('element, abundance offsets   OR   ',
     .            'isotope number, isotope name, factors')
1007     format ('Options: c = change abundance       ',
     .           'i = change isotopic ratio')
1008     format('         n = change # syntheses     q = rerun ',
     .          'syntheses        x = exit')

      end