Context Navigation

inifis.F @ 1217

Last change on this file since 1217 was 1216, checked in by emillour, 11 years ago

Generic model:
Major cleanup, in order to ease the use of LMDZ.GENERIC with (parallel) dynamics
in LMDZ.COMMON: (NB: this will break LMDZ.UNIVERSAL, which should be thrashed
in the near future)

Updated makegcm_* scripts (and makdim) and added the "-full" (to enforce full recomputation of the model) option
In dyn3d: converted control.h to module control_mod.F90 and converted iniadvtrac.F to module infotrac.F90
Added module mod_const_mpi.F90 in dyn3d (not used in serial mode)
Rearanged input/outputs routines everywhere to handle serial/MPI cases. physdem.F => phyredem.F90 , phyetat0.F => phyetat0.F90 ; all read/write routines for startfi files are gathered in module iostart.F90
added parallelism related routines init_phys_lmdz.F90, comgeomphy.F90, dimphy.F90, iniphysiq.F90, mod_grid_phy_lmdz.F90, mod_phys_lmdz_mpi_data.F90, mod_phys_lmdz_mpi_transfert.F90, mod_phys_lmdz_omp_data.F90, mod_phys_lmdz_omp_transfert.F90, mod_phys_lmdz_para.F90, mod_phys_lmdz_transfert_para.F90 in phymars and mod_const_mpi.F90 in dyn3d (for compliance with parallelism)
added created generic routines 'planetwide_maxval' and 'planetwide_minval', in module "planetwide_mod", that enable obtaining the max and min of a field over the whole planet. This should be further imroved with computation of means (possibly area weighed), etc.

File size: 23.3 KB

