1 | c======================================================================= |
---|
2 | SUBROUTINE tabfi(ngrid,nid,Lmodif,tab0,day_ini,lmax,p_rad, |
---|
3 | . p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time) |
---|
4 | c======================================================================= |
---|
5 | c |
---|
6 | c C. Hourdin 15/11/96 |
---|
7 | c |
---|
8 | c Object: Lecture du tab_cntrl physique dans un fichier |
---|
9 | c ------ et initialisation des constantes physiques |
---|
10 | c |
---|
11 | c Arguments: |
---|
12 | c ---------- |
---|
13 | c |
---|
14 | c Inputs: |
---|
15 | c ------ |
---|
16 | c |
---|
17 | c - nid: unitne logique du fichier ou on va lire le tab_cntrl |
---|
18 | c (ouvert dans le programme appellant) |
---|
19 | c |
---|
20 | c si nid=0: |
---|
21 | c pas de lecture du tab_cntrl mais |
---|
22 | c Valeurs par default des constantes physiques |
---|
23 | c |
---|
24 | c - tab0: Offset de tab_cntrl a partir duquel sont ranges |
---|
25 | c les parametres physiques (50 pour start_archive) |
---|
26 | c |
---|
27 | c - Lmodif: si on souhaite modifier les constantes Lmodif = 1 = TRUE |
---|
28 | c |
---|
29 | c |
---|
30 | c Outputs: |
---|
31 | c -------- |
---|
32 | c |
---|
33 | c - day_ini: tab_cntrl(tab0+3) (Dans les cas ou l'on souhaite |
---|
34 | c comparer avec le day_ini dynamique) |
---|
35 | c |
---|
36 | c - lmax: tab_cntrl(tab0+2) (pour test avec nlayer) |
---|
37 | c |
---|
38 | c - p_rad |
---|
39 | c - p_omeg ! |
---|
40 | c - p_g ! Constantes physiques ayant des |
---|
41 | c - p_mugaz ! homonymes dynamiques |
---|
42 | c - p_daysec ! |
---|
43 | c |
---|
44 | c======================================================================= |
---|
45 | ! to use 'getin' |
---|
46 | use ioipsl_getincom_p , only: getin_p |
---|
47 | |
---|
48 | use surfdat_h, only: albedice, emisice, iceradius, dtemisice, |
---|
49 | & emissiv |
---|
50 | use comsoil_h, only: volcapa |
---|
51 | use iostart, only: get_var |
---|
52 | use mod_phys_lmdz_para, only: is_parallel |
---|
53 | use planete_mod, only: year_day, periastr, apoastr, peri_day, |
---|
54 | & obliquit, z0, lmixmin, emin_turb |
---|
55 | |
---|
56 | implicit none |
---|
57 | |
---|
58 | #include "comcstfi.h" |
---|
59 | !#include "planete.h" |
---|
60 | #include "netcdf.inc" |
---|
61 | #include "callkeys.h" |
---|
62 | |
---|
63 | c----------------------------------------------------------------------- |
---|
64 | c Declarations |
---|
65 | c----------------------------------------------------------------------- |
---|
66 | |
---|
67 | c Arguments |
---|
68 | c --------- |
---|
69 | INTEGER,INTENT(IN) :: ngrid,nid,tab0 |
---|
70 | INTEGER*4,INTENT(OUT) :: day_ini |
---|
71 | INTEGER,INTENT(IN) :: Lmodif |
---|
72 | INTEGER,INTENT(OUT) :: lmax |
---|
73 | REAL,INTENT(OUT) :: p_rad,p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time |
---|
74 | |
---|
75 | c Variables |
---|
76 | c --------- |
---|
77 | INTEGER,PARAMETER :: length=100 |
---|
78 | REAL tab_cntrl(length) ! array in which are stored the run's parameters |
---|
79 | INTEGER ierr,nvarid |
---|
80 | INTEGER size |
---|
81 | CHARACTER modif*20 |
---|
82 | LOGICAL :: found |
---|
83 | |
---|
84 | write(*,*)"tabfi: nid=",nid," tab0=",tab0," Lmodif=",Lmodif |
---|
85 | |
---|
86 | IF (nid.eq.0) then |
---|
87 | c----------------------------------------------------------------------- |
---|
88 | c Initialization of various physical constants to defaut values (nid = 0 case) |
---|
89 | c----------------------------------------------------------------------- |
---|
90 | ! Ehouarn: Default Saturn values: |
---|
91 | tab_cntrl(:)=0 |
---|
92 | write(*,*) "Using default Saturn values..." |
---|
93 | ! these should be read in a def file I guess... |
---|
94 | lmax=0 ! not used anyways |
---|
95 | !day_ini=0 |
---|
96 | time=0 |
---|
97 | ! radius of the planet |
---|
98 | rad=60268000 |
---|
99 | call getin_p("radius",rad) |
---|
100 | ! Planetary rotation rate |
---|
101 | omeg=0.00016512100410182 |
---|
102 | call getin_p("omega",omeg) |
---|
103 | ! Gravity |
---|
104 | g=10.44 |
---|
105 | call getin_p("g",g) |
---|
106 | !mugaz=2.34 !EM: does not give cpp=11500) |
---|
107 | mugaz=2.53 ! with this value of mugaz, cpp=11500 |
---|
108 | call getin_p("mugaz",mugaz) |
---|
109 | ! recompute r=R/mugaz to be consitent |
---|
110 | r=8.3144621/(mugaz*1.e-3) |
---|
111 | ! kappa |
---|
112 | rcp=0.2857143 |
---|
113 | call getin_p("kappa",rcp) |
---|
114 | cpp=(8.3144621/(mugaz/1000.0))/rcp |
---|
115 | call getin_p("cpp",cpp) |
---|
116 | ! write(*,*) "tabfi: cpp=",cpp |
---|
117 | ! length (s) of a "standard" day |
---|
118 | daysec=38052 |
---|
119 | call getin_p("day_length",daysec) |
---|
120 | ! physics time step (s) ! not sure we need this here |
---|
121 | dtphys=19026 |
---|
122 | ! length of year, in standard days |
---|
123 | year_day=24430 |
---|
124 | ! Orbital parameters |
---|
125 | periastr=9.02151966094971 |
---|
126 | call getin_p("periastron",periastr) |
---|
127 | apoastr=10.054479598999 |
---|
128 | call getin_p("apoastron",apoastr) |
---|
129 | peri_day=19280 |
---|
130 | call getin_p("periastron_day",peri_day) |
---|
131 | obliquit=26.7299995422363 |
---|
132 | call getin_p("obliquity",obliquit) |
---|
133 | ! Other parameters some physical parametrizations need |
---|
134 | z0=1e-2 |
---|
135 | lmixmin=30 |
---|
136 | emin_turb=1.e-6 |
---|
137 | albedice(:)=0 |
---|
138 | emisice(:)=0 |
---|
139 | emissiv=0 |
---|
140 | iceradius(:)=1.e-6 |
---|
141 | dtemisice(:)=0 |
---|
142 | volcapa=1000000 |
---|
143 | c----------------------------------------------------------------------- |
---|
144 | c Save some constants for later use (as routine arguments) |
---|
145 | c----------------------------------------------------------------------- |
---|
146 | p_omeg = omeg |
---|
147 | p_g = g |
---|
148 | p_cpp = cpp |
---|
149 | p_mugaz = mugaz |
---|
150 | p_daysec = daysec |
---|
151 | p_rad=rad |
---|
152 | |
---|
153 | ELSE |
---|
154 | c----------------------------------------------------------------------- |
---|
155 | c Initialization of physical constants by reading array tab_cntrl(:) |
---|
156 | c which contains these parameters (nid != 0 case) |
---|
157 | c----------------------------------------------------------------------- |
---|
158 | c Read 'controle' array |
---|
159 | c |
---|
160 | |
---|
161 | call get_var("controle",tab_cntrl,found) |
---|
162 | if (.not.found) then |
---|
163 | write(*,*)"tabfi: Failed reading <controle> array" |
---|
164 | call abort |
---|
165 | else |
---|
166 | write(*,*)'tabfi: tab_cntrl',tab_cntrl |
---|
167 | endif |
---|
168 | c |
---|
169 | c Initialization of some physical constants |
---|
170 | c informations on physics grid |
---|
171 | ! if(ngrid.ne.tab_cntrl(tab0+1)) then |
---|
172 | ! print*,'tabfi: WARNING !!! tab_cntrl(tab0+1).ne.ngrid' |
---|
173 | ! print*,tab_cntrl(tab0+1),ngrid |
---|
174 | ! endif |
---|
175 | lmax = nint(tab_cntrl(tab0+2)) |
---|
176 | day_ini = tab_cntrl(tab0+3) |
---|
177 | time = tab_cntrl(tab0+4) |
---|
178 | write (*,*) 'IN tabfi day_ini=',day_ini |
---|
179 | c Informations about planet for dynamics and physics |
---|
180 | rad = tab_cntrl(tab0+5) |
---|
181 | omeg = tab_cntrl(tab0+6) |
---|
182 | g = tab_cntrl(tab0+7) |
---|
183 | mugaz = tab_cntrl(tab0+8) |
---|
184 | rcp = tab_cntrl(tab0+9) |
---|
185 | cpp=(8.314511/(mugaz/1000.0))/rcp |
---|
186 | daysec = tab_cntrl(tab0+10) |
---|
187 | dtphys = tab_cntrl(tab0+11) |
---|
188 | c Informations about planet for the physics only |
---|
189 | year_day = tab_cntrl(tab0+14) |
---|
190 | periastr = tab_cntrl(tab0+15) |
---|
191 | apoastr = tab_cntrl(tab0+16) |
---|
192 | peri_day = tab_cntrl(tab0+17) |
---|
193 | obliquit = tab_cntrl(tab0+18) |
---|
194 | c boundary layer and turbeulence |
---|
195 | z0 = tab_cntrl(tab0+19) |
---|
196 | lmixmin = tab_cntrl(tab0+20) |
---|
197 | emin_turb = tab_cntrl(tab0+21) |
---|
198 | c optical properties of polar caps and ground emissivity |
---|
199 | albedice(1)= tab_cntrl(tab0+22) |
---|
200 | albedice(2)= tab_cntrl(tab0+23) |
---|
201 | emisice(1) = tab_cntrl(tab0+24) |
---|
202 | emisice(2) = tab_cntrl(tab0+25) |
---|
203 | emissiv = tab_cntrl(tab0+26) |
---|
204 | iceradius(1)= tab_cntrl(tab0+31) ! mean scat radius of CO2 snow (north) |
---|
205 | iceradius(2)= tab_cntrl(tab0+32) ! mean scat radius of CO2 snow (south) |
---|
206 | dtemisice(1)= tab_cntrl(tab0+33) !time scale for snow metamorphism (north) |
---|
207 | dtemisice(2)= tab_cntrl(tab0+34) !time scale for snow metamorphism (south) |
---|
208 | c soil properties |
---|
209 | volcapa = tab_cntrl(tab0+35) ! volumetric heat capacity |
---|
210 | c----------------------------------------------------------------------- |
---|
211 | c Save some constants for later use (as routine arguments) |
---|
212 | c----------------------------------------------------------------------- |
---|
213 | p_omeg = omeg |
---|
214 | p_g = g |
---|
215 | p_cpp = cpp |
---|
216 | p_mugaz = mugaz |
---|
217 | p_daysec = daysec |
---|
218 | p_rad=rad |
---|
219 | |
---|
220 | ENDIF ! end of (nid = 0) |
---|
221 | |
---|
222 | c----------------------------------------------------------------------- |
---|
223 | c Write physical constants to output before modifying them |
---|
224 | c----------------------------------------------------------------------- |
---|
225 | |
---|
226 | 6 FORMAT(a20,e15.6,e15.6) |
---|
227 | 5 FORMAT(a20,f12.2,f12.2) |
---|
228 | |
---|
229 | write(*,*) '*****************************************************' |
---|
230 | write(*,*) 'Reading tab_cntrl when calling tabfi before changes' |
---|
231 | write(*,*) '*****************************************************' |
---|
232 | write(*,5) '(1) = ngrid?',tab_cntrl(tab0+1),float(ngrid) |
---|
233 | write(*,5) '(2) lmax',tab_cntrl(tab0+2),float(lmax) |
---|
234 | write(*,5) '(3) day_ini',tab_cntrl(tab0+3),float(day_ini) |
---|
235 | write(*,5) '(5) rad',tab_cntrl(tab0+5),rad |
---|
236 | write(*,5) '(10) daysec',tab_cntrl(tab0+10),daysec |
---|
237 | write(*,6) '(6) omeg',tab_cntrl(tab0+6),omeg |
---|
238 | write(*,5) '(7) g',tab_cntrl(tab0+7),g |
---|
239 | write(*,5) '(8) mugaz',tab_cntrl(tab0+8),mugaz |
---|
240 | write(*,5) '(9) rcp',tab_cntrl(tab0+9),rcp |
---|
241 | write(*,6) '(11) dtphys?',tab_cntrl(tab0+11),dtphys |
---|
242 | |
---|
243 | write(*,5) '(14) year_day',tab_cntrl(tab0+14),year_day |
---|
244 | write(*,5) '(15) periastr',tab_cntrl(tab0+15),periastr |
---|
245 | write(*,5) '(16) apoastr',tab_cntrl(tab0+16),apoastr |
---|
246 | write(*,5) '(17) peri_day',tab_cntrl(tab0+17),peri_day |
---|
247 | write(*,5) '(18) obliquit',tab_cntrl(tab0+18),obliquit |
---|
248 | |
---|
249 | write(*,6) '(19) z0',tab_cntrl(tab0+19),z0 |
---|
250 | write(*,6) '(21) emin_turb',tab_cntrl(tab0+21),emin_turb |
---|
251 | write(*,5) '(20) lmixmin',tab_cntrl(tab0+20),lmixmin |
---|
252 | |
---|
253 | write(*,5) '(26) emissiv',tab_cntrl(tab0+26),emissiv |
---|
254 | write(*,5) '(24) emisice(1)',tab_cntrl(tab0+24),emisice(1) |
---|
255 | write(*,5) '(25) emisice(2)',tab_cntrl(tab0+25),emisice(2) |
---|
256 | write(*,5) '(22) albedice(1)',tab_cntrl(tab0+22),albedice(1) |
---|
257 | write(*,5) '(23) albedice(2)',tab_cntrl(tab0+23),albedice(2) |
---|
258 | write(*,6) '(31) iceradius(1)',tab_cntrl(tab0+31),iceradius(1) |
---|
259 | write(*,6) '(32) iceradius(2)',tab_cntrl(tab0+32),iceradius(2) |
---|
260 | write(*,5) '(33) dtemisice(1)',tab_cntrl(tab0+33),dtemisice(1) |
---|
261 | write(*,5) '(34) dtemisice(2)',tab_cntrl(tab0+34),dtemisice(2) |
---|
262 | |
---|
263 | write(*,5) '(35) volcapa',tab_cntrl(tab0+35),volcapa |
---|
264 | |
---|
265 | write(*,*) |
---|
266 | write(*,*) 'Lmodif in tabfi!!!!!!!',Lmodif |
---|
267 | |
---|
268 | c----------------------------------------------------------------------- |
---|
269 | c Modifications... |
---|
270 | ! NB: Modifying controls should only be done by newstart, and in seq mode |
---|
271 | if ((Lmodif.eq.1).and.is_parallel) then |
---|
272 | write(*,*) "tabfi: Error modifying tab_control should", |
---|
273 | & " only happen in serial mode (eg: by newstart)" |
---|
274 | stop |
---|
275 | endif |
---|
276 | c----------------------------------------------------------------------- |
---|
277 | |
---|
278 | IF(Lmodif.eq.1) then |
---|
279 | |
---|
280 | write(*,*) |
---|
281 | write(*,*) 'Change values in tab_cntrl ? :' |
---|
282 | write(*,*) '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' |
---|
283 | write(*,*) '(Current values given above)' |
---|
284 | write(*,*) |
---|
285 | write(*,*) '(3) day_ini : Initial day (=0 at Ls=0)' |
---|
286 | write(*,*) '(19) z0 : surface roughness (m)' |
---|
287 | write(*,*) '(21) emin_turb : minimal energy (PBL)' |
---|
288 | write(*,*) '(20) lmixmin : mixing length (PBL)' |
---|
289 | write(*,*) '(26) emissiv : ground emissivity' |
---|
290 | write(*,*) '(24 et 25) emisice : CO2 ice max emissivity ' |
---|
291 | write(*,*) '(22 et 23) albedice : CO2 ice cap albedos' |
---|
292 | write(*,*) '(31 et 32) iceradius : mean scat radius of CO2 snow' |
---|
293 | write(*,*) '(33 et 34) dtemisice : time scale for snow', |
---|
294 | & 'metamorphism' |
---|
295 | write(*,*) '(35) volcapa : soil volumetric heat capacity' |
---|
296 | write(*,*) '(18) obliquit : planet obliquity (deg)' |
---|
297 | write(*,*) '(17) peri_day : periastron date (sols since Ls=0)' |
---|
298 | write(*,*) '(15) periastr : min. star-planet dist (UA)' |
---|
299 | write(*,*) '(16) apoastr : max. star-planet (UA)' |
---|
300 | write(*,*) '(14) year_day : length of year (in sols)' |
---|
301 | write(*,*) '(5) rad : radius of the planet (m)' |
---|
302 | write(*,*) '(6) omeg : planet rotation rate (rad/s)' |
---|
303 | write(*,*) '(7) g : gravity (m/s2)' |
---|
304 | write(*,*) '(8) mugaz : molecular mass ' |
---|
305 | write(*,*) ' of the atmosphere (g/mol)' |
---|
306 | write(*,*) '(9) rcp : r/Cp' |
---|
307 | write(*,*) '(8)+(9) calc_cpp_mugaz : r/Cp and mugaz ' |
---|
308 | write(*,*) ' computed from gases.def' |
---|
309 | write(*,*) '(10) daysec : length of a sol (s)' |
---|
310 | write(*,*) |
---|
311 | |
---|
312 | |
---|
313 | do while(modif(1:1).ne.'hello') |
---|
314 | write(*,*) |
---|
315 | write(*,*) |
---|
316 | write(*,*) 'Changes to perform ?' |
---|
317 | write(*,*) ' (enter keyword or return )' |
---|
318 | write(*,*) |
---|
319 | read(*,fmt='(a20)') modif |
---|
320 | if (modif(1:1) .eq. ' ') goto 999 |
---|
321 | |
---|
322 | write(*,*) |
---|
323 | write(*,*) modif(1:len_trim(modif)) , ' : ' |
---|
324 | |
---|
325 | if (modif(1:len_trim(modif)) .eq. 'day_ini') then |
---|
326 | write(*,*) 'current value:',day_ini |
---|
327 | write(*,*) 'enter new value:' |
---|
328 | 101 read(*,*,iostat=ierr) day_ini |
---|
329 | if(ierr.ne.0) goto 101 |
---|
330 | write(*,*) ' ' |
---|
331 | write(*,*) 'day_ini (new value):',day_ini |
---|
332 | |
---|
333 | else if (modif(1:len_trim(modif)) .eq. 'z0') then |
---|
334 | write(*,*) 'current value:',z0 |
---|
335 | write(*,*) 'enter new value:' |
---|
336 | 102 read(*,*,iostat=ierr) z0 |
---|
337 | if(ierr.ne.0) goto 102 |
---|
338 | write(*,*) ' ' |
---|
339 | write(*,*) ' z0 (new value):',z0 |
---|
340 | |
---|
341 | else if (modif(1:len_trim(modif)) .eq. 'emin_turb') then |
---|
342 | write(*,*) 'current value:',emin_turb |
---|
343 | write(*,*) 'enter new value:' |
---|
344 | 103 read(*,*,iostat=ierr) emin_turb |
---|
345 | if(ierr.ne.0) goto 103 |
---|
346 | write(*,*) ' ' |
---|
347 | write(*,*) ' emin_turb (new value):',emin_turb |
---|
348 | |
---|
349 | else if (modif(1:len_trim(modif)) .eq. 'lmixmin') then |
---|
350 | write(*,*) 'current value:',lmixmin |
---|
351 | write(*,*) 'enter new value:' |
---|
352 | 104 read(*,*,iostat=ierr) lmixmin |
---|
353 | if(ierr.ne.0) goto 104 |
---|
354 | write(*,*) ' ' |
---|
355 | write(*,*) ' lmixmin (new value):',lmixmin |
---|
356 | |
---|
357 | else if (modif(1:len_trim(modif)) .eq. 'emissiv') then |
---|
358 | write(*,*) 'current value:',emissiv |
---|
359 | write(*,*) 'enter new value:' |
---|
360 | 105 read(*,*,iostat=ierr) emissiv |
---|
361 | if(ierr.ne.0) goto 105 |
---|
362 | write(*,*) ' ' |
---|
363 | write(*,*) ' emissiv (new value):',emissiv |
---|
364 | |
---|
365 | else if (modif(1:len_trim(modif)) .eq. 'emisice') then |
---|
366 | write(*,*) 'current value emisice(1) North:',emisice(1) |
---|
367 | write(*,*) 'enter new value:' |
---|
368 | 106 read(*,*,iostat=ierr) emisice(1) |
---|
369 | if(ierr.ne.0) goto 106 |
---|
370 | write(*,*) |
---|
371 | write(*,*) ' emisice(1) (new value):',emisice(1) |
---|
372 | write(*,*) |
---|
373 | |
---|
374 | write(*,*) 'current value emisice(2) South:',emisice(2) |
---|
375 | write(*,*) 'enter new value:' |
---|
376 | 107 read(*,*,iostat=ierr) emisice(2) |
---|
377 | if(ierr.ne.0) goto 107 |
---|
378 | write(*,*) |
---|
379 | write(*,*) ' emisice(2) (new value):',emisice(2) |
---|
380 | |
---|
381 | else if (modif(1:len_trim(modif)) .eq. 'albedice') then |
---|
382 | write(*,*) 'current value albedice(1) North:',albedice(1) |
---|
383 | write(*,*) 'enter new value:' |
---|
384 | 108 read(*,*,iostat=ierr) albedice(1) |
---|
385 | if(ierr.ne.0) goto 108 |
---|
386 | write(*,*) |
---|
387 | write(*,*) ' albedice(1) (new value):',albedice(1) |
---|
388 | write(*,*) |
---|
389 | |
---|
390 | write(*,*) 'current value albedice(2) South:',albedice(2) |
---|
391 | write(*,*) 'enter new value:' |
---|
392 | 109 read(*,*,iostat=ierr) albedice(2) |
---|
393 | if(ierr.ne.0) goto 109 |
---|
394 | write(*,*) |
---|
395 | write(*,*) ' albedice(2) (new value):',albedice(2) |
---|
396 | |
---|
397 | else if (modif(1:len_trim(modif)) .eq. 'iceradius') then |
---|
398 | write(*,*) 'current value iceradius(1) North:',iceradius(1) |
---|
399 | write(*,*) 'enter new value:' |
---|
400 | 110 read(*,*,iostat=ierr) iceradius(1) |
---|
401 | if(ierr.ne.0) goto 110 |
---|
402 | write(*,*) |
---|
403 | write(*,*) ' iceradius(1) (new value):',iceradius(1) |
---|
404 | write(*,*) |
---|
405 | |
---|
406 | write(*,*) 'current value iceradius(2) South:',iceradius(2) |
---|
407 | write(*,*) 'enter new value:' |
---|
408 | 111 read(*,*,iostat=ierr) iceradius(2) |
---|
409 | if(ierr.ne.0) goto 111 |
---|
410 | write(*,*) |
---|
411 | write(*,*) ' iceradius(2) (new value):',iceradius(2) |
---|
412 | |
---|
413 | else if (modif(1:len_trim(modif)) .eq. 'dtemisice') then |
---|
414 | write(*,*) 'current value dtemisice(1) North:',dtemisice(1) |
---|
415 | write(*,*) 'enter new value:' |
---|
416 | 112 read(*,*,iostat=ierr) dtemisice(1) |
---|
417 | if(ierr.ne.0) goto 112 |
---|
418 | write(*,*) |
---|
419 | write(*,*) ' dtemisice(1) (new value):',dtemisice(1) |
---|
420 | write(*,*) |
---|
421 | |
---|
422 | write(*,*) 'current value dtemisice(2) South:',dtemisice(2) |
---|
423 | write(*,*) 'enter new value:' |
---|
424 | 113 read(*,*,iostat=ierr) dtemisice(2) |
---|
425 | if(ierr.ne.0) goto 113 |
---|
426 | write(*,*) |
---|
427 | write(*,*) ' dtemisice(2) (new value):',dtemisice(2) |
---|
428 | |
---|
429 | else if (modif(1:len_trim(modif)) .eq. 'obliquit') then |
---|
430 | write(*,*) 'current value:',obliquit |
---|
431 | write(*,*) 'obliquit should be 25.19 on current Mars' |
---|
432 | write(*,*) 'enter new value:' |
---|
433 | 115 read(*,*,iostat=ierr) obliquit |
---|
434 | if(ierr.ne.0) goto 115 |
---|
435 | write(*,*) |
---|
436 | write(*,*) ' obliquit (new value):',obliquit |
---|
437 | |
---|
438 | else if (modif(1:len_trim(modif)) .eq. 'peri_day') then |
---|
439 | write(*,*) 'current value:',peri_day |
---|
440 | write(*,*) 'peri_day should be 485 on current Mars' |
---|
441 | write(*,*) 'enter new value:' |
---|
442 | 116 read(*,*,iostat=ierr) peri_day |
---|
443 | if(ierr.ne.0) goto 116 |
---|
444 | write(*,*) |
---|
445 | write(*,*) ' peri_day (new value):',peri_day |
---|
446 | |
---|
447 | else if (modif(1:len_trim(modif)) .eq. 'periastr') then |
---|
448 | write(*,*) 'current value:',periastr |
---|
449 | write(*,*) 'periastr should be 206.66 on present-day Mars' |
---|
450 | write(*,*) 'enter new value:' |
---|
451 | 117 read(*,*,iostat=ierr) periastr |
---|
452 | if(ierr.ne.0) goto 117 |
---|
453 | write(*,*) |
---|
454 | write(*,*) ' periastr (new value):',periastr |
---|
455 | |
---|
456 | else if (modif(1:len_trim(modif)) .eq. 'apoastr') then |
---|
457 | write(*,*) 'current value:',apoastr |
---|
458 | write(*,*) 'apoastr should be 249.22 on present-day Mars' |
---|
459 | write(*,*) 'enter new value:' |
---|
460 | 118 read(*,*,iostat=ierr) apoastr |
---|
461 | if(ierr.ne.0) goto 118 |
---|
462 | write(*,*) |
---|
463 | write(*,*) ' apoastr (new value):',apoastr |
---|
464 | |
---|
465 | else if (modif(1:len_trim(modif)) .eq. 'volcapa') then |
---|
466 | write(*,*) 'current value:',volcapa |
---|
467 | write(*,*) 'enter new value:' |
---|
468 | 119 read(*,*,iostat=ierr) volcapa |
---|
469 | if(ierr.ne.0) goto 119 |
---|
470 | write(*,*) |
---|
471 | write(*,*) ' volcapa (new value):',volcapa |
---|
472 | |
---|
473 | else if (modif(1:len_trim(modif)).eq.'rad') then |
---|
474 | write(*,*) 'current value:',rad |
---|
475 | write(*,*) 'enter new value:' |
---|
476 | 120 read(*,*,iostat=ierr) rad |
---|
477 | if(ierr.ne.0) goto 120 |
---|
478 | write(*,*) |
---|
479 | write(*,*) ' rad (new value):',rad |
---|
480 | |
---|
481 | else if (modif(1:len_trim(modif)).eq.'omeg') then |
---|
482 | write(*,*) 'current value:',omeg |
---|
483 | write(*,*) 'enter new value:' |
---|
484 | 121 read(*,*,iostat=ierr) omeg |
---|
485 | if(ierr.ne.0) goto 121 |
---|
486 | write(*,*) |
---|
487 | write(*,*) ' omeg (new value):',omeg |
---|
488 | |
---|
489 | else if (modif(1:len_trim(modif)).eq.'g') then |
---|
490 | write(*,*) 'current value:',g |
---|
491 | write(*,*) 'enter new value:' |
---|
492 | 122 read(*,*,iostat=ierr) g |
---|
493 | if(ierr.ne.0) goto 122 |
---|
494 | write(*,*) |
---|
495 | write(*,*) ' g (new value):',g |
---|
496 | |
---|
497 | else if (modif(1:len_trim(modif)).eq.'mugaz') then |
---|
498 | write(*,*) 'current value:',mugaz |
---|
499 | write(*,*) 'enter new value:' |
---|
500 | 123 read(*,*,iostat=ierr) mugaz |
---|
501 | if(ierr.ne.0) goto 123 |
---|
502 | write(*,*) |
---|
503 | write(*,*) ' mugaz (new value):',mugaz |
---|
504 | r=8.314511/(mugaz/1000.0) |
---|
505 | write(*,*) ' R (new value):',r |
---|
506 | |
---|
507 | else if (modif(1:len_trim(modif)).eq.'rcp') then |
---|
508 | write(*,*) 'current value:',rcp |
---|
509 | write(*,*) 'enter new value:' |
---|
510 | 124 read(*,*,iostat=ierr) rcp |
---|
511 | if(ierr.ne.0) goto 124 |
---|
512 | write(*,*) |
---|
513 | write(*,*) ' rcp (new value):',rcp |
---|
514 | r=8.314511/(mugaz/1000.0) |
---|
515 | cpp=r/rcp |
---|
516 | write(*,*) ' cpp (new value):',cpp |
---|
517 | |
---|
518 | else if (modif(1:len_trim(modif)).eq.'calc_cpp_mugaz') then |
---|
519 | write(*,*) 'current value rcp, mugaz:',rcp,mugaz |
---|
520 | check_cpp_match=.false. |
---|
521 | force_cpp=.false. |
---|
522 | call su_gases |
---|
523 | call calc_cpp_mugaz |
---|
524 | write(*,*) |
---|
525 | write(*,*) ' cpp (new value):',cpp |
---|
526 | write(*,*) ' mugaz (new value):',mugaz |
---|
527 | r=8.314511/(mugaz/1000.0) |
---|
528 | rcp=r/cpp |
---|
529 | write(*,*) ' rcp (new value):',rcp |
---|
530 | |
---|
531 | else if (modif(1:len_trim(modif)).eq.'daysec') then |
---|
532 | write(*,*) 'current value:',daysec |
---|
533 | write(*,*) 'enter new value:' |
---|
534 | 125 read(*,*,iostat=ierr) daysec |
---|
535 | if(ierr.ne.0) goto 125 |
---|
536 | write(*,*) |
---|
537 | write(*,*) ' daysec (new value):',daysec |
---|
538 | |
---|
539 | ! added by RW! |
---|
540 | else if (modif(1:len_trim(modif)).eq.'year_day') then |
---|
541 | write(*,*) 'current value:',year_day |
---|
542 | write(*,*) 'enter new value:' |
---|
543 | 126 read(*,*,iostat=ierr) year_day |
---|
544 | if(ierr.ne.0) goto 126 |
---|
545 | write(*,*) |
---|
546 | write(*,*) ' year_day (new value):',year_day |
---|
547 | |
---|
548 | endif |
---|
549 | enddo ! of do while(modif(1:1).ne.'hello') |
---|
550 | |
---|
551 | 999 continue |
---|
552 | |
---|
553 | c----------------------------------------------------------------------- |
---|
554 | c Write values of physical constants after modifications |
---|
555 | c----------------------------------------------------------------------- |
---|
556 | |
---|
557 | write(*,*) '*****************************************************' |
---|
558 | write(*,*) 'Reading tab_cntrl when calling tabfi AFTER changes' |
---|
559 | write(*,*) '*****************************************************' |
---|
560 | write(*,5) '(1) = ngrid?',tab_cntrl(tab0+1),float(ngrid) |
---|
561 | write(*,5) '(2) lmax',tab_cntrl(tab0+2),float(lmax) |
---|
562 | write(*,5) '(3) day_ini',tab_cntrl(tab0+3),float(day_ini) |
---|
563 | write(*,5) '(5) rad',tab_cntrl(tab0+5),rad |
---|
564 | write(*,5) '(10) daysec',tab_cntrl(tab0+10),daysec |
---|
565 | write(*,6) '(6) omeg',tab_cntrl(tab0+6),omeg |
---|
566 | write(*,5) '(7) g',tab_cntrl(tab0+7),g |
---|
567 | write(*,5) '(8) mugaz',tab_cntrl(tab0+8),mugaz |
---|
568 | write(*,5) '(9) rcp',tab_cntrl(tab0+9),rcp |
---|
569 | write(*,6) '(11) dtphys?',tab_cntrl(tab0+11),dtphys |
---|
570 | |
---|
571 | write(*,5) '(14) year_day',tab_cntrl(tab0+14),year_day |
---|
572 | write(*,5) '(15) periastr',tab_cntrl(tab0+15),periastr |
---|
573 | write(*,5) '(16) apoastr',tab_cntrl(tab0+16),apoastr |
---|
574 | write(*,5) '(17) peri_day',tab_cntrl(tab0+17),peri_day |
---|
575 | write(*,5) '(18) obliquit',tab_cntrl(tab0+18),obliquit |
---|
576 | |
---|
577 | write(*,6) '(19) z0',tab_cntrl(tab0+19),z0 |
---|
578 | write(*,6) '(21) emin_turb',tab_cntrl(tab0+21),emin_turb |
---|
579 | write(*,5) '(20) lmixmin',tab_cntrl(tab0+20),lmixmin |
---|
580 | |
---|
581 | write(*,5) '(26) emissiv',tab_cntrl(tab0+26),emissiv |
---|
582 | write(*,5) '(24) emisice(1)',tab_cntrl(tab0+24),emisice(1) |
---|
583 | write(*,5) '(25) emisice(2)',tab_cntrl(tab0+25),emisice(2) |
---|
584 | write(*,5) '(22) albedice(1)',tab_cntrl(tab0+22),albedice(1) |
---|
585 | write(*,5) '(23) albedice(2)',tab_cntrl(tab0+23),albedice(2) |
---|
586 | write(*,6) '(31) iceradius(1)',tab_cntrl(tab0+31),iceradius(1) |
---|
587 | write(*,6) '(32) iceradius(2)',tab_cntrl(tab0+32),iceradius(2) |
---|
588 | write(*,5) '(33) dtemisice(1)',tab_cntrl(tab0+33),dtemisice(1) |
---|
589 | write(*,5) '(34) dtemisice(2)',tab_cntrl(tab0+34),dtemisice(2) |
---|
590 | |
---|
591 | write(*,5) '(35) volcapa',tab_cntrl(tab0+35),volcapa |
---|
592 | |
---|
593 | write(*,*) |
---|
594 | write(*,*) |
---|
595 | |
---|
596 | ENDIF ! of if (Lmodif == 1) |
---|
597 | |
---|
598 | c----------------------------------------------------------------------- |
---|
599 | c Save some constants for later use (as routine arguments) |
---|
600 | c----------------------------------------------------------------------- |
---|
601 | p_omeg = omeg |
---|
602 | p_g = g |
---|
603 | p_cpp = cpp |
---|
604 | p_mugaz = mugaz |
---|
605 | p_daysec = daysec |
---|
606 | p_rad=rad |
---|
607 | |
---|
608 | |
---|
609 | end |
---|