Context Navigation

inifis.F @ 3436

Last change on this file since 3436 was 3373, checked in by afalco, 6 months ago
Pluto.old PCM: Corrected some runtime issues with gfortran. AF
File size: 32.1 KB

Line
1	SUBROUTINE inifis(ngrid,nlayer,
2	$ day_ini,pdaysec,ptimestep,
3	$ plat,plon,parea,
4	$ prad,pg,pr,pcpp)
5	!
6	!=======================================================================
7	!
8	! purpose:
9	! -------
10	!
11	! Initialisation for the physical parametrisations of the LMD
12	! Pluto atmospheric general circulation modele.
13	! used by 1D Model only
14	!
15	! author: Frederic Hourdin 15 / 10 /93
16	! -------
17	! modified: Sebastien Lebonnois 11/06/2003 (new callphys.def)
18	! Ehouarn Millour (oct. 2008) tracers are now identified
19	! by their names and may not be contiguously
20	! stored in the q(:,:,:,:) array
21	! E.M. (june 2009) use getin routine to load parameters
22	!
23	!
24	! arguments:
25	! ----------
26	!
27	! input:
28	! ------
29	!
30	! ngrid Size of the horizontal grid.
31	! All internal loops are performed on that grid.
32	! nlayer Number of vertical layers.
33	! pdayref Day of reference for the simulation
34	! pday Number of days counted from the North. Spring
35	! equinoxe.
36	!
37	!=======================================================================
38	!
39	!-----------------------------------------------------------------------
40	! declarations:
41	! -------------
42	! to use 'getin'
43	USE ioipsl_getincom
44	use planet_h
45	use comgeomfi_h
46	use datafile_mod
47	IMPLICIT NONE
48	#include "dimensions.h"
49	#include "dimphys.h"
50	!#include "planete.h"
51	#include "comcstfi.h"
52	#include "comsaison.h"
53	#include "comdiurn.h"
54	#include "callkeys.h"
55	#include "surfdat.h"
56
57
58	REAL prad,pg,pr,pcpp,pdaysec,ptimestep
59
60	INTEGER ngrid,nlayer
61	REAL plat(ngrid),plon(ngrid),parea(ngridmx)
62	integer day_ini
63	INTEGER ig,ierr
64
65	EXTERNAL iniorbit,orbite
66	EXTERNAL SSUM
67	REAL SSUM
68
69	CHARACTER ch1*12
70	CHARACTER ch80*80
71
72	logical chem, h2o
73	logical :: parameter, doubleq=.false.
74
75	rad=prad
76	daysec=pdaysec
77	dtphys=ptimestep
78	cpp=pcpp
79	g=pg
80	r=pr
81	rcp=r/cpp
82
83
84	avocado = 6.02214179e23 ! added by RW
85
86
87	! --------------------------------------------------------
88	! The usual Tests
89	! --------------
90	IF (nlayer.NE.nlayermx) THEN
91	PRINT*,'STOP in inifis'
92	PRINT*,'Probleme de dimensions :'
93	PRINT*,'nlayer = ',nlayer
94	PRINT*,'nlayermx = ',nlayermx
95	STOP
96	ENDIF
97
98	IF (ngrid.NE.ngridmx) THEN
99	PRINT*,'STOP in inifis'
100	PRINT*,'Probleme de dimensions :'
101	PRINT*,'ngrid = ',ngrid
102	PRINT*,'ngridmx = ',ngridmx
103	STOP
104	ENDIF
105
106	! --------------------------------------------------------------
107	! Reading the "callphys.def" file controlling some key options
108	! --------------------------------------------------------------
109
110	! check that 'callphys.def' file is around
111	OPEN(99,file='callphys.def',status='old',form='formatted'
112	& ,iostat=ierr)
113	CLOSE(99)
114
115	IF(ierr.EQ.0) THEN
116	PRINT*
117	PRINT*
118	PRINT*,'--------------------------------------------'
119	PRINT*,' inifis: Parametres pour la physique (callphys.def)'
120	PRINT*,'--------------------------------------------'
121
122	!***************************************************************
123	!! Directories / Files
124
125	write(,) "Directory where external input files are:"
126	datadir="../../../datagcm" ! default path
127	call getin("datadir",datadir)
128	write(,) " datafile = ",trim(datadir)
129
130	write(,) "Haze optical properties datafile"
131	hazeprop="../../../datagcm" ! default path
132	call getin("hazeprop",hazeprop)
133	write(,) " hazeprop = ",trim(hazeprop)
134	write(,) "Haze radii datafile"
135	hazerad_file="hazerad.txt" ! default file
136	call getin("hazerad_file",hazerad_file)
137	write(,) " hazerad_file = ",trim(hazerad_file)
138	write(,) "Haze mmr datafile"
139	hazemmr_file="None" ! default file
140	call getin("hazemmr_file",hazemmr_file)
141	write(,) " hazemmr_file = ",trim(hazemmr_file)
142	write(,) "Haze dens datafile"
143	hazedens_file="None" ! default file
144	call getin("hazedens_file",hazedens_file)
145	write(,) " hazedens_file = ",trim(hazedens_file)
146
147	!***************************************************************
148	!! Basics
149
150	write(,) "Run with or without tracer transport ?"
151	tracer=.false. ! default value
152	call getin("tracer",tracer)
153	write(,) " tracer = ",tracer
154
155	write(,) "Diurnal cycle ?"
156	write(,) "(if diurnal=false, diurnal averaged solar heating)"
157	diurnal=.true. ! default value
158	call getin("diurnal",diurnal)
159	write(,) " diurnal = ",diurnal
160
161	write(,) "Seasonal cycle ?"
162	write(,) "(if season=false, Ls stays constant, to value ",
163	& "set in 'start'"
164	season=.true. ! default value
165	call getin("season",season)
166	write(,) " season = ",season
167
168	write(,) "Run 1D or 3D version (false) / fast version
169	& volatile transport (true) ?"
170	fast=.false. ! default value
171	call getin("fast",fast)
172	write(,) " fast = ",fast
173
174	write(,) "1D options:"
175	globmean1d=.false. ! default value
176	call getin("globmean1d",globmean1d)
177	write(,) " globmean1d = ",globmean1d
178	kmixmin=1 ! default value
179	call getin("kmixmin",kmixmin)
180	write(,) " kmixmin = ",kmixmin
181	kmix_proffix=.false. ! default value
182	call getin("kmix_proffix",kmix_proffix)
183	write(,) " kmix_proffix = ",kmix_proffix
184
185	write(,) "Write some extra output to the screen ?"
186	lwrite=.false. ! default value
187	call getin("lwrite",lwrite)
188	write(,) " lwrite = ",lwrite
189
190	write(,) "Save statistics in file stats.nc ?"
191	callstats=.true. ! default value
192	call getin("callstats",callstats)
193	write(,) " callstats = ",callstats
194
195	write(,) "Test energy conservation of model physics ?"
196	enertest=.false. ! default value
197	call getin("enertest",enertest)
198	write(,) " enertest = ",enertest
199	conservn2=.false. ! default value
200	call getin("conservn2",conservn2)
201	write(,) " conservn2 = ",conservn2
202
203	write(,) "Run with Triton orbit ?"
204	triton=.false. ! default value
205	call getin("triton",triton)
206	write(,) " triton = ",triton
207	convergeps=.false. ! default value
208	call getin("convergeps",convergeps)
209	write(,) " convergeps = ",convergeps
210
211	write(,) "KBO runs (eris, makemake) ?"
212	kbo=.false. ! default value
213	call getin("kbo",kbo)
214	write(,) " kbo = ",kbo
215
216	write(,) "Specific paleo run ?"
217	paleo=.false. ! default value
218	call getin("paleo",paleo)
219	write(,) " paleo = ",paleo
220
221	write(,) "paleoclimate step (Earth years) "
222	paleoyears=10000. ! default value
223	call getin("paleoyears",paleoyears)
224	write(,) " paleoyears = ",paleoyears
225
226	!***************************************************************
227	!! Radiative transfer options
228
229	write(,) "call radiative transfer ?"
230	callrad=.true. ! default value
231	call getin("callrad",callrad)
232	write(,) " callrad = ",callrad
233
234	write(,)"rad.transfer is computed every iradia",
235	& " physical timestep"
236	iradia=1 ! default value
237	call getin("iradia",iradia)
238	write(,)" iradia = ",iradia
239
240	write(,) "call correlated-k radiative transfer ?"
241	corrk=.true. ! default value
242	call getin("corrk",corrk)
243	write(,) " corrk = ",corrk
244
245	write(,) "call specific cooling for radiative transfer ?"
246	cooling=.false. ! default value
247	call getin("cooling",cooling)
248	write(,) " cooling = ",cooling
249
250	write(,) "fixed gaseous CH4 mixing ratio for rad transfer ?"
251	ch4fix=.false. ! default value
252	call getin("ch4fix",ch4fix)
253	write(,) " ch4fix = ",ch4fix
254	write(,) "fixed gaseous CH4 mixing ratio for rad transfer ?"
255	vmrch4fix=0.5 ! default value
256	call getin("vmrch4fix",vmrch4fix)
257	write(,) " vmrch4fix = ",vmrch4fix
258	vmrch4_proffix=.false. ! default value
259	call getin("vmrch4_proffix",vmrch4_proffix)
260	write(,) " vmrch4_proffix = ",vmrch4_proffix
261
262	write(,) "NLTE correction for methane heating rates?"
263	nlte=.false. ! default value
264	call getin("nlte",nlte)
265	write(,) " nlte = ",nlte
266	strobel=.false. ! default value
267	call getin("strobel",strobel)
268	write(,) " strobel = ",strobel
269
270	write(,) "call gaseous absorption in the visible bands?"
271	callgasvis=.false. ! default value
272	call getin("callgasvis",callgasvis)
273	write(,) " callgasvis = ",callgasvis
274
275	write(,)"Output mean OLR in 1D?"
276	meanOLR=.false.
277	call getin("meanOLR",meanOLR)
278	write(,)" meanOLR = ",meanOLR
279
280	write(,)"Output spectral OLR in 3D?"
281	specOLR=.false.
282	call getin("specOLR",specOLR)
283	write(,)" specOLR = ",specOLR
284
285	!***************************************************************
286	!! From generic model - to be cleaned
287
288	write(,)"Which star?"
289	startype=1 ! default value = Sol
290	call getin("startype",startype)
291	write(,)" startype = ",startype
292	write(,)"Value of stellar flux at 1 AU?"
293	Fat1AU=1366.0 ! default value = Sol today
294	call getin("Fat1AU",Fat1AU)
295	write(,)" Fat1AU = ",Fat1AU
296	write(,)"Default planetary temperature?"
297	tplanet=100.0
298	call getin("tplanet",tplanet)
299	write(,)" tplanet = ",tplanet
300
301	!***************************************************************
302	!! TRACERS options
303
304	write(,) "call N2 condensation ?"
305	N2cond=.true. ! default value
306	call getin("N2cond",N2cond)
307	write(,) " N2cond = ",N2cond
308
309	write(,) "call N2 cloud condensation ?"
310	condensn2=.false. ! default value
311	call getin("condensn2",condensn2)
312	write(,) "condensn2 = ",condensn2
313
314	write(,) "call no N2 frost formation?"
315	no_n2frost=.false. ! default value
316	call getin("no_n2frost",no_n2frost)
317	write(,) "no_n2frost = ",no_n2frost
318
319	write(,) "N2 condensation subtimestep?"
320	nbsub=20 ! default value
321	call getin("nbsub",nbsub)
322	write(,) " nbsub = ",nbsub
323
324	write(,) "Gravitationnal sedimentation ?"
325	sedimentation=.true. ! default value
326	call getin("sedimentation",sedimentation)
327	write(,) " sedimentation = ",sedimentation
328
329	write(,) "Compute glacial flow ?"
330	glaflow=.false. ! default value
331	call getin("glaflow",glaflow)
332	write(,) " glaflow = ",glaflow
333
334	write(,) "Compute methane cycle ?"
335	methane=.false. ! default value
336	call getin("methane",methane)
337	write(,) " methane = ",methane
338	condmetsurf=.true. ! default value
339	call getin("condmetsurf",condmetsurf)
340	write(,) " condmetsurf = ",condmetsurf
341
342	write(,) "Compute methane clouds ?"
343	metcloud=.false. ! default value
344	call getin("metcloud",metcloud)
345	write(,) " metcloud = ",metcloud
346
347	write(,) "Compute CO cycle ?"
348	carbox=.false. ! default value
349	call getin("carbox",carbox)
350	write(,) " carbox = ",carbox
351	condcosurf=.true. ! default value
352	call getin("condcosurf",condcosurf)
353	write(,) " condcosurf = ",condcosurf
354
355	write(,) "Compute CO clouds ?"
356	monoxcloud=.false. ! default value
357	call getin("monoxcloud",monoxcloud)
358	write(,) " monoxcloud = ",monoxcloud
359
360	write(,)"atmospheric redistribution (s):"
361	tau_n2=1. ! default value
362	call getin("tau_n2",tau_n2)
363	write(,)" tau_n2 = ",tau_n2
364	tau_ch4=1.E7 ! default value
365	call getin("tau_ch4",tau_ch4)
366	write(,)" tau_ch4 = ",tau_ch4
367	tau_co=1. ! default value
368	call getin("tau_co",tau_co)
369	write(,)" tau_co = ",tau_co
370	tau_prechaze=1. ! default value
371	call getin("tau_prechaze",tau_prechaze)
372	write(,)" tau_prechaze = ",tau_prechaze
373
374	write(,) "Day fraction for limited cold trap in SP?"
375	dayfrac=0. ! default value
376	call getin("dayfrac",dayfrac)
377	write(,) " dayfrac = ",dayfrac
378	thresh_non2=0. ! default value
379	call getin("thresh_non2",thresh_non2)
380	write(,) " thresh_non2 = ",thresh_non2
381
382	!***************************************************************
383	!! Haze options
384
385	write(,)"Production of haze ?"
386	haze=.false. ! default value
387	call getin("haze",haze)
388	write(,)" haze = ",haze
389	hazeconservch4=.false. ! conservch4, by default value ch4 is photolyzed
390	call getin("hazeconservch4",hazeconservch4)
391	write(,)"hazconservch4 = ",hazeconservch4
392	write(,)"Production of haze (fast version) ?"
393	fasthaze=.false. ! default value
394	call getin("fasthaze",fasthaze)
395	write(,)"fasthaze = ",fasthaze
396
397	write(,)"Add source of haze ?"
398	source_haze=.false. ! default value
399	call getin("source_haze",source_haze)
400	write(,)" source_haze = ",source_haze
401	mode_hs=0 ! mode haze source default value
402	call getin("mode_hs",mode_hs)
403	write(,)" mode_hs",mode_hs
404	kfix=1 ! default value
405	call getin("kfix",kfix)
406	write(,)"levels kfix",kfix
407	fracsource=0.1 ! default value
408	call getin("fracsource",fracsource)
409	write(,)" fracsource",fracsource
410	latsource=30. ! default value
411	call getin("latsource",latsource)
412	write(,)" latsource",latsource
413	lonsource=180. ! default value
414	call getin("lonsource",lonsource)
415	write(,)" lonsource",lonsource
416
417	write(,)"Radiatively active haze ?"
418	aerohaze=.false. ! default value
419	call getin("aerohaze",aerohaze)
420	write(,)"aerohaze = ",aerohaze
421
422	write(,)"Haze monomer radius ?"
423	rad_haze=10.e-9 ! default value
424	call getin("rad_haze",rad_haze)
425	write(,)"rad_haze = ",rad_haze
426
427	write(,)"fractal particle ?"
428	fractal=.false. ! default value
429	call getin("fractal",fractal)
430	write(,)"fractal = ",fractal
431	nb_monomer=10 ! default value
432	call getin("nb_monomer",nb_monomer)
433	write(,)"nb_monomer = ",nb_monomer
434
435	write(,)"fixed radius profile from txt file ?"
436	haze_radproffix=.false. ! default value
437	call getin("haze_radproffix",haze_radproffix)
438	write(,)"haze_radproffix = ",haze_radproffix
439	write(,)"fixed MMR profile from txt file ?"
440	haze_proffix=.false. ! default value
441	call getin("haze_proffix",haze_proffix)
442	write(,)"haze_proffix = ",haze_proffix
443
444	write(,)"Number mix ratio of haze particles for co clouds:"
445	Nmix_co=100000. ! default value
446	call getin("Nmix_co",Nmix_co)
447	write(,)" Nmix_co = ",Nmix_co
448
449	write(,)"Number mix ratio of haze particles for ch4 clouds:"
450	Nmix_ch4=100000. ! default value
451	call getin("Nmix_ch4",Nmix_ch4)
452	write(,)" Nmix_ch4 = ",Nmix_ch4
453
454	!***************************************************************
455	!! Surface properties
456
457	!********* N2 *******************************
458
459	write(,) "Mode for changing N2 albedo"
460	mode_n2=0 ! default value
461	call getin("mode_n2",mode_n2)
462	write(,) " mode_n2 = ",mode_n2
463	thres_n2ice=1. ! default value
464	call getin("thres_n2ice",thres_n2ice)
465	write(,) " thres_n2ice = ",thres_n2ice
466
467	write(,) "Diff of N2 albedo with thickness"
468	deltab=0. ! default value
469	call getin("deltab",deltab)
470	write(,) " deltab = ",deltab
471
472	write(,) "albedo N2 beta "
473	alb_n2b=0.7 ! default value
474	call getin("alb_n2b",alb_n2b)
475	write(,) " alb_n2b = ",alb_n2b
476
477	write(,) "albedo N2 alpha "
478	alb_n2a=0.7 ! default value
479	call getin("alb_n2a",alb_n2a)
480	write(,) " alb_n2a = ",alb_n2a
481
482	write(,) "emis N2 beta "
483	emis_n2b=0.7 ! default value
484	call getin("emis_n2b",emis_n2b)
485	write(,) " emis_n2b = ",emis_n2b
486
487	write(,) "emis N2 alpha "
488	emis_n2a=0.7 ! default value
489	call getin("emis_n2a",emis_n2a)
490	write(,) " emis_n2a = ",emis_n2a
491
492	!********* CH4 *******************************
493
494	write(,) "Mode for changing CH4 albedo"
495	mode_ch4=0 ! default value
496	call getin("mode_ch4",mode_ch4)
497	write(,) " mode_ch4 = ",mode_ch4
498	feedback_met=0 ! default value
499	call getin("feedback_met",feedback_met)
500	write(,) " feedback_met = ",feedback_met
501	thres_ch4ice=1. ! default value, mm
502	call getin("thres_ch4ice",thres_ch4ice)
503	write(,) " thres_ch4ice = ",thres_ch4ice
504	fdch4_latn=-1.
505	fdch4_lats=0.
506	fdch4_lone=0.
507	fdch4_lonw=-1.
508	fdch4_depalb=0.5
509	fdch4_finalb=0.95
510	fdch4_maxalb=0.99
511	fdch4_ampl=200.
512	fdch4_maxice=100.
513	call getin(" fdch4_latn",fdch4_latn)
514	call getin(" fdch4_lats",fdch4_lats)
515	call getin(" fdch4_lone",fdch4_lone)
516	call getin(" fdch4_lonw",fdch4_lonw)
517	call getin(" fdch4_depalb",fdch4_depalb)
518	call getin(" fdch4_finalb",fdch4_finalb)
519	call getin(" fdch4_maxalb",fdch4_maxalb)
520	call getin(" fdch4_maxice",fdch4_maxice)
521	call getin(" fdch4_ampl",fdch4_ampl)
522	write(,) " Values for albedo feedback = ",fdch4_latn,
523	& fdch4_lats,fdch4_lone,fdch4_lonw,fdch4_depalb,
524	& fdch4_finalb,fdch4_maxalb,fdch4_maxice,fdch4_ampl
525
526	write(,) "Latitude for diff albedo methane"
527	metlateq=25. ! default value
528	call getin("metlateq",metlateq)
529	write(,) " metlateq = ",metlateq
530
531	write(,) "Ls1 and Ls2 for change of ch4 albedo"
532	metls1=-1. ! default value
533	metls2=-2. ! default value
534	call getin("metls1",metls1)
535	call getin("metls2",metls2)
536
537	write(,) "albedo CH4 "
538	alb_ch4=0.5 ! default value
539	call getin("alb_ch4",alb_ch4)
540	write(,) " alb_ch4 = ",alb_ch4
541
542	write(,) "albedo equatorial CH4 "
543	alb_ch4_eq=alb_ch4 ! default value
544	call getin("alb_ch4_eq",alb_ch4_eq)
545	write(,) " alb_ch4_eq = ",alb_ch4_eq
546
547	write(,) "albedo south hemis CH4 "
548	alb_ch4_s=alb_ch4 ! default value
549	call getin("alb_ch4_s",alb_ch4_s)
550	write(,) " alb_ch4_s = ",alb_ch4_s
551
552	write(,) "emis CH4 "
553	emis_ch4=0.5 ! default value
554	call getin("emis_ch4",emis_ch4)
555	write(,) " emis_ch4 = ",emis_ch4
556
557	write(,) "CH4 lag for n2 sublimation limitation"
558	ch4lag=.false. ! default value
559	latlag=-90. ! default value
560	vmrlag=1. ! default value
561	call getin("ch4lag",ch4lag)
562	call getin("latlag",latlag)
563	call getin("vmrlag",vmrlag)
564	write(,) " ch4lag = ",ch4lag
565	write(,) " latlag = ",latlag
566	write(,) " vmrlag = ",vmrlag
567
568	write(,) "Max temperature for surface ?"
569	tsurfmax=.false. ! default value
570	albmin_ch4=0.3 ! default value
571	call getin("tsurfmax",tsurfmax)
572	call getin("albmin_ch4",albmin_ch4)
573	write(,) " tsurfmax = ",tsurfmax
574	write(,) " albmin_ch4 = ",albmin_ch4
575
576	!********* CO *******************************
577
578	write(,) "albedo CO "
579	alb_co=0.4 ! default value
580	call getin("alb_co",alb_co)
581	write(,) " alb_co = ",alb_co
582	thres_coice=1. ! default value, mm
583	call getin("thres_coice",thres_coice)
584	write(,) " thres_coice = ",thres_coice
585
586	write(,) "emis CO "
587	emis_co=0.4 ! default value
588	call getin("emis_co",emis_co)
589	write(,) " emis_co = ",emis_co
590
591	!********* THOLINS *******************************
592	write(,) "Mode for changing tholins albedo/emis"
593	mode_tholins=0 ! default value
594	call getin("mode_tholins",mode_tholins)
595	write(,) " mode_tholins = ",mode_tholins
596
597	write(,) "albedo tho "
598	alb_tho=0.1 ! default value
599	call getin("alb_tho",alb_tho)
600	write(,) " alb_tho = ",alb_tho
601
602	write(,) "albedo tho eq"
603	alb_tho_eq=0.1 ! default value
604	call getin("alb_tho_eq",alb_tho_eq)
605	write(,) " alb_tho_eq = ",alb_tho_eq
606
607	write(,) "emis tholins "
608	emis_tho=1. ! default value
609	call getin("emis_tho",emis_tho)
610	write(,) " emis_tho = ",emis_tho
611
612	write(,) "emis tholins eq"
613	emis_tho_eq=1. ! default value
614	call getin("emis_tho_eq",emis_tho_eq)
615	write(,) " emis_tho_eq = ",emis_tho_eq
616
617	write(,) "Latitude for diff albedo tholins"
618	tholateq=25. ! default value
619	call getin("tholateq",tholateq)
620	write(,) " tholateq = ",tholateq
621	tholatn=-1.
622	tholats=0.
623	tholone=0.
624	tholonw=-1.
625	alb_tho_spe=0.1 ! default value
626	emis_tho_spe=1. ! default value
627	call getin(" tholatn",tholatn)
628	call getin(" tholats",tholats)
629	call getin(" tholone",tholone)
630	call getin(" tholonw",tholonw)
631	write(,) " Values for special tholins albedo = ",tholatn,
632	& tholats,tholone,tholonw,alb_tho_spe,emis_tho_spe
633
634	write(,) "Specific albedo"
635	spelon1=-180. ! default value
636	spelon2=180. ! default value
637	spelat1=-90. ! default value
638	spelat2=90. ! default value
639	specalb=.false. ! default value
640	write(,) "albedo/emis spe "
641	albspe=0.1 ! default value
642	emispe=1. ! default value
643	call getin("spelon1",spelon1)
644	call getin("spelon2",spelon2)
645	call getin("spelat1",spelat1)
646	call getin("spelat2",spelat2)
647	call getin("specalb",specalb)
648	call getin("albspe",albspe)
649	call getin("emispe",emispe)
650
651	write(,) " specific = ",specalb
652
653	!******************** TI ***************************
654
655	write(,) "Change TI with time"
656	changeti=.false. ! default value
657	call getin("changeti",changeti)
658	write(,) " changeti = ",changeti
659	changetid=.false. ! default value for diurnal TI
660	call getin("changetid",changetid)
661	write(,) " changetid = ",changetid
662
663	write(,) "IT N2 "
664	ITN2=800. ! default value
665	call getin("ITN2",ITN2)
666	write(,) " ITN2 = ",ITN2
667	ITN2d=20. ! default value
668	call getin("ITN2d",ITN2d)
669	write(,) " ITN2d = ",ITN2d
670
671	write(,) "IT CH4"
672	ITCH4=800. ! default value
673	call getin("ITCH4",ITCH4)
674	write(,) " ITCH4 = ",ITCH4
675	ITCH4d=20. ! default value
676	call getin("ITCH4d",ITCH4d)
677	write(,) " ITCH4d = ",ITCH4d
678
679	write(,) "IT H2O"
680	ITH2O=800. ! default value
681	call getin("ITH2O",ITH2O)
682	write(,) " ITH2O = ",ITH2O
683	ITH2Od=20. ! default value
684	call getin("ITH2Od",ITH2Od)
685	write(,) " ITH2Od = ",ITH2Od
686
687	!***************************************************************
688	!! Other Physics Options
689
690	write(,)"Geothermal flux (W) at the bottom layer"
691	fluxgeo=0. ! default value
692	call getin("fluxgeo",fluxgeo)
693	write(,)" fluxgeo = ",fluxgeo
694
695	write(,)"Assymetry flux (W) at the bottom layer"
696	assymflux=.false. ! default value
697	call getin("assymflux",assymflux)
698	write(,)" assymflux = ",assymflux
699	write(,)"Geothermal flux (W) for assymetry"
700	fluxgeo2=fluxgeo ! default value
701	call getin("fluxgeo2",fluxgeo2)
702	write(,)" fluxgeo2 = ",fluxgeo2
703	write(,)"Warm patch of flux"
704	patchflux=0 ! default value
705	call getin("patchflux",patchflux)
706	write(,)" patchflux = ",patchflux
707
708	write(,) "call turbulent vertical diffusion ?"
709	calldifv=.false. ! default value
710	call getin("calldifv",calldifv)
711	write(,) " calldifv = ",calldifv
712	vertdiff=.true. ! default value
713	call getin("vertdiff",vertdiff)
714	write(,) " vertdiff = ",vertdiff
715
716	write(,) "call convective adjustment ?"
717	calladj=.false. ! default value
718	call getin("calladj",calladj)
719	write(,) " calladj = ",calladj
720
721	write(,)"call thermal conduction in the soil ?"
722	callsoil=.true. ! default value
723	call getin("callsoil",callsoil)
724	write(,) " callsoil = ",callsoil
725
726	write(,) "photochemistry: include chemical species"
727	photochem=.false. ! default value
728	call getin("photochem",photochem)
729	write(,) " photochem = ",photochem
730
731	write(,) "call thermal conduction ?"
732	callconduct=.false. ! default value
733	call getin("callconduct",callconduct)
734	write(,) " callconduct = ",callconduct
735
736	write(,) "call phitop ?"
737	phitop=0. ! default value
738	call getin("phitop",phitop)
739	write(,) " phitop = ",phitop
740
741	write(,) "call molecular viscosity ?"
742	callmolvis=.false. ! default value
743	call getin("callmolvis",callmolvis)
744	write(,) " callmolvis = ",callmolvis
745
746	write(,) "call molecular diffusion ?"
747	callmoldiff=.false. ! default value
748	call getin("callmoldiff",callmoldiff)
749	write(,) " callmoldiff = ",callmoldiff
750
751	c write(,) "call thermosphere ?"
752	c callthermos=.false. ! default value
753	c call getin("callthermos",callthermos)
754	c write(,) " callthermos = ",callthermos
755
756
757	c write(,) "call EUV heating ?",
758	c & " (only if callthermos=.true.)"
759	c calleuv=.false. ! default value
760	c call getin("calleuv",calleuv)
761	c write(,) " calleuv = ",calleuv
762
763	c write(,) "call thermospheric photochemistry ?",
764	c & " (only if callthermos=.true.)"
765	c thermochem=.false. ! default value
766	c call getin("thermochem",thermochem)
767	c write(,) " thermochem = ",thermochem
768
769	c if (.not.callthermos) then
770	c if (callconduct) then
771	c print*,'if callconduct is set, callthermos must be true'
772	c stop
773	c endif
774	c if (calleuv) then
775	c print*,'if calleuv is set, callthermos must be true'
776	c stop
777	c endif
778
779	!*************************************************************
780	!*************************************************************
781	! Tests of incompatibility:
782	!*************************************************************
783	!*************************************************************
784
785	if ((.not.tracer).and.(haze)) then
786	print*,'if haze are on, tracers must be on!'
787	stop
788	endif
789	if (callmolvis.and..not.callconduct) then
790	print*,'if callmolvis is set, callconduc must be true'
791	stop
792	endif
793	if (paleo.and..not.fast) then
794	print*,'if paleo is set, fast must be true'
795	stop
796	endif
797	if (fasthaze.and..not.fast) then
798	print*,'if fasthaze is set, fast must be true'
799	stop
800	endif
801	if (fasthaze.and..not.methane) then
802	print*,'if fasthaze is set, methane must be true'
803	stop
804	endif
805	if ((haze.or.fasthaze).and..not.methane) then
806	print*,'if haze/fasthaze is set, methane must be true'
807	stop
808	endif
809	if ((haze_proffix.or.haze_radproffix).and..not.aerohaze) then
810	print*,'for now, haze/rad proffix only works w aerohaze=T'
811	stop
812	endif
813	c if (callmoldiff) then
814	c print*,'if callmoldiff is set, callthermos must be true'
815	c stop
816	c endif
817	c if (thermochem) then
818	c print*,'if thermochem is set, callthermos must be true'
819	c stop
820	c endif
821	c endif
822
823	PRINT*,'--------------------------------------------'
824	PRINT*
825	PRINT*
826	ELSE
827	write(,)
828	write(,) 'Cannot read file callphys.def. Is it here ?'
829	stop
830	ENDIF
831
832	8000 FORMAT(t5,a12,l8)
833	8001 FORMAT(t5,a12,i8)
834
835	PRINT*
836	PRINT*,'inifis: daysec',daysec
837	PRINT*
838	PRINT*,'inifis: The radiative transfer is computed:'
839	PRINT*,' each ',iradia,' physical time-step'
840	!PRINT,' or each ',iradiadtphys,' seconds'
841	PRINT*
842
843	! --------------------------------------------------------------
844	! Managing the Longwave radiative transfer
845	! --------------------------------------------------------------
846
847	! In most cases, the run just use the following values :
848	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
849	callemis=.true.
850	! ilwd=10*int(daysec/dtphys) ! bug before 22/10/01
851	ilwd=10*int(daysec/dtphys)
852	ilwn=2
853	linear=.true.
854	ncouche=3
855	alphan=0.4
856	ilwb=2
857	semi=0
858
859	c$$$! BUT people working hard on the LW may want to read them in 'radia.def'
860	c$$$! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
861	c$$$
862	c$$$ OPEN(99,file='radia.def',status='old',form='formatted'
863	c$$$ . ,iostat=ierr)
864	c$$$ IF(ierr.EQ.0) THEN
865	c$$$ write(,) 'inifis: Reading radia.def !!!'
866	c$$$ READ(99,fmt='(a)') ch1
867	c$$$ READ(99,*) callemis
868	c$$$ WRITE(*,8000) ch1,callemis
869	c$$$
870	c$$$ READ(99,fmt='(a)') ch1
871	c$$$ READ(99,*) iradia
872	c$$$ WRITE(*,8001) ch1,iradia
873	c$$$
874	c$$$ READ(99,fmt='(a)') ch1
875	c$$$ READ(99,*) ilwd
876	c$$$ WRITE(*,8001) ch1,ilwd
877	c$$$
878	c$$$ READ(99,fmt='(a)') ch1
879	c$$$ READ(99,*) ilwn
880	c$$$ WRITE(*,8001) ch1,ilwn
881	c$$$
882	c$$$ READ(99,fmt='(a)') ch1
883	c$$$ READ(99,*) linear
884	c$$$ WRITE(*,8000) ch1,linear
885	c$$$
886	c$$$ READ(99,fmt='(a)') ch1
887	c$$$ READ(99,*) ncouche
888	c$$$ WRITE(*,8001) ch1,ncouche
889	c$$$
890	c$$$ READ(99,fmt='(a)') ch1
891	c$$$ READ(99,*) alphan
892	c$$$ WRITE(,) ch1,alphan
893	c$$$
894	c$$$ READ(99,fmt='(a)') ch1
895	c$$$ READ(99,*) ilwb
896	c$$$ WRITE(*,8001) ch1,ilwb
897	c$$$
898	c$$$
899	c$$$ READ(99,fmt='(a)') ch1
900	c$$$ READ(99,'(l1)') callg2d
901	c$$$ WRITE(*,8000) ch1,callg2d
902	c$$$
903	c$$$ READ(99,fmt='(a)') ch1
904	c$$$ READ(99,*) semi
905	c$$$ WRITE(,) ch1,semi
906	c$$$ end if
907	c$$$ CLOSE(99)
908
909	!-----------------------------------------------------------------------
910	! Some more initialization:
911	! ------------------------
912
913	CALL SCOPY(ngrid,plon,1,long,1)
914	CALL SCOPY(ngrid,plat,1,lati,1)
915	CALL SCOPY(ngrid,parea,1,area,1)
916	totarea=SSUM(ngridmx,area,1)
917	!! Computing realarea
918	DO ig=2,ngridmx-1
919	realarea(ig)=area(ig)
920	ENDDO
921	realarea(1)=area(1)*iim
922	realarea(ngridmx)=area(ngridmx)*iim
923
924	DO ig=1,ngrid
925	sinlat(ig)=sin(plat(ig))
926	coslat(ig)=cos(plat(ig))
927	sinlon(ig)=sin(plon(ig))
928	coslon(ig)=cos(plon(ig))
929	ENDDO
930
931	pi=2.*asin(1.) ! NB: pi is a common in comcstfi.h
932
933	! managing the tracers, and tests:
934	! -------------------------------
935
936	if(tracer) then
937
938	! when photochem is used, nqchem_min is the rank
939	! of the first chemical species
940
941	! Ehouarn: nqchem_min is now meaningless and no longer used
942	! nqchem_min = 1
943	if (photochem .or. callthermos) then
944	chem = .true.
945	end if
946
947	endif ! of if (tracer)
948
949	RETURN
950	END

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format