Line
1	SUBROUTINE inifis(ngrid,nlayer,nq,
2	$ day_ini,pdaysec,ptimestep,
3	$ plat,plon,parea,
4	$ prad,pg,pr,pcpp)
5
6	use radinc_h, only : naerkind
7	use datafile_mod, only: datadir
8	use comdiurn_h, only: sinlat, coslat, sinlon, coslon
9	use comgeomfi_h, only: long, lati, area, totarea
10	use comsoil_h, only: ini_comsoil_h
11	use control_mod, only: ecritphy
12
13	!=======================================================================
14	!
15	! purpose:
16	! -------
17	!
18	! Initialisation for the physical parametrisations of the LMD
19	! Generic Model.
20	!
21	! author: Frederic Hourdin 15 / 10 /93
22	! -------
23	! modified: Sebastien Lebonnois 11/06/2003 (new callphys.def)
24	! Ehouarn Millour (oct. 2008) tracers are now identified
25	! by their names and may not be contiguously
26	! stored in the q(:,:,:,:) array
27	! E.M. (june 2009) use getin routine to load parameters
28	!
29	!
30	! arguments:
31	! ----------
32	!
33	! input:
34	! ------
35	!
36	! ngrid Size of the horizontal grid.
37	! All internal loops are performed on that grid.
38	! nlayer Number of vertical layers.
39	! pdayref Day of reference for the simulation
40	! pday Number of days counted from the North. Spring
41	! equinoxe.
42	!
43	!=======================================================================
44	!
45	!-----------------------------------------------------------------------
46	! declarations:
47	! -------------
48	use datafile_mod, only: datadir
49	! to use 'getin'
50	USE ioipsl_getincom
51	IMPLICIT NONE
52	#include "dimensions.h"
53	#include "dimphys.h"
54	#include "planete.h"
55	#include "comcstfi.h"
56	#include "callkeys.h"
57
58
59
60	REAL,INTENT(IN) :: prad,pg,pr,pcpp,pdaysec,ptimestep
61
62	INTEGER,INTENT(IN) :: ngrid,nlayer,nq
63	REAL,INTENT(IN) :: plat(ngrid),plon(ngrid),parea(ngrid)
64	integer day_ini
65	INTEGER ig,ierr
66
67	EXTERNAL iniorbit,orbite
68	EXTERNAL SSUM
69	REAL SSUM
70
71	CHARACTER ch1*12
72	CHARACTER ch80*80
73
74	logical chem, h2o
75	logical :: parameter, doubleq=.false.
76
77	real psurf,pN2 ! added by RW for Gliese 581d N2+CO2
78
79	rad=prad
80	daysec=pdaysec
81	dtphys=ptimestep
82	cpp=pcpp
83	g=pg
84	r=pr
85	rcp=r/cpp
86
87	avocado = 6.02214179e23 ! added by RW
88
89	cursor = 1 ! added by AS in dimphys. 1 for sequential runs.
90
91	! --------------------------------------------------------
92	! The usual Tests
93	! --------------
94	IF (nlayer.NE.nlayermx) THEN
95	PRINT*,'STOP in inifis'
96	PRINT*,'Probleme de dimensions :'
97	PRINT*,'nlayer = ',nlayer
98	PRINT*,'nlayermx = ',nlayermx
99	STOP
100	ENDIF
101
102	! read in 'ecritphy' (frequency of calls to physics, in dynamical steps)
103	! (also done in dyn3d/defrun_new but not in LMDZ.COMMON)
104	call getin("ecritphy",ecritphy)
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	write(,) "Directory where external input files are:"
123	! default 'datadir' is set in "datadir_mod"
124	call getin("datadir",datadir) ! default path
125	write(,) " datadir = ",trim(datadir)
126
127	write(,) "Run with or without tracer transport ?"
128	tracer=.false. ! default value
129	call getin("tracer",tracer)
130	write(,) " tracer = ",tracer
131
132	write(,) "Run with or without atm mass update ",
133	& " due to tracer evaporation/condensation?"
134	mass_redistrib=.false. ! default value
135	call getin("mass_redistrib",mass_redistrib)
136	write(,) " mass_redistrib = ",mass_redistrib
137
138	write(,) "Diurnal cycle ?"
139	write(,) "(if diurnal=false, diurnal averaged solar heating)"
140	diurnal=.true. ! default value
141	call getin("diurnal",diurnal)
142	write(,) " diurnal = ",diurnal
143
144	write(,) "Seasonal cycle ?"
145	write(,) "(if season=false, Ls stays constant, to value ",
146	& "set in 'start'"
147	season=.true. ! default value
148	call getin("season",season)
149	write(,) " season = ",season
150
151	write(,) "Tidally resonant rotation ?"
152	tlocked=.false. ! default value
153	call getin("tlocked",tlocked)
154	write(,) "tlocked = ",tlocked
155
156	write(,) "Saturn ring shadowing ?"
157	rings_shadow = .false.
158	call getin("rings_shadow", rings_shadow)
159	write(,) "rings_shadow = ", rings_shadow
160
161	write(,) "Compute latitude-dependent gravity field?"
162	oblate = .false.
163	call getin("oblate", oblate)
164	write(,) "oblate = ", oblate
165
166	write(,) "Flattening of the planet (a-b)/a "
167	flatten = 0.0
168	call getin("flatten", flatten)
169	write(,) "flatten = ", flatten
170
171
172	write(,) "Needed if oblate=.true.: J2"
173	J2 = 0.0
174	call getin("J2", J2)
175	write(,) "J2 = ", J2
176
177	write(,) "Needed if oblate=.true.: Planet mass (*1e24 kg)"
178	MassPlanet = 0.0
179	call getin("MassPlanet", MassPlanet)
180	write(,) "MassPlanet = ", MassPlanet
181
182	write(,) "Needed if oblate=.true.: Planet mean radius (m)"
183	Rmean = 0.0
184	call getin("Rmean", Rmean)
185	write(,) "Rmean = ", Rmean
186
187	! Test of incompatibility:
188	! if tlocked, then diurnal should be false
189	if (tlocked.and.diurnal) then
190	print*,'If diurnal=true, we should turn off tlocked.'
191	stop
192	endif
193
194	write(,) "Tidal resonance ratio ?"
195	nres=0 ! default value
196	call getin("nres",nres)
197	write(,) "nres = ",nres
198
199	write(,) "Write some extra output to the screen ?"
200	lwrite=.false. ! default value
201	call getin("lwrite",lwrite)
202	write(,) " lwrite = ",lwrite
203
204	write(,) "Save statistics in file stats.nc ?"
205	callstats=.true. ! default value
206	call getin("callstats",callstats)
207	write(,) " callstats = ",callstats
208
209	write(,) "Test energy conservation of model physics ?"
210	enertest=.false. ! default value
211	call getin("enertest",enertest)
212	write(,) " enertest = ",enertest
213
214	write(,) "Check to see if cpp values used match gases.def ?"
215	check_cpp_match=.true. ! default value
216	call getin("check_cpp_match",check_cpp_match)
217	write(,) " check_cpp_match = ",check_cpp_match
218
219	write(,) "call radiative transfer ?"
220	callrad=.true. ! default value
221	call getin("callrad",callrad)
222	write(,) " callrad = ",callrad
223
224	write(,) "call correlated-k radiative transfer ?"
225	corrk=.true. ! default value
226	call getin("corrk",corrk)
227	write(,) " corrk = ",corrk
228
229	write(,) "prohibit calculations outside corrk T grid?"
230	strictboundcorrk=.true. ! default value
231	call getin("strictboundcorrk",strictboundcorrk)
232	write(,) "strictboundcorrk = ",strictboundcorrk
233
234	write(,) "call gaseous absorption in the visible bands?",
235	& "(matters only if callrad=T)"
236	callgasvis=.false. ! default value
237	call getin("callgasvis",callgasvis)
238	write(,) " callgasvis = ",callgasvis
239
240	write(,) "call continuum opacities in radiative transfer ?",
241	& "(matters only if callrad=T)"
242	continuum=.true. ! default value
243	call getin("continuum",continuum)
244	write(,) " continuum = ",continuum
245
246	write(,) "use analytic function for H2O continuum ?"
247	H2Ocont_simple=.false. ! default value
248	call getin("H2Ocont_simple",H2Ocont_simple)
249	write(,) " H2Ocont_simple = ",H2Ocont_simple
250
251	write(,) "call turbulent vertical diffusion ?"
252	calldifv=.true. ! default value
253	call getin("calldifv",calldifv)
254	write(,) " calldifv = ",calldifv
255
256	write(,) "use turbdiff instead of vdifc ?"
257	UseTurbDiff=.true. ! default value
258	call getin("UseTurbDiff",UseTurbDiff)
259	write(,) " UseTurbDiff = ",UseTurbDiff
260
261	write(,) "call convective adjustment ?"
262	calladj=.true. ! default value
263	call getin("calladj",calladj)
264	write(,) " calladj = ",calladj
265
266	write(,) "call CO2 condensation ?"
267	co2cond=.false. ! default value
268	call getin("co2cond",co2cond)
269	write(,) " co2cond = ",co2cond
270	! Test of incompatibility
271	if (co2cond.and.(.not.tracer)) then
272	print*,'We need a CO2 ice tracer to condense CO2'
273	call abort
274	endif
275
276	write(,) "CO2 supersaturation level ?"
277	co2supsat=1.0 ! default value
278	call getin("co2supsat",co2supsat)
279	write(,) " co2supsat = ",co2supsat
280
281	write(,) "Radiative timescale for Newtonian cooling ?"
282	tau_relax=30. ! default value
283	call getin("tau_relax",tau_relax)
284	write(,) " tau_relax = ",tau_relax
285	tau_relax=tau_relax243600 ! convert Earth days --> seconds
286
287	write(,)"call thermal conduction in the soil ?"
288	callsoil=.true. ! default value
289	call getin("callsoil",callsoil)
290	write(,) " callsoil = ",callsoil
291
292	write(,)"Rad transfer is computed every iradia",
293	& " physical timestep"
294	iradia=1 ! default value
295	call getin("iradia",iradia)
296	write(,)" iradia = ",iradia
297
298	write(,)"Rayleigh scattering ?"
299	rayleigh=.false.
300	call getin("rayleigh",rayleigh)
301	write(,)" rayleigh = ",rayleigh
302
303	write(,) "Use blackbody for stellar spectrum ?"
304	stelbbody=.false. ! default value
305	call getin("stelbbody",stelbbody)
306	write(,) " stelbbody = ",stelbbody
307
308	write(,) "Stellar blackbody temperature ?"
309	stelTbb=5800.0 ! default value
310	call getin("stelTbb",stelTbb)
311	write(,) " stelTbb = ",stelTbb
312
313	write(,)"Output mean OLR in 1D?"
314	meanOLR=.false.
315	call getin("meanOLR",meanOLR)
316	write(,)" meanOLR = ",meanOLR
317
318	write(,)"Output spectral OLR in 3D?"
319	specOLR=.false.
320	call getin("specOLR",specOLR)
321	write(,)" specOLR = ",specOLR
322
323	write(,)"Operate in kastprof mode?"
324	kastprof=.false.
325	call getin("kastprof",kastprof)
326	write(,)" kastprof = ",kastprof
327
328	write(,)"Uniform absorption in radiative transfer?"
329	graybody=.false.
330	call getin("graybody",graybody)
331	write(,)" graybody = ",graybody
332
333	! Test of incompatibility:
334	! if kastprof used, we must be in 1D
335	if (kastprof.and.(ngrid.gt.1)) then
336	print*,'kastprof can only be used in 1D!'
337	call abort
338	endif
339
340	write(,)"Stratospheric temperature for kastprof mode?"
341	Tstrat=167.0
342	call getin("Tstrat",Tstrat)
343	write(,)" Tstrat = ",Tstrat
344
345	write(,)"Remove lower boundary?"
346	nosurf=.false.
347	call getin("nosurf",nosurf)
348	write(,)" nosurf = ",nosurf
349
350	! Tests of incompatibility:
351	if (nosurf.and.callsoil) then
352	print*,'nosurf not compatible with soil scheme!'
353	print*,'... got to make a choice!'
354	call abort
355	endif
356
357	write(,)"Add an internal heat flux?",
358	. "... matters only if callsoil=F"
359	intheat=0.
360	call getin("intheat",intheat)
361	write(,)" intheat = ",intheat
362
363	write(,)"Use Newtonian cooling for radiative transfer?"
364	newtonian=.false.
365	call getin("newtonian",newtonian)
366	write(,)" newtonian = ",newtonian
367
368	! Tests of incompatibility:
369	if (newtonian.and.corrk) then
370	print*,'newtonian not compatible with correlated-k!'
371	call abort
372	endif
373	if (newtonian.and.calladj) then
374	print*,'newtonian not compatible with adjustment!'
375	call abort
376	endif
377	if (newtonian.and.calldifv) then
378	print*,'newtonian not compatible with a boundary layer!'
379	call abort
380	endif
381
382	write(,)"Test physics timescale in 1D?"
383	testradtimes=.false.
384	call getin("testradtimes",testradtimes)
385	write(,)" testradtimes = ",testradtimes
386
387	! Test of incompatibility:
388	! if testradtimes used, we must be in 1D
389	if (testradtimes.and.(ngrid.gt.1)) then
390	print*,'testradtimes can only be used in 1D!'
391	call abort
392	endif
393
394	write(,)"Default planetary temperature?"
395	tplanet=215.0
396	call getin("tplanet",tplanet)
397	write(,)" tplanet = ",tplanet
398
399	write(,)"Which star?"
400	startype=1 ! default value = Sol
401	call getin("startype",startype)
402	write(,)" startype = ",startype
403
404	write(,)"Value of stellar flux at 1 AU?"
405	Fat1AU=1356.0 ! default value = Sol today
406	call getin("Fat1AU",Fat1AU)
407	write(,)" Fat1AU = ",Fat1AU
408
409
410	! TRACERS:
411
412	write(,)"Varying H2O cloud fraction?"
413	CLFvarying=.false. ! default value
414	call getin("CLFvarying",CLFvarying)
415	write(,)" CLFvarying = ",CLFvarying
416
417	write(,)"Value of fixed H2O cloud fraction?"
418	CLFfixval=1.0 ! default value
419	call getin("CLFfixval",CLFfixval)
420	write(,)" CLFfixval = ",CLFfixval
421
422	write(,)"fixed radii for Cloud particles?"
423	radfixed=.false. ! default value
424	call getin("radfixed",radfixed)
425	write(,)" radfixed = ",radfixed
426
427	if(kastprof)then
428	radfixed=.true.
429	endif
430
431	write(,)"Number mixing ratio of CO2 ice particles:"
432	Nmix_co2=1.e6 ! default value
433	call getin("Nmix_co2",Nmix_co2)
434	write(,)" Nmix_co2 = ",Nmix_co2
435
436	! write(,)"Number of radiatively active aerosols:"
437	! naerkind=0. ! default value
438	! call getin("naerkind",naerkind)
439	! write(,)" naerkind = ",naerkind
440
441	write(,)"Opacity of dust (if used):"
442	dusttau=0. ! default value
443	call getin("dusttau",dusttau)
444	write(,)" dusttau = ",dusttau
445
446	write(,)"Radiatively active CO2 aerosols?"
447	aeroco2=.false. ! default value
448	call getin("aeroco2",aeroco2)
449	write(,)" aeroco2 = ",aeroco2
450
451	write(,)"Fixed CO2 aerosol distribution?"
452	aerofixco2=.false. ! default value
453	call getin("aerofixco2",aerofixco2)
454	write(,)" aerofixco2 = ",aerofixco2
455
456	write(,)"Radiatively active water ice?"
457	aeroh2o=.false. ! default value
458	call getin("aeroh2o",aeroh2o)
459	write(,)" aeroh2o = ",aeroh2o
460
461	write(,)"Fixed H2O aerosol distribution?"
462	aerofixh2o=.false. ! default value
463	call getin("aerofixh2o",aerofixh2o)
464	write(,)" aerofixh2o = ",aerofixh2o
465
466	write(,)"Radiatively active sulfuric acid aersols?"
467	aeroh2so4=.false. ! default value
468	call getin("aeroh2so4",aeroh2so4)
469	write(,)" aeroh2so4 = ",aeroh2so4
470
471	!=================================
472
473	write(,)"Radiatively active two-layer aersols?"
474	aeroback2lay=.false. ! default value
475	call getin("aeroback2lay",aeroback2lay)
476	write(,)" aeroback2lay = ",aeroback2lay
477
478	write(,)"TWOLAY AEROSOL: total optical depth ",
479	& "in the tropospheric layer (visible)"
480	obs_tau_col_tropo=8.D0
481	call getin("obs_tau_col_tropo",obs_tau_col_tropo)
482	write(,)" obs_tau_col_tropo = ",obs_tau_col_tropo
483
484	write(,)"TWOLAY AEROSOL: total optical depth ",
485	& "in the stratospheric layer (visible)"
486	obs_tau_col_strato=0.08D0
487	call getin("obs_tau_col_strato",obs_tau_col_strato)
488	write(,)" obs_tau_col_strato = ",obs_tau_col_strato
489
490	write(,)"TWOLAY AEROSOL: pres_bottom_tropo? in pa"
491	pres_bottom_tropo=66000.0
492	call getin("pres_bottom_tropo",pres_bottom_tropo)
493	write(,)" pres_bottom_tropo = ",pres_bottom_tropo
494
495	write(,)"TWOLAY AEROSOL: pres_top_tropo? in pa"
496	pres_top_tropo=18000.0
497	call getin("pres_top_tropo",pres_top_tropo)
498	write(,)" pres_top_tropo = ",pres_top_tropo
499
500	write(,)"TWOLAY AEROSOL: pres_bottom_strato? in pa"
501	pres_bottom_strato=2000.0
502	call getin("pres_bottom_strato",pres_bottom_strato)
503	write(,)" pres_bottom_strato = ",pres_bottom_strato
504
505	write(,)"TWOLAY AEROSOL: pres_top_strato? in pa"
506	pres_top_strato=100.0
507	call getin("pres_top_strato",pres_top_strato)
508	write(,)" pres_top_strato = ",pres_top_strato
509
510	write(,)"TWOLAY AEROSOL: particle size in the ",
511	& "tropospheric layer, in meters"
512	size_tropo=2.e-6
513	call getin("size_tropo",size_tropo)
514	write(,)" size_tropo = ",size_tropo
515
516	write(,)"TWOLAY AEROSOL: particle size in the ",
517	& "stratospheric layer, in meters"
518	size_strato=1.e-7
519	call getin("size_strato",size_strato)
520	write(,)" size_strato = ",size_strato
521
522	!=================================
523
524	write(,)"Cloud pressure level (with kastprof only):"
525	cloudlvl=0. ! default value
526	call getin("cloudlvl",cloudlvl)
527	write(,)" cloudlvl = ",cloudlvl
528
529	write(,)"Is the variable gas species radiatively active?"
530	Tstrat=167.0
531	varactive=.false.
532	call getin("varactive",varactive)
533	write(,)" varactive = ",varactive
534
535	write(,)"Is the variable gas species distribution set?"
536	varfixed=.false.
537	call getin("varfixed",varfixed)
538	write(,)" varfixed = ",varfixed
539
540	write(,)"What is the saturation % of the variable species?"
541	satval=0.8
542	call getin("satval",satval)
543	write(,)" satval = ",satval
544
545
546	! Test of incompatibility:
547	! if varactive, then varfixed should be false
548	if (varactive.and.varfixed) then
549	print*,'if varactive, varfixed must be OFF!'
550	stop
551	endif
552
553	write(,) "Gravitationnal sedimentation ?"
554	sedimentation=.false. ! default value
555	call getin("sedimentation",sedimentation)
556	write(,) " sedimentation = ",sedimentation
557
558	write(,) "Compute water cycle ?"
559	water=.false. ! default value
560	call getin("water",water)
561	write(,) " water = ",water
562
563	! Test of incompatibility:
564	! if water is true, there should be at least a tracer
565	if (water.and.(.not.tracer)) then
566	print*,'if water is ON, tracer must be ON too!'
567	stop
568	endif
569
570	write(,) "Include water condensation ?"
571	watercond=.false. ! default value
572	call getin("watercond",watercond)
573	write(,) " watercond = ",watercond
574
575	! Test of incompatibility:
576	! if watercond is used, then water should be used too
577	if (watercond.and.(.not.water)) then
578	print*,'if watercond is used, water should be used too'
579	stop
580	endif
581
582	write(,) "Include water precipitation ?"
583	waterrain=.false. ! default value
584	call getin("waterrain",waterrain)
585	write(,) " waterrain = ",waterrain
586
587	write(,) "Include surface hydrology ?"
588	hydrology=.false. ! default value
589	call getin("hydrology",hydrology)
590	write(,) " hydrology = ",hydrology
591
592	write(,) "Evolve surface water sources ?"
593	sourceevol=.false. ! default value
594	call getin("sourceevol",sourceevol)
595	write(,) " sourceevol = ",sourceevol
596
597	write(,) "Ice evolution timestep ?"
598	icetstep=100.0 ! default value
599	call getin("icetstep",icetstep)
600	write(,) " icetstep = ",icetstep
601
602	write(,) "Snow albedo ?"
603	albedosnow=0.5 ! default value
604	call getin("albedosnow",albedosnow)
605	write(,) " albedosnow = ",albedosnow
606
607	write(,) "Maximum ice thickness ?"
608	maxicethick=2.0 ! default value
609	call getin("maxicethick",maxicethick)
610	write(,) " maxicethick = ",maxicethick
611
612	write(,) "Freezing point of seawater ?"
613	Tsaldiff=-1.8 ! default value
614	call getin("Tsaldiff",Tsaldiff)
615	write(,) " Tsaldiff = ",Tsaldiff
616
617	write(,) "Does user want to force cpp and mugaz?"
618	force_cpp=.false. ! default value
619	call getin("force_cpp",force_cpp)
620	write(,) " force_cpp = ",force_cpp
621
622	if (force_cpp) then
623	mugaz = -99999.
624	PRINT *,'MEAN MOLECULAR MASS in g mol-1 ?'
625	call getin("mugaz",mugaz)
626	IF (mugaz.eq.-99999.) THEN
627	PRINT *, "mugaz must be set if force_cpp = T"
628	STOP
629	ELSE
630	write(,) "mugaz=",mugaz
631	ENDIF
632	cpp = -99999.
633	PRINT *,'SPECIFIC HEAT CAPACITY in J K-1 kg-1 ?'
634	call getin("cpp",cpp)
635	IF (cpp.eq.-99999.) THEN
636	PRINT *, "cpp must be set if force_cpp = T"
637	STOP
638	ELSE
639	write(,) "cpp=",cpp
640	ENDIF
641	! else
642	! mugaz=8.314*1000./pr
643	endif
644	call su_gases
645	call calc_cpp_mugaz
646
647	PRINT*,'--------------------------------------------'
648	PRINT*
649	PRINT*
650	ELSE
651	write(,)
652	write(,) 'Cannot read file callphys.def. Is it here ?'
653	stop
654	ENDIF
655
656	8000 FORMAT(t5,a12,l8)
657	8001 FORMAT(t5,a12,i8)
658
659	PRINT*
660	PRINT*,'inifis: daysec',daysec
661	PRINT*
662	PRINT*,'inifis: The radiative transfer is computed:'
663	PRINT*,' each ',iradia,' physical time-step'
664	PRINT,' or each ',iradiadtphys,' seconds'
665	PRINT*
666
667
668	!-----------------------------------------------------------------------
669	! Some more initialization:
670	! ------------------------
671
672	! ALLOCATE ARRAYS IN comgeomfi_h
673	IF (.not. ALLOCATED(lati)) ALLOCATE(lati(ngrid))
674	IF (.not. ALLOCATED(long)) ALLOCATE(long(ngrid))
675	IF (.not. ALLOCATED(area)) ALLOCATE(area(ngrid))
676
677	CALL SCOPY(ngrid,plon,1,long,1)
678	CALL SCOPY(ngrid,plat,1,lati,1)
679	CALL SCOPY(ngrid,parea,1,area,1)
680	totarea=SSUM(ngrid,area,1)
681
682	!! those are defined in comdiurn_h.F90
683	IF (.not.ALLOCATED(sinlat)) ALLOCATE(sinlat(ngrid))
684	IF (.not.ALLOCATED(coslat)) ALLOCATE(coslat(ngrid))
685	IF (.not.ALLOCATED(sinlon)) ALLOCATE(sinlon(ngrid))
686	IF (.not.ALLOCATED(coslon)) ALLOCATE(coslon(ngrid))
687
688	DO ig=1,ngrid
689	sinlat(ig)=sin(plat(ig))
690	coslat(ig)=cos(plat(ig))
691	sinlon(ig)=sin(plon(ig))
692	coslon(ig)=cos(plon(ig))
693	ENDDO
694
695	pi=2.*asin(1.) ! NB: pi is a common in comcstfi.h
696
697	! allocate "comsoil_h" arrays
698	call ini_comsoil_h(ngrid)
699
700	END

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

Download in other formats:

Original Format