Context Navigation

conf_phys.F @ 3337

Last change on this file since 3337 was 3333, checked in by llange, 21 months ago
Mars PCM Fixing a bug in vdif_cd: a "residual", used as criterion to enter an iterative loop, was wrongly initialized. Hence, for some points, the algorithm does not go into the loop, and a wrong value of Cd, Ch was computed. Also some cleaning/small fixing with save variables with OMP. LL
File size: 49.2 KB

Line
1	MODULE conf_phys_mod
2
3	IMPLICIT NONE
4
5	CONTAINS
6
7	SUBROUTINE conf_phys(ngrid,nlayer,nq)
8
9	!=======================================================================
10	!
11	! purpose:
12	! -------
13	!
14	! Initialisation for the physical parametrisations
15	! flags (i.e. run-time options) of the Mars PCM
16	!-----------------------------------------------------------------------
17
18	USE ioipsl_getin_p_mod, ONLY : getin_p
19	use tracer_mod, only : nuice_sed, ccn_factor, nuiceco2_sed,
20	& nuice_ref,nuiceco2_ref
21	use surfdat_h, only: albedo_h2o_cap,albedo_h2o_frost,
22	& frost_albedo_threshold, inert_h2o_ice,
23	& frost_metam_threshold,old_wsublimation_scheme
24	use time_phylmdz_mod, only: ecritphy,day_step,iphysiq,ecritstart,
25	& daysec,dtphys
26	use dimradmars_mod, only: naerkind, name_iaer,
27	& ini_scatterers,tauvis
28	use datafile_mod, only: datadir
29	use wstats_mod, only: callstats
30	use writediagsoil_mod, only: diagsoil
31	use calchim_mod, only: ichemistry
32	use co2condens_mod, only: scavco2cond
33	use dust_param_mod, only: dustbin, doubleq, submicron, active,
34	& lifting, freedust, callddevil,
35	& dustscaling_mode,
36	& reff_driven_IRtoVIS_scenario
37	use aeropacity_mod, only: iddist, topdustref
38	USE mod_phys_lmdz_transfert_para, ONLY: bcast
39	USE paleoclimate_mod,ONLY: paleoclimate,albedo_perennialco2,
40	& lag_layer, include_waterbuoyancy
41	use microphys_h, only: mteta
42	use comsoil_h, only: adsorption_soil, choice_ads,ads_const_D,
43	& ads_massive_ice
44	use nonoro_gwd_mix_mod, only: calljliu_gwimix
45
46	IMPLICIT NONE
47
48	include "callkeys.h"
49
50	INTEGER,INTENT(IN) :: ngrid ! number of atmospheric columns
51	INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers
52	INTEGER,INTENT(IN) :: nq ! number of tracers
53
54	INTEGER ierr,j
55	character(len=20),parameter :: modname="conf_phys"
56
57	CHARACTER ch1*12
58	#ifndef MESOSCALE
59	! read in some parameters from "run.def" for physics,
60	! or shared between dynamics and physics.
61	ecritphy=240 ! default value
62	call getin_p("ecritphy",ecritphy) ! frequency of outputs in physics,
63	! in dynamical steps
64	day_step=960 ! default value
65	call getin_p("day_step",day_step) ! number of dynamical steps per day
66	iphysiq=20 ! default value
67	call getin_p("iphysiq",iphysiq) ! call physics every iphysiq dyn step
68	ecritstart=0 ! default value
69	call getin_p("ecritstart",ecritstart) ! write a restart every ecristart steps
70	#endif
71
72	! --------------------------------------------------------------
73	! Reading the "callphys.def" file controlling some key options
74	! --------------------------------------------------------------
75
76	!$OMP MASTER
77	! check that 'callphys.def' file is around
78	OPEN(99,file='callphys.def',status='old',form='formatted'
79	& ,iostat=ierr)
80	CLOSE(99)
81	!$OMP END MASTER
82	call bcast(ierr)
83	! ierr=0
84
85	IF(ierr.EQ.0) THEN
86	PRINT*
87	PRINT*
88	PRINT*,'--------------------------------------------'
89	PRINT*,' conf_phys: Parameters for the physics (callphys.def)'
90	PRINT*,'--------------------------------------------'
91
92	write(,) "Directory where external input files are:"
93	! default path is set in datafile_mod
94	call getin_p("datadir",datadir)
95	write(,) " datadir = ",trim(datadir)
96
97	write(,) "Initialize physics with startfi.nc file ?"
98	startphy_file=.true.
99	call getin_p("startphy_file",startphy_file)
100	write(,) "startphy_file", startphy_file
101
102	write(,) "Diurnal cycle ?"
103	write(,) "(if diurnal=False, diurnal averaged solar heating)"
104	diurnal=.true. ! default value
105	call getin_p("diurnal",diurnal)
106	write(,) " diurnal = ",diurnal
107
108	write(,) "Seasonal cycle ?"
109	write(,) "(if season=False, Ls stays constant, to value ",
110	& "set in 'start'"
111	season=.true. ! default value
112	call getin_p("season",season)
113	write(,) " season = ",season
114
115	write(,) "Write some extra output to the screen ?"
116	lwrite=.false. ! default value
117	call getin_p("lwrite",lwrite)
118	write(,) " lwrite = ",lwrite
119
120	write(,) "Save statistics in file stats.nc ?"
121	#ifdef MESOSCALE
122	callstats=.false. ! default value
123	#else
124	callstats=.true. ! default value
125	#endif
126	call getin_p("callstats",callstats)
127	write(,) " callstats = ",callstats
128
129	write(,) "Write sub-surface fields in file diagsoil.nc ?"
130	diagsoil=.false. ! default value
131	call getin_p("diagsoil",diagsoil)
132	write(,) " diagsoil = ",diagsoil
133
134	write(,) "Save EOF profiles in file 'profiles' for ",
135	& "Climate Database?"
136	calleofdump=.false. ! default value
137	call getin_p("calleofdump",calleofdump)
138	write(,) " calleofdump = ",calleofdump
139
140	write(,) "Dust scenario: 1=constant dust (read from startfi",
141	& " or set as tauvis); 2=Viking scenario; =3 MGS scenario,",
142	& "=6 cold (low dust) scenario; =7 warm (high dust) scenario ",
143	& "=24,25 ... 30 :Mars Year 24, ... or 30 from TES assimilation"
144	iaervar=3 ! default value
145	call getin_p("iaervar",iaervar)
146	write(,) " iaervar = ",iaervar
147
148	write(,) "Reference (visible) dust opacity at 610 Pa ",
149	& "(matters only if iaervar=1)"
150	! NB: default value of tauvis is set/read in startfi.nc file
151	call getin_p("tauvis",tauvis)
152	write(,) " tauvis = ",tauvis
153
154	write(,) "Dust vertical distribution:"
155	write(,) "(=1 top set by topdustref parameter;",
156	& " =2 Viking scenario; =3 MGS scenario)"
157	iddist=3 ! default value
158	call getin_p("iddist",iddist)
159	write(,) " iddist = ",iddist
160
161	write(,) "Dust top altitude (km). (Matters only if iddist=1)"
162	topdustref= 90.0 ! default value
163	call getin_p("topdustref",topdustref)
164	write(,) " topdustref = ",topdustref
165
166	write(,) "Prescribed surface thermal flux (H/(rho*cp),K m/s)"
167	tke_heat_flux=0. ! default value
168	call getin_p("tke_heat_flux",tke_heat_flux)
169	write(,) " tke_heat_flux = ",tke_heat_flux
170	write(,) " 0 means the usual schemes are computing"
171
172	write(,) "call radiative transfer ?"
173	callrad=.true. ! default value
174	call getin_p("callrad",callrad)
175	write(,) " callrad = ",callrad
176
177	write(,) "call slope insolation scheme ?",
178	& "(matters only if callrad=T)"
179	#ifdef MESOSCALE
180	callslope=.true. ! default value
181	#else
182	callslope=.false. ! default value (not supported yet)
183	#endif
184	call getin_p("callslope",callslope)
185	write(,) " callslope = ",callslope
186
187	write(,) "call NLTE radiative schemes ?",
188	& "(matters only if callrad=T)"
189	callnlte=.false. ! default value
190	call getin_p("callnlte",callnlte)
191	write(,) " callnlte = ",callnlte
192
193	nltemodel=0 !default value
194	write(,) "NLTE model?"
195	write(,) "0 -> old model, static O"
196	write(,) "1 -> old model, dynamic O"
197	write(,) "2 -> new model"
198	write(,) "(matters only if callnlte=T)"
199	call getin_p("nltemodel",nltemodel)
200	write(,) " nltemodel = ",nltemodel
201
202	write(,) "call CO2 NIR absorption ?",
203	& "(matters only if callrad=T)"
204	callnirco2=.false. ! default value
205	call getin_p("callnirco2",callnirco2)
206	write(,) " callnirco2 = ",callnirco2
207
208	write(,) "New NIR NLTE correction ?",
209	$ "0-> old model (no correction)",
210	$ "1-> new correction",
211	$ "(matters only if callnirco2=T)"
212	#ifdef MESOSCALE
213	nircorr=0 !default value. this is OK below 60 km.
214	#else
215	nircorr=0 !default value
216	#endif
217	call getin_p("nircorr",nircorr)
218	write(,) " nircorr = ",nircorr
219
220	write(,) "call turbulent vertical diffusion ?"
221	calldifv=.true. ! default value
222	call getin_p("calldifv",calldifv)
223	write(,) " calldifv = ",calldifv
224
225	write(,) "call thermals ?"
226	calltherm=.false. ! default value
227	call getin_p("calltherm",calltherm)
228	write(,) " calltherm = ",calltherm
229
230	write(,) "call convective adjustment ?"
231	calladj=.true. ! default value
232	call getin_p("calladj",calladj)
233	write(,) " calladj = ",calladj
234
235	if (calltherm .and. calladj) then
236	print*,'!!! PLEASE NOTE !!!'
237	print*,'convective adjustment is on'
238	print*,'but since thermal plume model is on'
239	print*,'convadj is only activated above the PBL'
240	endif
241
242	write(,) "used latest version of yamada scheme?"
243	callyamada4=.true. ! default value
244	call getin_p("callyamada4",callyamada4)
245	write(,) " callyamada4 = ",callyamada4
246
247	write(,) "used latest version of ATKE scheme?"
248	callatke =.false. ! default value
249	call getin_p("callatke",callatke)
250	write(,) " callatke = ",callatke
251
252	if(callyamada4.and.callatke) then
253	print*,"You can't use both yamada and atke"
254	print*,"Choose either Yamada or ATKE"
255	call abort_physic(modname,
256	& "Can't use both yamada and ATKE",1)
257	endif
258
259	if (calltherm .and. .not.callyamada4) then
260	print*,'!!!! WARNING WARNING WARNING !!!!'
261	print*,'if calltherm=T we strongly advise that '
262	print*,'you set the flag callyamada4 to T '
263	print*,'!!!! WARNING WARNING WARNING !!!!'
264	endif
265
266	write(,) "call Richardson-based surface layer ?"
267	callrichsl=.false. ! default value
268	call getin_p("callrichsl",callrichsl)
269	write(,) " callrichsl = ",callrichsl
270
271	if (calltherm .and. .not.callrichsl) then
272	print*,'WARNING WARNING WARNING'
273	print*,'if calltherm=T we strongly advise that '
274	print*,'you use the new surface layer scheme '
275	print*,'by setting callrichsl=T '
276	endif
277
278	if (calladj .and. callrichsl .and. (.not. calltherm)) then
279	print*,'You should not be calling the convective adjustment
280	& scheme with the Richardson surface-layer and without the thermals
281	&. This approach is not
282	& physically consistent and can lead to unrealistic friction
283	& values.'
284	print*,'If you want to use the Ri. surface-layer, either
285	& activate thermals OR de-activate the convective adjustment.'
286	call abort_physic(modname,
287	& "Richardson layer must be used with thermals",1)
288	endif
289
290	write(,) "call CO2 condensation ?"
291	callcond=.true. ! default value
292	call getin_p("callcond",callcond)
293	write(,) " callcond = ",callcond
294
295	write(,)"call thermal conduction in the soil ?"
296	callsoil=.true. ! default value
297	call getin_p("callsoil",callsoil)
298	write(,) " callsoil = ",callsoil
299
300
301	write(,)"call Lott's gravity wave/subgrid topography ",
302	& "scheme ?"
303	calllott=.true. ! default value
304	call getin_p("calllott",calllott)
305	write(,)" calllott = ",calllott
306
307	write(,)"call Lott's non-oro GWs parameterisation ",
308	& "scheme ?"
309	calllott_nonoro=.false. ! default value
310	call getin_p("calllott_nonoro",calllott_nonoro)
311	write(,)" calllott_nonoro = ",calllott_nonoro
312
313	write(,)"call jliu's nonoro GWs-induced mixing ",
314	& "scheme ?"
315	calljliu_gwimix=.false. ! default value
316	call getin_p("calljliu_gwimix",calljliu_gwimix)
317	write(,)" calljliu_gwimix = ",calljliu_gwimix
318
319
320	! rocket dust storm injection scheme
321	write(,)"call rocket dust storm parametrization"
322	rdstorm=.false. ! default value
323	call getin_p("rdstorm",rdstorm)
324	write(,)" rdstorm = ",rdstorm
325	! rocket dust storm detrainment coefficient
326	coeff_detrainment=0.02 ! default value
327	call getin_p("coeff_detrainment",coeff_detrainment)
328	write(,)" coeff_detrainment = ",coeff_detrainment
329
330	! entrainment by mountain top dust flows scheme
331	write(,)"call mountain top dust flows parametrization"
332	topflows=.false. ! default value
333	call getin_p("topflows",topflows)
334	write(,)" topflows = ",topflows
335
336	! latent heat release from ground water ice sublimation/condensation
337	write(,)"latent heat release during sublimation",
338	& " /condensation of ground water ice"
339	latentheat_surfwater=.true. ! default value
340	call getin_p("latentheat_surfwater",latentheat_surfwater)
341	write(,)" latentheat_surfwater = ",latentheat_surfwater
342
343	write(,)"rad.transfer is computed every iradia",
344	& " physical timestep"
345	iradia=1 ! default value
346	call getin_p("iradia",iradia)
347	write(,)" iradia = ",iradia
348
349
350	write(,)"Output of the exchange coefficient mattrix ?",
351	& "(for diagnostics only)"
352	callg2d=.false. ! default value
353	call getin_p("callg2d",callg2d)
354	write(,)" callg2d = ",callg2d
355
356	write(,)"Rayleigh scattering : (should be .false. for now)"
357	rayleigh=.false.
358	call getin_p("rayleigh",rayleigh)
359	write(,)" rayleigh = ",rayleigh
360
361	! PALEOCLIMATE
362
363	write(,)"Using lag layer??"
364	lag_layer=.false.
365	call getin_p("lag_layer",lag_layer)
366	write(,) "lag_layer = ", lag_layer
367
368	write(,)"Is it paleoclimate run?"
369	paleoclimate=.false. ! default value
370	call getin_p("paleoclimate",paleoclimate)
371	write(,)"paleoclimate = ",paleoclimate
372
373	write(,)"Albedo for perennial CO2 ice?"
374	albedo_perennialco2 = 0.85 ! default value
375	call getin_p("albedo_perennialco2",albedo_perennialco2)
376	write(,)"albedo_perennialco2 = ",albedo_perennialco2
377
378	write(,)"Include water buoyancy effect??"
379	include_waterbuoyancy=.false.
380	call getin_p("include_waterbuoyancy",include_waterbuoyancy)
381	write(,) "include_waterbuoyancy = ", include_waterbuoyancy
382
383	if (include_waterbuoyancy.and.(.not.paleoclimate)) then
384	print*,'You should not include the water buoyancy effect
385	& on current climate (model not tunned with this effect'
386	call abort_physic(modname,
387	& "include_waterbuoyancy must be used with paleoclimate",1)
388	endif
389
390	! TRACERS:
391
392	! dustbin
393	write(,)"Transported dust ? (if >0, use 'dustbin' dust bins)"
394	dustbin=0 ! default value
395	call getin_p("dustbin",dustbin)
396	write(,)" dustbin = ",dustbin
397	! active
398	write(,)"Radiatively active dust ? (matters if dustbin>0)"
399	active=.false. ! default value
400	call getin_p("active",active)
401	write(,)" active = ",active
402
403	! Test of incompatibility:
404	! if active is used, then dustbin should be > 0
405
406	if (active.and.(dustbin.lt.1)) then
407	print*,'if active is used, then dustbin should > 0'
408	call abort_physic(modname,
409	& "active option requires dustbin < 0",1)
410	endif
411	! doubleq
412	write(,)"use mass and number mixing ratios to predict",
413	& " dust size ?"
414	doubleq=.false. ! default value
415	call getin_p("doubleq",doubleq)
416	write(,)" doubleq = ",doubleq
417	! submicron
418	submicron=.false. ! default value
419	call getin_p("submicron",submicron)
420	write(,)" submicron = ",submicron
421
422	! Test of incompatibility:
423	! if doubleq is used, then dustbin should be 2
424
425	if (doubleq.and.(dustbin.ne.2)) then
426	print*,'if doubleq is used, then dustbin should be 2'
427	call abort_physic(modname,
428	& "doubleq option requires dustbin = 2",1)
429	endif
430	if (doubleq.and.submicron.and.(nq.LT.3)) then
431	print*,'If doubleq is used with a submicron tracer,'
432	print*,' then the number of tracers has to be'
433	print*,' larger than 3.'
434	call abort_physic(modname,
435	& "submicron option requires dustbin > 2",1)
436	endif
437
438	! lifting
439	write(,)"dust lifted by GCM surface winds ?"
440	lifting=.false. ! default value
441	call getin_p("lifting",lifting)
442	write(,)" lifting = ",lifting
443
444	! Test of incompatibility:
445	! if lifting is used, then dustbin should be > 0
446
447	if (lifting.and.(dustbin.lt.1)) then
448	print*,'if lifting is used, then dustbin should > 0'
449	call abort_physic(modname,
450	& "lifting option requires dustbin > 0",1)
451	endif
452
453	! dust injection scheme
454	dustinjection=0 ! default: no injection scheme
455	call getin_p("dustinjection",dustinjection)
456	write(,)" dustinjection = ",dustinjection
457	! dust injection scheme coefficient
458	coeff_injection=0.25 ! default value
459	call getin_p("coeff_injection",coeff_injection)
460	write(,)" coeff_in,jection = ",coeff_injection
461	! timing for dust injection
462	ti_injection=0. ! default value
463	tf_injection=24. ! default value
464	call getin_p("ti_injection",ti_injection)
465	write(,)" ti_injection = ",ti_injection
466	call getin_p("tf_injection",tf_injection)
467	write(,)" tf_injection = ",tf_injection
468
469	! free evolving dust
470	! freedust=true just says that there is no lifting and no dust opacity scaling.
471	write(,)"dust lifted by GCM surface winds ?"
472	freedust=.false. ! default value
473	call getin_p("freedust",freedust)
474	write(,)" freedust = ",freedust
475	if (freedust.and..not.doubleq) then
476	print*,'freedust should be used with doubleq !'
477	call abort_physic(modname,
478	& "freedust option requires doubleq",1)
479	endif
480
481	! dust rescaling mode (if any)
482	! =0 --> freedust, tauscaling=1
483	! =1 --> GCM5.3-like, tauscaling computed wrt tau_pref_scenario
484	! =2 --> tauscaling=1 but dust_rad_adjust (rescaling for the radiative transfer only)
485	if (freedust) then
486	dustscaling_mode=0
487	else
488	dustscaling_mode=1 ! GCMv5.3 style
489	endif
490	call getin_p("dustscaling_mode",dustscaling_mode)
491	write(,) "dustscaling_mode=",dustscaling_mode
492	if ((dustscaling_mode.eq.1).and.freedust) then
493	print*,'freedust and dustscaling_mode=1 are incompatible !'
494	call abort_physic(modname,
495	& "freedust and dustscaling_mode=1 are incompatible",1)
496	endif
497
498	#ifndef MESOSCALE
499	! this test is valid in GCM case
500	! ... not in mesoscale case, for we want to activate mesoscale lifting
501	if (freedust.and.dustinjection.eq.0)then
502	if(lifting) then
503	print*,'if freedust is used and dustinjection = 0,
504	& then lifting should not be used'
505	call abort_physic(modname,
506	& "freedust option with dustinjection = 0"//
507	& " requires lifting to be false",1)
508	endif
509	endif
510	#endif
511	if (dustinjection.eq.1)then
512	if(.not.lifting) then
513	print*,"if dustinjection=1, then lifting should be true"
514	call abort_physic(modname,
515	& "dustinjection=1 requires lifting",1)
516	endif
517	if(.not.freedust) then
518	print*,"if dustinjection=1, then freedust should be true"
519	call abort_physic(modname,
520	& "dustinjection=1 requires freedust",1)
521	endif
522	endif
523	! rocket dust storm and entrainment by top flows
524	! Test of incompatibility:
525	! if rdstorm or topflows is used, then doubleq should be true
526	if ((rdstorm.or.topflows).and..not.doubleq) then
527	print*,'if rdstorm or topflows is used, then doubleq
528	& should be used !'
529	call abort_physic(modname,
530	& "rdstorm or topflows requires doubleq",1)
531	endif
532	if ((rdstorm.or.topflows).and..not.active) then
533	print*,'if rdstorm or topflows is used, then active
534	& should be used !'
535	call abort_physic(modname,
536	& "rdstorm or topflows requires active",1)
537	endif
538	if (rdstorm.and..not.lifting) then
539	print*,'if rdstorm is used, then lifting
540	& should be used !'
541	call abort_physic(modname,
542	& "rdstorm requires lifting",1)
543	endif
544	if ((rdstorm.or.topflows).and..not.freedust) then
545	print*,'if rdstorm or topflows is used, then freedust
546	& should be used !'
547	call abort_physic(modname,
548	& "rdstorm or topflows requires freedust",1)
549	endif
550	if (rdstorm.and.(dustinjection.eq.0)) then
551	print*,'if rdstorm is used, then dustinjection
552	& should be used !'
553	call abort_physic(modname,
554	& "rdstorm requires dustinjection",1)
555	endif
556	! Dust IR opacity
557	write(,)" Wavelength for infrared opacity of dust ?"
558	write(,)" Choices are:"
559	write(,)" tes --- > 9.3 microns [default]"
560	write(,)" mcs --- > 21.6 microns"
561	!
562	! WARNING WARNING WARNING WARNING WARNING WARNING
563	!
564	! BEFORE ADDING A NEW VALUE, BE SURE THAT THE
565	! CORRESPONDING WAVELENGTH IS IN THE LOOKUP TABLE,
566	! OR AT LEAST NO TO FAR, TO AVOID FALLACIOUS INTERPOLATIONS.
567	!
568	dustiropacity="tes" !default value
569	call getin_p("dustiropacity",dustiropacity)
570	write(,)" dustiropacity = ",trim(dustiropacity)
571	select case (trim(dustiropacity))
572	case ("tes")
573	dustrefir = 9.3E-6
574	case ("mcs")
575	dustrefir = 21.6E-6
576	case default
577	write(,) trim(dustiropacity),
578	& " is not a valid option for dustiropacity"
579	call abort_physic(modname,
580	& "invalid dustiropacity option value",1)
581	end select
582	! Dust scenario IR to VIS conversion
583	write(,)"Use an IR to VIS conversion coefficient"
584	write(,)"for the dust scenario, that is dependent"
585	write(,)"on the GCM dust effective radius,"
586	write(,)"instead of a fixed 2.6 coefficient ?"
587	reff_driven_IRtoVIS_scenario=.false. !default value
588	call getin_p("reff_driven_IRtoVIS_scenario",
589	& reff_driven_IRtoVIS_scenario)
590	write(,)" reff_driven_IRtoVIS_scenario = ",
591	& reff_driven_IRtoVIS_scenario
592	! Test of incompatibility:
593	! if reff_driven_IRtoVIS_scenario=.true.,
594	! dustrefir must be 9.3E-6 = scenarios' wavelength
595	if (reff_driven_IRtoVIS_scenario .and.
596	& (dustrefir.ne.9.3E-6)) then
597	print*,'if reff_driven_IRtoVIS_scenario is used, then '//
598	& 'the GCM IR reference wavelength should be the one '//
599	& 'of the scenarios (dustiropacity=tes)'
600	call abort_physic(modname,
601	& "reff_driven_IRtoVIS_scenario requires tes wavelength",1)
602	endif
603
604	! callddevil
605	write(,)" dust lifted by dust devils ?"
606	callddevil=.false. !default value
607	call getin_p("callddevil",callddevil)
608	write(,)" callddevil = ",callddevil
609
610	! Test of incompatibility:
611	! if dustdevil is used, then dustbin should be > 0
612	if (callddevil.and.(dustbin.lt.1)) then
613	print*,'if dustdevil is used, then dustbin should > 0'
614	call abort_physic(modname,
615	& "callddevil requires dustbin > 0",1)
616	endif
617
618	! sedimentation
619	write(,) "Gravitationnal sedimentation ?"
620	sedimentation=.true. ! default value
621	call getin_p("sedimentation",sedimentation)
622	write(,) " sedimentation = ",sedimentation
623	! activice
624	write(,) "Radiatively active transported atmospheric ",
625	& "water ice ?"
626	activice=.false. ! default value
627	call getin_p("activice",activice)
628	write(,) " activice = ",activice
629	! water
630	write(,) "Compute water cycle ?"
631	water=.false. ! default value
632	call getin_p("water",water)
633	write(,) " water = ",water
634	! hdo
635	write(,) "Compute hdo cycle ?"
636	hdo=.false. ! default value
637	call getin_p("hdo",hdo)
638	write(,) " hdo = ",hdo
639
640	write(,) "Use fractionation for hdo?"
641	hdofrac=.true. ! default value
642	call getin_p("hdofrac",hdofrac)
643	write(,) " hdofrac = ",hdofrac
644
645	! Activeco2ice
646	write(,) "Radiatively active transported atmospheric ",
647	& "Co2 ice ?"
648	activeco2ice=.false. ! default value
649	call getin_p("activeco2ice",activeco2ice)
650	write(,) " activeco2ice = ",activeco2ice
651
652	! sub-grid cloud fraction: fixed
653	write(,) "Fixed cloud fraction?"
654	CLFfixval=1.0 ! default value
655	call getin_p("CLFfixval",CLFfixval)
656	write(,) "CLFfixval=",CLFfixval
657	! sub-grid cloud fraction: varying
658	write(,) "Use partial nebulosity?"
659	CLFvarying=.false. ! default value
660	call getin_p("CLFvarying",CLFvarying)
661	write(,)"CLFvarying=",CLFvarying
662
663	!CO2 clouds scheme?
664	write(,) "Compute CO2 clouds (implies microphysical scheme)?"
665	co2clouds=.false. ! default value
666	call getin_p("co2clouds",co2clouds)
667	write(,) " co2clouds = ",co2clouds
668	!Can water ice particles serve as CCN for CO2clouds
669	write(,) "Use water ice as CO2 clouds CCN ?"
670	co2useh2o=.false. ! default value
671	call getin_p("co2useh2o",co2useh2o)
672	write(,) " co2useh2o = ",co2useh2o
673	!Do we allow a supply of meteoritic paricles to serve as CO2 ice CCN?
674	write(,) "Supply meteoritic particle for CO2 clouds ?"
675	meteo_flux=.false. !Default value
676	call getin_p("meteo_flux",meteo_flux)
677	write(,) " meteo_flux = ",meteo_flux
678	!Do we allow a sub-grid temperature distribution for the CO2 microphysics
679	write(,) "sub-grid temperature distribution for CO2 clouds?"
680	CLFvaryingCO2=.false. !Default value
681	call getin_p("CLFvaryingCO2",CLFvaryingCO2)
682	write(,) " CLFvaryingCO2 = ",CLFvaryingCO2
683	!Amplitude of the sub-grid temperature distribution for the CO2 microphysics
684	write(,) "sub-grid temperature amplitude for CO2 clouds?"
685	spantCO2=0 !Default value
686	call getin_p("spantCO2",spantCO2)
687	write(,) " spantCO2 = ",spantCO2
688	!Do you want to filter the sub-grid T distribution by a Saturation index?
689	write(,) "filter sub-grid temperature by Saturation index?"
690	satindexco2=.true.
691	call getin_p("satindexco2",satindexco2)
692	write(,) " satindexco2 = ",satindexco2
693
694
695	! thermal inertia feedback
696	write(,) "Activate the thermal inertia feedback ?"
697	surfaceice_tifeedback=.false. ! default value
698	call getin_p("surfaceice_tifeedback",surfaceice_tifeedback)
699	write(,) " surfaceice_tifeedback = ",surfaceice_tifeedback
700
701	! Test of incompatibility:
702
703	if (surfaceice_tifeedback.and..not.water) then
704	print*,'if surfaceice_tifeedback is used,'
705	print*,'water should be used too'
706	call abort_physic(modname,
707	& "surfaceice_tifeedback requires water",1)
708	endif
709
710	if (surfaceice_tifeedback.and..not.callsoil) then
711	print*,'if surfaceice_tifeedback is used,'
712	print*,'callsoil should be used too'
713	call abort_physic(modname,
714	& "surfaceice_tifeedback requires callsoil",1)
715	endif
716
717	if (activice.and..not.water) then
718	print*,'if activice is used, water should be used too'
719	call abort_physic(modname,
720	& "activeice requires water",1)
721	endif
722
723	if (hdo.and..not.water) then
724	print*,'if hdo is used, water should be used too'
725	call abort_physic(modname,
726	& "hd2 requires tracer",1)
727	endif
728
729
730	if (activeco2ice.and..not.co2clouds) then
731	print*,'if activeco2ice is used, co2clouds should be used too'
732	call abort_physic(modname,
733	& "activeco2ice requires co2clouds",1)
734	endif
735
736	! water ice clouds effective variance distribution for sedimentaion
737	write(,) "Sed effective variance for water ice clouds ?"
738	nuice_sed=0.45
739	call getin_p("nuice_sed",nuice_sed)
740	write(,) "water_param nueff Sedimentation:", nuice_sed
741
742	write(,) "Sed effective variance for CO2 clouds ?"
743	nuiceco2_sed=0.45
744	call getin_p("nuiceco2_sed",nuiceco2_sed)
745	write(,) "CO2 nueff Sedimentation:", nuiceco2_sed
746
747	write(,) "REF effective variance for CO2 clouds ?"
748	nuiceco2_ref=0.45
749	call getin_p("nuiceco2_ref",nuiceco2_ref)
750	write(,) "CO2 nueff Sedimentation:", nuiceco2_ref
751
752	write(,) "REF effective variance for water clouds ?"
753	nuice_ref=0.1
754	call getin_p("nuice_ref",nuice_ref)
755	write(,) "H2O nueff Sedimentation:", nuice_ref
756
757
758	! ccn factor if no scavenging
759	write(,) "water param CCN reduc. factor ?"
760	ccn_factor = 4.5
761	call getin_p("ccn_factor",ccn_factor)
762	write(,)" ccn_factor = ",ccn_factor
763	write(,)"Careful: only used when microphys=F, otherwise"
764	write(,)"the contact parameter is used instead;"
765
766	! microphys
767	write(,)"Microphysical scheme for water-ice clouds?"
768	microphys=.false. ! default value
769	call getin_p("microphys",microphys)
770	write(,)" microphys = ",microphys
771
772	! supersat
773	write(,)"Allow super-saturation of water vapor?"
774	supersat=.true. ! default value
775	call getin_p("supersat",supersat)
776	write(,)"supersat = ",supersat
777
778	! microphysical parameter contact
779	write(,) "water contact parameter ?"
780	mteta = 0.95 ! default value
781	temp_dependent_m = .false. ! default value
782	call getin_p("temp_dependent_m",temp_dependent_m)
783	if (temp_dependent_m) then !(JN 2023)
784	print*,'You have chosen a temperature-dependent water'
785	print*,'contact parameter ! From Maattanen et al. 2014'
786	else if (.not.temp_dependent_m) then
787	print*,'Water contact parameter is constant'
788	call getin_p("mteta",mteta)
789	write(,) "mteta = ", mteta
790	endif
791
792	! Adaptative timestep for cloud microphysics (JN 2023)
793	write(,)"Adaptative timestep for cloud",
794	& " microphysics ? (default is false)"
795	cloud_adapt_ts=.false. ! default value
796	call getin_p("cloud_adapt_ts",cloud_adapt_ts)
797	write(,)"cloud_adapt_ts= ",cloud_adapt_ts
798
799	! Test of incompatibility:
800	! If you want the adaptative timestep, you should use the
801	! temperature dependent contact parameter (otherwise the
802	! global water cycle is nuts !)
803	! However one can use the temperature dependent contact parameter
804	! without the adaptative timestep (MCD6.1 configuration)
805
806	if (cloud_adapt_ts.and..not.temp_dependent_m) then
807	print*,'Water cycle v6 : if cloud_adapt_ts is used'
808	print*,'then temp_dependent_m must be used!'
809	print*,'Otherwise the water cycle is unrealistic.'
810	call abort_physic(modname,
811	& "cloud_adapt_ts requires temp_dependent_m",1)
812	endif
813	! scavenging
814	write(,)"Dust scavenging by H2O/CO2 snowfall ?"
815	scavenging=.false. ! default value
816	call getin_p("scavenging",scavenging)
817	write(,)" scavenging = ",scavenging
818
819
820	! Test of incompatibility:
821	! if scavenging is used, then dustbin should be > 0
822
823	if ((microphys.and..not.doubleq).or.
824	& (microphys.and..not.water)) then
825	print*,'if microphys is used, then doubleq,'
826	print*,'and water must be used!'
827	call abort_physic(modname,
828	& "microphys requires water and doubleq",1)
829	endif
830	if (microphys.and..not.scavenging) then
831	print*,''
832	print*,'----------------WARNING-----------------'
833	print*,'microphys is used without scavenging !!!'
834	print*,'----------------WARNING-----------------'
835	print*,''
836	endif
837
838	if ((scavenging.and..not.microphys).or.
839	& (scavenging.and.(dustbin.lt.1)))then
840	print*,'if scavenging is used, then microphys'
841	print*,'must be used!'
842	call abort_physic(modname,
843	& "scavenging requires microphys",1)
844	endif
845
846	! Instantaneous scavenging by CO2
847	! -> expected to be replaced by scavenging with microphysics (flag scavenging) one day
848	write(,)"Dust scavenging by instantaneous CO2 snowfall ?"
849	scavco2cond=.false. ! default value
850	call getin_p("scavco2cond",scavco2cond)
851	write(,)" scavco2cond = ",scavco2cond
852	! Test of incompatibility:
853	! if scavco2cond is used, then dustbin should be > 0
854	if (scavco2cond.and.(dustbin.lt.1))then
855	print*,'if scavco2cond is used, then dustbin should be > 0'
856	call abort_physic(modname,
857	& "scavco2cond requires dustbin > 0",1)
858	endif
859	! if co2clouds is used, then there is no need for scavco2cond
860	if (co2clouds.and.scavco2cond) then
861	print*,''
862	print*,'----------------WARNING-----------------'
863	print*,' microphys scavenging is used so '
864	print*,' no need for scavco2cond !!! '
865	print*,'----------------WARNING-----------------'
866	print*,''
867	call abort_physic(modname,
868	& "incompatible co2cloud and scavco2cond options",1)
869	endif
870
871	! Test of incompatibility:
872
873	write(,) "Permanent water caps at poles ?",
874	& " .true. is RECOMMENDED"
875	write(,) "(with .true., North cap is a source of water ",
876	& "and South pole is a cold trap)"
877	caps=.true. ! default value
878	call getin_p("caps",caps)
879	write(,) " caps = ",caps
880
881	! JN : now separated between albedo_h2o_cap and
882	! albedo_h2o_frost. Retrocompatible with old
883	! callphys.def with albedo_h2o_ice
884	write(,) "water ice albedo ? Old settings use ",
885	& "albedo_h2o_ice, new settings use ",
886	& "albedo_h2o_cap and albedo_h2o_frost "
887	albedo_h2o_cap=0.35
888	albedo_h2o_frost=0.35
889	call getin_p("albedo_h2o_ice",albedo_h2o_cap)
890	albedo_h2o_frost=albedo_h2o_cap
891	call getin_p("albedo_h2o_cap",albedo_h2o_cap)
892	write(,) " albedo_h2o_cap = ",albedo_h2o_cap
893	call getin_p("albedo_h2o_frost",albedo_h2o_frost)
894	write(,) " albedo_h2o_frost = ",albedo_h2o_frost
895
896	! Northern polar cap albedo (JN 2021)
897	write(,)"Watercaptag albedo is unchanged by water frost",
898	& " deposition (default is false)"
899	cst_cap_albedo=.false. ! default value
900	call getin_p("cst_cap_albedo",cst_cap_albedo)
901	write(,)"cst_cap_albedo = ",cst_cap_albedo
902
903	! Watercap evolution & refill (with discriminated albedos) (JN 2021)
904	write(,)"Watercap is replenished by water frost",
905	& " accumulation (default is false)"
906	refill_watercap=.false. ! default value
907	call getin_p("refill_watercap",refill_watercap)
908	write(,)"refill_watercap= ",refill_watercap
909	! frost thickness threshold for refill_watercap (ice metamorphism)
910	write(,) "frost thickness threshold for metamorphism ?",
911	& "ie converted into watercap",
912	& "only if refill_watercap is .true."
913	frost_metam_threshold=0.05 ! (i.e 0.05 kg.m-2)
914	call getin_p("frost_metam_threshold",
915	& frost_metam_threshold)
916	write(,) " frost_metam_threshold = ",
917	& frost_metam_threshold
918
919
920	! inert_h2o_ice
921	write(,) "water ice thermal inertia ?"
922	inert_h2o_ice=2400 ! (J.m^-2.K^-1.s^-1/2)
923	call getin_p("inert_h2o_ice",inert_h2o_ice)
924	write(,) " inert_h2o_ice = ",inert_h2o_ice
925	! frost_albedo_threshold
926	write(,) "frost thickness threshold for albedo ?"
927	frost_albedo_threshold=0.005 ! 5.4 mic (i.e 0.005 kg.m-2)
928	call getin_p("frost_albedo_threshold",
929	& frost_albedo_threshold)
930	write(,) " frost_albedo_threshold = ",
931	& frost_albedo_threshold
932
933	! TMP: old_wsublimation_scheme
934	write(,) "Old water sublimation scheme?"
935	old_wsublimation_scheme = .true.
936	call getin_p("old_wsublimation_scheme",old_wsublimation_scheme)
937	write(,) "old_wsublimation_scheme",old_wsublimation_scheme
938
939	! call Titus crocus line -- DEFAULT IS NONE
940	write(,) "Titus crocus line ?"
941	tituscap=.false. ! default value
942	call getin_p("tituscap",tituscap)
943	write(,) "tituscap",tituscap
944
945	! Chemistry:
946	write(,) "photochemistry: include chemical species"
947	photochem=.false. ! default value
948	call getin_p("photochem",photochem)
949	write(,) " photochem = ",photochem
950
951	write(,) "Compute chemistry (if photochem is .true.)",
952	& "every ichemistry physics step (default: ichemistry=1)"
953	ichemistry=1
954	call getin_p("ichemistry",ichemistry)
955	write(,) " ichemistry = ",ichemistry
956
957
958	! SCATTERERS
959	write(,) "how many scatterers?"
960	naerkind=1 ! default value
961	call getin_p("naerkind",naerkind)
962	write(,)" naerkind = ",naerkind
963
964	! Test of incompatibility
965	c Logical tests for radiatively active water-ice clouds:
966	IF ( (activice.AND.(.NOT.water)).OR.
967	& (activice.AND.(naerkind.LT.2)) ) THEN
968	WRITE(,) 'If activice is TRUE, water has to be set'
969	WRITE(,) 'to TRUE, and "naerkind" must be at least'
970	WRITE(,) 'equal to 2.'
971	call abort_physic(modname,
972	& "radiatively active dust and water"//
973	& " require naerkind > 1",1)
974	ENDIF
975
976	!------------------------------------------
977	!------------------------------------------
978	! once naerkind is known allocate arrays
979	! -- we do it here and not in phys_var_init
980	! -- because we need to know naerkind first
981	CALL ini_scatterers(ngrid,nlayer)
982	!------------------------------------------
983	!------------------------------------------
984
985
986	c Please name the different scatterers here ----------------
987	name_iaer(1) = "dust_conrath" !! default choice is good old Conrath profile
988	IF (doubleq.AND.active) name_iaer(1) = "dust_doubleq" !! two-moment scheme
989
990	if (nq.gt.1) then
991	! trick to avoid problems compiling with 1 tracer
992	! and picky compilers who know name_iaer(2) is out of bounds
993	j=2
994	IF (rdstorm.AND..NOT.activice.AND..NOT.topflows) then
995	name_iaer(j) = "stormdust_doubleq" !! storm dust two-moment scheme
996	j = j+1
997	END IF
998
999	IF (rdstorm.AND.water.AND.activice.AND..NOT.topflows) then
1000	name_iaer(j) = "stormdust_doubleq"
1001	j = j+1
1002	END IF
1003
1004	IF (topflows.AND..NOT.activice.AND..NOT.rdstorm) then
1005	name_iaer(j) = "topdust_doubleq" !! storm dust two-moment scheme
1006	j = j+1
1007	END IF
1008
1009	IF (topflows.AND.water.AND.activice.AND..NOT.rdstorm) then
1010	name_iaer(j) = "topdust_doubleq"
1011	j = j+1
1012	END IF
1013
1014	IF (rdstorm.AND.topflows.AND..NOT.activice) THEN
1015	name_iaer(j) = "stormdust_doubleq"
1016	name_iaer(j+1) = "topdust_doubleq"
1017	j = j+2
1018	ENDIF
1019
1020	IF (rdstorm.AND.topflows.AND.water.AND.activice) THEN
1021	name_iaer(j) = "stormdust_doubleq"
1022	name_iaer(j+1) = "topdust_doubleq"
1023	j = j+2
1024	ENDIF
1025
1026	IF (water.AND.activice) then
1027	name_iaer(j) = "h2o_ice" !! radiatively-active clouds
1028	j = j+1
1029	END IF
1030
1031	IF (co2clouds.AND.activeco2ice) then
1032	name_iaer(j) = "co2_ice" !! radiatively-active co2 clouds
1033	j = j+1
1034	ENDIF
1035
1036	IF (submicron.AND.active) then
1037	name_iaer(j) = "dust_submicron" !! JBM experimental stuff
1038	j = j+1
1039	ENDIF
1040	endif ! of if (nq.gt.1)
1041	c ----------------------------------------------------------
1042
1043	! Adsorption
1044	adsorption_soil = .false. ! default value
1045	call getin_p("adsorption_soil",adsorption_soil)
1046	if (adsorption_soil .and. (.not. water)) then
1047	write(,)"Adsorption can be run only if water = True"
1048	call abort_physic(modname,
1049	& "Adsorption must be used with water = true",1)
1050	endif
1051
1052	choice_ads = 0 ! default value (no adsorption)
1053	if (adsorption_soil) choice_ads = 1
1054	call getin_p("choice_ads",choice_ads)
1055
1056	ads_const_D = .false.
1057	call getin_p("ads_const_D",ads_const_D)
1058
1059	ads_massive_ice = .false.
1060	call getin_p("ads_massive_ice",ads_massive_ice)
1061
1062	poreice_tifeedback = .false.
1063	call getin_p("poreice_tifeedback",poreice_tifeedback)
1064	write(,) 'poreice_tifeedback=',poreice_tifeedback
1065	if (poreice_tifeedback .and. (.not. adsorption_soil)) then
1066	write(,)"Pore ice TI feedback can be run only if
1067	& adsorption_soil = True"
1068	call abort_physic(modname,
1069	& "Pore ice TI feedback be used with adsorption_soil
1070	& = true",1)
1071	endif
1072	c ----------------------------------------------------------
1073
1074	! THERMOSPHERE
1075
1076	write(,) "call thermosphere ?"
1077	callthermos=.false. ! default value
1078	call getin_p("callthermos",callthermos)
1079	write(,) " callthermos = ",callthermos
1080
1081
1082	write(,) " water included without cycle ",
1083	& "(only if water=.false.)"
1084	thermoswater=.false. ! default value
1085	call getin_p("thermoswater",thermoswater)
1086	write(,) " thermoswater = ",thermoswater
1087
1088	write(,) "call thermal conduction ?",
1089	& " (only if callthermos=.true.)"
1090	callconduct=.false. ! default value
1091	call getin_p("callconduct",callconduct)
1092	write(,) " callconduct = ",callconduct
1093
1094	write(,) "call EUV heating ?",
1095	& " (only if callthermos=.true.)"
1096	calleuv=.false. ! default value
1097	call getin_p("calleuv",calleuv)
1098	write(,) " calleuv = ",calleuv
1099
1100	write(,) "call molecular viscosity ?",
1101	& " (only if callthermos=.true.)"
1102	callmolvis=.false. ! default value
1103	call getin_p("callmolvis",callmolvis)
1104	write(,) " callmolvis = ",callmolvis
1105
1106	write(,) "call molecular diffusion ?",
1107	& " (only if callthermos=.true.)"
1108	callmoldiff=.false. ! default value
1109	call getin_p("callmoldiff",callmoldiff)
1110	write(,) " callmoldiff = ",callmoldiff
1111
1112
1113	write(,) "call thermospheric photochemistry ?",
1114	& " (only if callthermos=.true.)"
1115	thermochem=.false. ! default value
1116	call getin_p("thermochem",thermochem)
1117	write(,) " thermochem = ",thermochem
1118
1119	write(,) "Method to include solar variability"
1120	write(,) "0-> fixed value of E10.7 (fixed_euv_value); ",
1121	& "1-> daily evolution of E10.7 (for given solvaryear)"
1122	solvarmod=1
1123	call getin_p("solvarmod",solvarmod)
1124	write(,) " solvarmod = ",solvarmod
1125
1126	write(,) "Fixed euv (for solvarmod==0) 10.7 value?"
1127	write(,) " (min=80 , ave=140, max=320)"
1128	fixed_euv_value=140 ! default value
1129	call getin_p("fixed_euv_value",fixed_euv_value)
1130	write(,) " fixed_euv_value = ",fixed_euv_value
1131
1132	write(,) "Solar variability as observed for MY: "
1133	write(,) "Only if solvarmod=1"
1134	solvaryear=24
1135	call getin_p("solvaryear",solvaryear)
1136	write(,) " solvaryear = ",solvaryear
1137
1138	write(,) "UV heating efficiency:",
1139	& "measured values between 0.19 and 0.23 (Fox et al. 1996)",
1140	& "lower values may be used to compensate low 15 um cooling"
1141	euveff=0.21 !default value
1142	call getin_p("euveff",euveff)
1143	write(,) " euveff = ", euveff
1144
1145
1146	if (.not.callthermos) then
1147	if (thermoswater) then
1148	print*,'if thermoswater is set, callthermos must be true'
1149	call abort_physic(modname,
1150	& "thermoswater requires callthermos",1)
1151	endif
1152	if (callconduct) then
1153	print*,'if callconduct is set, callthermos must be true'
1154	call abort_physic(modname,
1155	& "callconduct requires callthermos",1)
1156	endif
1157	if (calleuv) then
1158	print*,'if calleuv is set, callthermos must be true'
1159	call abort_physic(modname,
1160	& "calleuv requires callthermos",1)
1161	endif
1162	if (callmolvis) then
1163	print*,'if callmolvis is set, callthermos must be true'
1164	call abort_physic(modname,
1165	& "callmolvis requires callthermos",1)
1166	endif
1167	if (callmoldiff) then
1168	print*,'if callmoldiff is set, callthermos must be true'
1169	call abort_physic(modname,
1170	& "callmoldiff requires callthermos",1)
1171	endif
1172	if (thermochem) then
1173	print*,'if thermochem is set, callthermos must be true'
1174	call abort_physic(modname,
1175	& "thermochem requires callthermos",1)
1176	endif
1177	endif
1178
1179	! Test of incompatibility:
1180	! if photochem is used, then water should be used too
1181
1182	if (photochem.and..not.water) then
1183	print*,'if photochem is used, water should be used too'
1184	call abort_physic(modname,
1185	& "photochem requires water",1)
1186	endif
1187
1188	! if callthermos is used, then thermoswater should be used too
1189	! (if water not used already)
1190
1191	if (callthermos .and. .not.water) then
1192	if (callthermos .and. .not.thermoswater) then
1193	print*,'if callthermos is used, water or thermoswater
1194	& should be used too'
1195	call abort_physic(modname,
1196	& "callthermos requires water or thermoswater",1)
1197	endif
1198	endif
1199
1200	PRINT*,'--------------------------------------------'
1201	PRINT*
1202	PRINT*
1203	ELSE
1204	write(,)
1205	write(,) 'Cannot read file callphys.def. Is it here ?'
1206	call abort_physic(modname,
1207	& "missing callphys.def file",1)
1208	ENDIF
1209
1210	8000 FORMAT(t5,a12,l8)
1211	8001 FORMAT(t5,a12,i8)
1212
1213	PRINT*
1214	PRINT*,'conf_phys: daysec',daysec
1215	PRINT*
1216	PRINT*,'conf_phys: The radiative transfer is computed:'
1217	PRINT*,' each ',iradia,' physical time-step'
1218	PRINT,' or each ',iradiadtphys,' seconds'
1219	PRINT*
1220	! --------------------------------------------------------------
1221	! Managing the Longwave radiative transfer
1222	! --------------------------------------------------------------
1223
1224	! In most cases, the run just use the following values :
1225	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1226	callemis=.true.
1227	! ilwd=10*int(daysec/dtphys) ! bug before 22/10/01
1228	ilwd=1
1229	ilwn=1 !2
1230	ilwb=1 !2
1231	linear=.true.
1232	ncouche=3
1233	alphan=0.4
1234	semi=0
1235
1236	! BUT people working hard on the LW may want to read them in 'radia.def'
1237	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1238	!$OMP MASTER
1239	OPEN(99,file='radia.def',status='old',form='formatted'
1240	. ,iostat=ierr)
1241	IF(ierr.EQ.0) THEN
1242	write(,) 'conf_phys: Reading radia.def !!!'
1243	READ(99,fmt='(a)') ch1
1244	READ(99,*) callemis
1245	WRITE(*,8000) ch1,callemis
1246
1247	READ(99,fmt='(a)') ch1
1248	READ(99,*) iradia
1249	WRITE(*,8001) ch1,iradia
1250
1251	READ(99,fmt='(a)') ch1
1252	READ(99,*) ilwd
1253	WRITE(*,8001) ch1,ilwd
1254
1255	READ(99,fmt='(a)') ch1
1256	READ(99,*) ilwn
1257	WRITE(*,8001) ch1,ilwn
1258
1259	READ(99,fmt='(a)') ch1
1260	READ(99,*) linear
1261	WRITE(*,8000) ch1,linear
1262
1263	READ(99,fmt='(a)') ch1
1264	READ(99,*) ncouche
1265	WRITE(*,8001) ch1,ncouche
1266
1267	READ(99,fmt='(a)') ch1
1268	READ(99,*) alphan
1269	WRITE(,) ch1,alphan
1270
1271	READ(99,fmt='(a)') ch1
1272	READ(99,*) ilwb
1273	WRITE(*,8001) ch1,ilwb
1274
1275
1276	READ(99,fmt='(a)') ch1
1277	READ(99,'(l1)') callg2d
1278	WRITE(*,8000) ch1,callg2d
1279
1280	READ(99,fmt='(a)') ch1
1281	READ(99,*) semi
1282	WRITE(,) ch1,semi
1283	end if
1284	CLOSE(99)
1285	!$OMP END MASTER
1286	call bcast(ch1)
1287	call bcast(callemis)
1288	call bcast(iradia)
1289	call bcast(ilwd)
1290	call bcast(ilwn)
1291	call bcast(linear)
1292	call bcast(ncouche)
1293	call bcast(alphan)
1294	call bcast(ilwb)
1295	call bcast(callg2d)
1296	call bcast(semi)
1297
1298	END SUBROUTINE conf_phys
1299
1300	END MODULE conf_phys_mod

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

Download in other formats:

Original Format