Context Navigation

meso_inifis.F @ 190

Last change on this file since 190 was 162, checked in by aslmd, 14 years ago
MESOSCALE: changements de la revision precedente (thermiques) impactes aux routines specifiques meso_physiq.F, meso_inifis.F et meso_callkeys.h. transparent pour l'utilisateur en attendant que les thermiques soient validees. compilation OK. il faudra ajouter des liens symboliques dans MESOSCALE/LMDZ.MARS une fois valide.
File size: 25.8 KB

Rev	Line
[42]	1	SUBROUTINE meso_inifis(ngrid,nlayer,nq,
	2	$ wdt,
	3	$ wday_ini,wdaysec,
	4	$ wappel_phys,wecri_phys,
	5	$ plat,plon,parea,
	6	$ prad,pg,pr,pcpp,
	7	$ womeg,wmugaz,
	8	$ wyear_day,wperiheli,waphelie,wperi_day,wobliquit,
	9	$ wz0,wemin_turb,wlmixmin,
	10	$ wemissiv,wemissiceN,wemissiceS,walbediceN,walbediceS,
	11	$ wiceradiusN,wiceradiusS,wdtemisiceN,wdtemisiceS,
	12	$ walbedodat,wphisfi,wvolcapa,
	13	$ wzmea,wzstd,wzsig,wzgam,wzthe,
	14	$ wtheta,wpsi)
	15
	16	c=======================================================================
	17	c
	18	c CAREFUL: THIS IS A VERSION TO BE USED WITH WRF !!!
	19	c
	20	c ... CHECK THE ****WRF lines
	21	c
	22	c=======================================================================
	23	!
	24	!=======================================================================
	25	!
	26	! purpose:
	27	! -------
	28	!
	29	! Initialisation for the physical parametrisations of the LMD
	30	! martian atmospheric general circulation modele.
	31	!
	32	! author: Frederic Hourdin 15 / 10 /93
	33	! -------
	34	! modified: Sebastien Lebonnois 11/06/2003 (new callphys.def)
	35	! Ehouarn Millour (oct. 2008) tracers are now identified
	36	! by their names and may not be contiguously
	37	! stored in the q(:,:,:,:) array
	38	! E.M. (june 2009) use getin routine to load parameters
	39	! adapted to the WRF use - Aymeric Spiga - Jan 2009 - Jan 2007
	40	!
	41	!
	42	! arguments:
	43	! ----------
	44	!
	45	! input:
	46	! ------
	47	!
	48	! ngrid Size of the horizontal grid.
	49	! All internal loops are performed on that grid.
	50	! nlayer Number of vertical layers.
	51	! pdayref Day of reference for the simulation
	52	! pday Number of days counted from the North. Spring
	53	! equinoxe.
	54	!
	55	!=======================================================================
	56	!
	57	!-----------------------------------------------------------------------
	58	! declarations:
	59	! -------------
	60	! to use 'getin'
	61	!! **WRF a new stuff to be checked
	62	USE ioipsl_getincom
	63	IMPLICIT NONE
	64	#include "dimensions.h"
	65	#include "dimphys.h"
	66	#include "planete.h"
	67	#include "comcstfi.h"
	68	#include "comsaison.h"
	69	#include "comdiurn.h"
	70	#include "comgeomfi.h"
	71	#include "callkeys.h"
	72	#include "surfdat.h"
	73	#include "dimradmars.h" !!new
	74	#include "yomaer.h" !!new -- prob. pour remplir les proprietes
	75	#include "datafile.h"
	76	#include "meso_slope.h"
	77	#include "comsoil.h" !!**WRF -- needed to fill volcapa
	78
	79	REAL prad,pg,pr,pcpp,pdaysec
	80
	81	INTEGER ngrid,nlayer,nq
	82
	83	REAL womeg,wmugaz,wdaysec
	84	REAL wyear_day,wperiheli,waphelie,wperi_day,wobliquit
	85	REAL wz0,wemin_turb,wlmixmin
	86	REAL wemissiv,wemissiceN,wemissiceS,walbediceN,walbediceS
	87
	88	REAL wiceradiusN,wiceradiusS,wdtemisiceN,wdtemisiceS
	89	REAL walbedodat(ngrid),wphisfi(ngrid)
	90	REAL wzmea(ngrid),wzstd(ngrid),wzsig(ngrid)
	91	REAL wzgam(ngrid),wzthe(ngrid)
	92	REAL wtheta(ngrid),wpsi(ngrid)
	93	REAL wvolcapa
	94	REAL wdt
	95
	96	REAL plat(ngrid),plon(ngrid),parea(ngridmx)
	97	integer wday_ini
	98	INTEGER ig,ierr
	99
	100	INTEGER wecri_phys, wappel_phys
	101
	102	! EXTERNAL iniorbit,orbite
	103	EXTERNAL SSUM
	104	REAL SSUM
	105
	106	CHARACTER ch1*12
	107	CHARACTER ch80*80
	108
[86]	109	#ifdef MESOSCALE
	110
[42]	111	! logical chem, h2o
	112
	113	! chem = .false.
	114	! h2o = .false.
	115
	116	c ****WRF
	117	c
	118	c------------------------------------------------------
	119	c Fill some parameters in the 'include' files
	120	c >> Do part of the job previously done by phyetat0.F
	121	c >> Complete list of parameters is found in tabfi.F
	122	c------------------------------------------------------
	123	c
	124	c Values are defined in the module_model_constants.F WRF routine
	125	c
	126	! in 'comcstfi.h'
	127	rad=prad
	128	cpp=pcpp
	129	g=pg
	130	r=pr
	131	rcp=r/cpp
	132	daysec=wdaysec
	133	omeg=womeg
	134	mugaz=wmugaz
	135	print*,"check: rad,cpp,g,r,rcp,daysec,omeg,mugaz"
	136	print*,rad,cpp,g,r,rcp,daysec,omeg,mugaz
	137
	138	! in 'planet.h'
	139	year_day=wyear_day
	140	periheli=wperiheli
	141	aphelie=waphelie
	142	peri_day=wperi_day
	143	obliquit=wobliquit
	144	z0=wz0
	145	emin_turb=wemin_turb
	146	lmixmin=wlmixmin
	147	print*,"check: year_day,periheli,aphelie,peri_day,obliquit"
	148	print*,year_day,periheli,aphelie,peri_day,obliquit
	149	print*,"check: z0,emin_turb,lmixmin"
	150	print*,z0,emin_turb,lmixmin
	151
	152	! in 'surfdat.h'
	153	emissiv=wemissiv
	154	emisice(1)=wemissiceN
	155	emisice(2)=wemissiceS
	156	albedice(1)=walbediceN
	157	albedice(2)=walbediceS
	158	iceradius(1)=wiceradiusN
	159	iceradius(2)=wiceradiusS
	160	dtemisice(1)=wdtemisiceN
	161	dtemisice(2)=wdtemisiceS
	162	print*,"check: emissiv,emisice,albedice,iceradius,dtemisice"
	163	print*,emissiv,emisice,albedice,iceradius,dtemisice
	164
	165	c
	166	c Values are defined in the WPS processing
	167	c
	168	albedodat(:)=walbedodat(:)
	169	!!!!! ***WRF inertiedat was moved, new physics !!
	170	!inertiedat(:)=winertiedat(:)
	171	phisfi(:)=wphisfi(:)
	172	print*,"check: albedodat(1),phisfi(1)"
	173	print*,albedodat(1),phisfi(1)
	174	print*,"check: albedodat(end),phisfi(end)"
	175	print*,albedodat(ngrid),phisfi(ngrid)
	176
	177	! NB: usually, gravity wave scheme is useless in mesoscale modeling
	178	! NB: we however keep the option for coarse grid case ...
	179	zmea(:)=wzmea(:)
	180	zstd(:)=wzstd(:)
	181	zsig(:)=wzsig(:)
	182	zgam(:)=wzgam(:)
	183	zthe(:)=wzthe(:)
	184	print*,"check: gw param"
	185	print*,zmea(1),zmea(ngrid)
	186	print*,zstd(1),zstd(ngrid)
	187	print*,zsig(1),zsig(ngrid)
	188	print*,zgam(1),zgam(ngrid)
	189	print*,zthe(1),zthe(ngrid)
	190
	191	!
	192	! in meso_slope.h
	193	!
	194	theta_sl(:)=wtheta(:)
	195	psi_sl(:)=wpsi(:)
	196	print*,"check: theta_sl(1),psi_sl(1)"
	197	print*,theta_sl(1),psi_sl(1)
	198	print*,"check: theta_sl(end),psi_sl(end)"
	199	print*,theta_sl(ngrid),psi_sl(ngrid)
	200
	201	!
	202	! in comsoil.h
	203	!
	204	volcapa=wvolcapa
	205	print*,"check: volcapa"
	206	print*,volcapa
	207
	208
	209	c ****WRF
	210
	211
	212	! rad=prad
	213	! daysec=pdaysec
	214	! dtphys=ptimestep
	215	! cpp=pcpp
	216	! g=pg
	217	! r=pr
	218	! rcp=r/cpp
	219
	220	! --------------------------------------------------------
	221	! The usual Tests
	222	! --------------
	223	IF (nlayer.NE.nlayermx) THEN
	224	PRINT*,'STOP in inifis'
	225	PRINT*,'Probleme de dimensions :'
	226	PRINT*,'nlayer = ',nlayer
	227	PRINT*,'nlayermx = ',nlayermx
	228	STOP
	229	ENDIF
	230
	231	IF (ngrid.NE.ngridmx) THEN
	232	PRINT*,'STOP in inifis'
	233	PRINT*,'Probleme de dimensions :'
	234	PRINT*,'ngrid = ',ngrid
	235	PRINT*,'ngridmx = ',ngridmx
	236	STOP
	237	ENDIF
	238
	239	! --------------------------------------------------------------
	240	! Reading the "callphys.def" file controlling some key options
	241	! --------------------------------------------------------------
	242
	243	! check that 'callphys.def' file is around
	244	OPEN(99,file='callphys.def',status='old',form='formatted'
	245	& ,iostat=ierr)
	246	CLOSE(99)
	247
	248	IF(ierr.EQ.0) THEN
	249	PRINT*
	250	PRINT*
	251	PRINT*,'--------------------------------------------'
	252	PRINT*,' inifis: Parameters for the physics (callphys.def)'
	253	PRINT*,'--------------------------------------------'
	254
	255	write(,) "Directory where external input files are:"
[162]	256	!datafile="/u/forget/WWW/datagcm/datafile"
[42]	257	! NB: default 'datafile' is set in datafile.h
	258	call getin("datadir",datafile)
	259	write(,) " datafile = ",trim(datafile)
	260
	261	write(,) "Run with or without tracer transport ?"
	262	tracer=.false. ! default value
	263	call getin("tracer",tracer)
	264	write(,) " tracer = ",tracer
	265
	266	write(,) "Diurnal cycle ?"
	267	write(,) "(if diurnal=False, diurnal averaged solar heating)"
	268	diurnal=.true. ! default value
	269	call getin("diurnal",diurnal)
	270	write(,) " diurnal = ",diurnal
	271
	272	write(,) "Seasonal cycle ?"
	273	write(,) "(if season=False, Ls stays constant, to value ",
	274	& "set in 'start'"
	275	season=.true. ! default value
	276	call getin("season",season)
	277	write(,) " season = ",season
	278
	279	write(,) "Write some extra output to the screen ?"
	280	lwrite=.false. ! default value
	281	call getin("lwrite",lwrite)
	282	write(,) " lwrite = ",lwrite
	283
	284	write(,) "Save statistics in file stats.nc ?"
	285	!callstats=.true. ! default value
	286	callstats=.false. ! default value
	287	call getin("callstats",callstats)
	288	write(,) " callstats = ",callstats
	289
	290	write(,) "Save EOF profiles in file 'profiles' for ",
	291	& "Climate Database?"
	292	calleofdump=.false. ! default value
	293	call getin("calleofdump",calleofdump)
	294	write(,) " calleofdump = ",calleofdump
	295
	296	write(,) "Dust scenario: 1=constant dust (read from startfi",
	297	& " or set as tauvis); 2=Viking scenario; =3 MGS scenario,",
	298	& "=4 Mars Year 24 from TES assimilation, ",
	299	& "=24,25 or 26 :Mars Year 24,25 or 26 from TES assimilation"
	300	iaervar=3 ! default value
	301	call getin("iaervar",iaervar)
	302	write(,) " iaervar = ",iaervar
	303
	304	write(,) "Reference (visible) dust opacity at 700 Pa ",
	305	& "(matters only if iaervar=1)"
	306	! NB: default value of tauvis is set/read in startfi.nc file
	307	call getin("tauvis",tauvis)
	308	write(,) " tauvis = ",tauvis
	309
	310	write(,) "Dust vertical distribution:"
	311	write(,) "(=1 top set by topdustref parameter;",
	312	& " =2 Viking scenario; =3 MGS scenario)"
	313	iddist=3 ! default value
	314	call getin("iddist",iddist)
	315	write(,) " iddist = ",iddist
	316
	317	write(,) "Dust top altitude (km). (Matters only if iddist=1)"
	318	topdustref= 90.0 ! default value
	319	call getin("topdustref",topdustref)
	320	write(,) " topdustref = ",topdustref
	321
	322	write(,) "call radiative transfer ?"
	323	callrad=.true. ! default value
	324	call getin("callrad",callrad)
	325	write(,) " callrad = ",callrad
	326
	327	write(,) "call NLTE radiative schemes ?",
	328	& "(matters only if callrad=T)"
	329	callnlte=.false. ! default value
	330	call getin("callnlte",callnlte)
	331	write(,) " callnlte = ",callnlte
	332
	333	write(,) "call CO2 NIR absorption ?",
	334	& "(matters only if callrad=T)"
	335	callnirco2=.false. ! default value
	336	call getin("callnirco2",callnirco2)
	337	write(,) " callnirco2 = ",callnirco2
	338
	339	write(,) "call turbulent vertical diffusion ?"
	340	calldifv=.true. ! default value
	341	call getin("calldifv",calldifv)
	342	write(,) " calldifv = ",calldifv
	343
[162]	344	write(,) "call thermals ?"
	345	calltherm=.false. ! default value
	346	call getin("calltherm",calltherm)
	347	write(,) " calltherm = ",calltherm
	348
	349	write(,) "output thermal diagnostics ?"
	350	outptherm=.false. ! default value
	351	call getin("outptherm",outptherm)
	352	write(,) " outptherm = ",outptherm
	353
[42]	354	write(,) "call convective adjustment ?"
	355	calladj=.true. ! default value
	356	call getin("calladj",calladj)
	357	write(,) " calladj = ",calladj
	358
[162]	359	if (calltherm .and. (.not. calladj)) then
	360	print*,'Convadj has to be activated when using thermals'
	361	stop
	362	endif
[42]	363
	364	write(,) "call CO2 condensation ?"
	365	callcond=.true. ! default value
	366	call getin("callcond",callcond)
	367	write(,) " callcond = ",callcond
	368
	369	write(,)"call thermal conduction in the soil ?"
	370	callsoil=.true. ! default value
	371	call getin("callsoil",callsoil)
	372	write(,) " callsoil = ",callsoil
	373
	374
	375	write(,)"call Lott's gravity wave/subgrid topography ",
	376	& "scheme ?"
	377	calllott=.true. ! default value
	378	call getin("calllott",calllott)
	379	write(,)" calllott = ",calllott
	380
	381
	382	write(,)"rad.transfer is computed every iradia",
	383	& " physical timestep"
	384	iradia=1 ! default value
	385	call getin("iradia",iradia)
	386	write(,)" iradia = ",iradia
	387
	388
	389	write(,)"Output of the exchange coefficient mattrix ?",
	390	& "(for diagnostics only)"
	391	callg2d=.false. ! default value
	392	call getin("callg2d",callg2d)
	393	write(,)" callg2d = ",callg2d
	394
	395	write(,)"Rayleigh scattering : (should be .false. for now)"
	396	rayleigh=.false.
	397	call getin("rayleigh",rayleigh)
	398	write(,)" rayleigh = ",rayleigh
	399
	400
	401	! TRACERS:
	402
	403	write(,)"Transported dust ? (if >0, use 'dustbin' dust bins)"
	404	dustbin=0 ! default value
	405	call getin("dustbin",dustbin)
	406	write(,)" dustbin = ",dustbin
	407
	408	write(,)"Radiatively active dust ? (matters if dustbin>0)"
	409	active=.false. ! default value
	410	call getin("active",active)
	411	write(,)" active = ",active
	412
	413	! Test of incompatibility:
	414	! if active is used, then dustbin should be > 0
	415
	416	if (active.and.(dustbin.lt.1)) then
	417	print*,'if active is used, then dustbin should > 0'
	418	stop
	419	endif
	420
	421	write(,)"use mass and number mixing ratios to predict",
	422	& " dust size ?"
	423	doubleq=.false. ! default value
	424	call getin("doubleq",doubleq)
	425	write(,)" doubleq = ",doubleq
	426
	427	submicron=.false. ! default value
	428	call getin("submicron",submicron)
	429	write(,)" submicron = ",submicron
	430
	431	! Test of incompatibility:
	432	! if doubleq is used, then dustbin should be 2
	433
	434	if (doubleq.and.(dustbin.ne.2)) then
	435	print*,'if doubleq is used, then dustbin should be 2'
	436	stop
	437	endif
	438	if (doubleq.and.submicron.and.(nqmx.LT.3)) then
	439	print*,'If doubleq is used with a submicron tracer,'
	440	print*,' then the number of tracers has to be'
	441	print*,' larger than 3.'
	442	stop
	443	endif
	444
	445	write(,)"dust lifted by GCM surface winds ?"
	446	lifting=.false. ! default value
	447	call getin("lifting",lifting)
	448	write(,)" lifting = ",lifting
	449
	450	! Test of incompatibility:
	451	! if lifting is used, then dustbin should be > 0
	452
	453	if (lifting.and.(dustbin.lt.1)) then
	454	print*,'if lifting is used, then dustbin should > 0'
	455	stop
	456	endif
	457
	458	write(,)" dust lifted by dust devils ?"
	459	callddevil=.false. !default value
	460	call getin("callddevil",callddevil)
	461	write(,)" callddevil = ",callddevil
	462
	463
	464	! Test of incompatibility:
	465	! if dustdevil is used, then dustbin should be > 0
	466
	467	if (callddevil.and.(dustbin.lt.1)) then
	468	print*,'if dustdevil is used, then dustbin should > 0'
	469	stop
	470	endif
	471
	472	write(,)"Dust scavenging by CO2 snowfall ?"
	473	scavenging=.false. ! default value
	474	call getin("scavenging",scavenging)
	475	write(,)" scavenging = ",scavenging
	476
	477
	478	! Test of incompatibility:
	479	! if scavenging is used, then dustbin should be > 0
	480
	481	if (scavenging.and.(dustbin.lt.1)) then
	482	print*,'if scavenging is used, then dustbin should > 0'
	483	stop
	484	endif
	485
	486	write(,) "Gravitationnal sedimentation ?"
	487	sedimentation=.true. ! default value
	488	call getin("sedimentation",sedimentation)
	489	write(,) " sedimentation = ",sedimentation
	490
	491	write(,) "Radiatively active transported atmospheric ",
	492	& "water ice ?"
	493	activice=.false. ! default value
	494	call getin("activice",activice)
	495	write(,) " activice = ",activice
	496
	497	write(,) "Compute water cycle ?"
	498	water=.false. ! default value
	499	call getin("water",water)
	500	write(,) " water = ",water
	501
	502	! Test of incompatibility:
	503
	504	if (activice.and..not.water) then
	505	print*,'if activice is used, water should be used too'
	506	stop
	507	endif
	508
	509	if (water.and..not.tracer) then
	510	print*,'if water is used, tracer should be used too'
	511	stop
	512	endif
	513
	514	! Test of incompatibility:
	515
	516	write(,) "Permanent water caps at poles ?",
	517	& " .true. is RECOMMENDED"
	518	write(,) "(with .true., North cap is a source of water ",
	519	& "and South pole is a cold trap)"
	520	caps=.true. ! default value
	521	call getin("caps",caps)
	522	write(,) " caps = ",caps
	523
	524	write(,) "photochemistry: include chemical species"
	525	photochem=.false. ! default value
	526	call getin("photochem",photochem)
	527	write(,) " photochem = ",photochem
	528
	529
	530	! THERMOSPHERE
	531
	532	write(,) "call thermosphere ?"
	533	callthermos=.false. ! default value
	534	call getin("callthermos",callthermos)
	535	write(,) " callthermos = ",callthermos
	536
	537
	538	write(,) " water included without cycle ",
	539	& "(only if water=.false.)"
	540	thermoswater=.false. ! default value
	541	call getin("thermoswater",thermoswater)
	542	write(,) " thermoswater = ",thermoswater
	543
	544	write(,) "call thermal conduction ?",
	545	& " (only if callthermos=.true.)"
	546	callconduct=.false. ! default value
	547	call getin("callconduct",callconduct)
	548	write(,) " callconduct = ",callconduct
	549
	550	write(,) "call EUV heating ?",
	551	& " (only if callthermos=.true.)"
	552	calleuv=.false. ! default value
	553	call getin("calleuv",calleuv)
	554	write(,) " calleuv = ",calleuv
	555
	556	write(,) "call molecular viscosity ?",
	557	& " (only if callthermos=.true.)"
	558	callmolvis=.false. ! default value
	559	call getin("callmolvis",callmolvis)
	560	write(,) " callmolvis = ",callmolvis
	561
	562	write(,) "call molecular diffusion ?",
	563	& " (only if callthermos=.true.)"
	564	callmoldiff=.false. ! default value
	565	call getin("callmoldiff",callmoldiff)
	566	write(,) " callmoldiff = ",callmoldiff
	567
	568
	569	write(,) "call thermospheric photochemistry ?",
	570	& " (only if callthermos=.true.)"
	571	thermochem=.false. ! default value
	572	call getin("thermochem",thermochem)
	573	write(,) " thermochem = ",thermochem
	574
	575	write(,) "date for solar flux calculation:",
	576	& " (1985 < date < 2002)"
	577	write(,) "(Solar min=1996.4 ave=1993.4 max=1990.6)"
	578	solarcondate=1993.4 ! default value
	579	call getin("solarcondate",solarcondate)
	580	write(,) " solarcondate = ",solarcondate
	581
	582
	583	if (.not.callthermos) then
	584	if (thermoswater) then
	585	print*,'if thermoswater is set, callthermos must be true'
	586	stop
	587	endif
	588	if (callconduct) then
	589	print*,'if callconduct is set, callthermos must be true'
	590	stop
	591	endif
	592	if (calleuv) then
	593	print*,'if calleuv is set, callthermos must be true'
	594	stop
	595	endif
	596	if (callmolvis) then
	597	print*,'if callmolvis is set, callthermos must be true'
	598	stop
	599	endif
	600	if (callmoldiff) then
	601	print*,'if callmoldiff is set, callthermos must be true'
	602	stop
	603	endif
	604	if (thermochem) then
	605	print*,'if thermochem is set, callthermos must be true'
	606	stop
	607	endif
	608	endif
	609
	610	! Test of incompatibility:
	611	! if photochem is used, then water should be used too
	612
	613	if (photochem.and..not.water) then
	614	print*,'if photochem is used, water should be used too'
	615	stop
	616	endif
	617
	618	! if callthermos is used, then thermoswater should be used too
	619	! (if water not used already)
	620
	621	if (callthermos .and. .not.water) then
	622	if (callthermos .and. .not.thermoswater) then
	623	print*,'if callthermos is used, water or thermoswater
	624	& should be used too'
	625	stop
	626	endif
	627	endif
	628
	629	PRINT*,'--------------------------------------------'
	630	PRINT*
	631	PRINT*
	632	ELSE
	633	write(,)
	634	write(,) 'Cannot read file callphys.def. Is it here ?'
	635	stop
	636	ENDIF
	637
	638	8000 FORMAT(t5,a12,l8)
	639	8001 FORMAT(t5,a12,i8)
	640
	641
	642
	643	c ****WRF
	644	c*****************************************************
	645	c Since it comes from WRF settings, we have to
	646	c fill dtphys in the include file
	647	c It must be set now, because it is used afterwards
	648	c
	649	c Opportunity is taken to fill ecri_phys as well
	650	c*****************************************************
	651	dtphys=wdt*float(wappel_phys)
	652	print*,'Physical timestep (s) ',dtphys
	653	ecri_phys=wecri_phys
	654	print*,'Physical frequency for writing ',ecri_phys
	655	c
	656	c ****WRF
	657
	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	! Managing the Longwave radiative transfer
	668	! --------------------------------------------------------------
	669
	670	! In most cases, the run just use the following values :
	671	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	672	callemis=.true.
	673	! ilwd=10*int(daysec/dtphys) ! bug before 22/10/01
	674	ilwd=1
	675	ilwn=1 !2
	676	ilwb=1 !2
	677	linear=.true.
	678	ncouche=3
	679	alphan=0.4
	680	semi=0
	681
	682	! BUT people working hard on the LW may want to read them in 'radia.def'
	683	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
	684
	685	OPEN(99,file='radia.def',status='old',form='formatted'
	686	. ,iostat=ierr)
	687	IF(ierr.EQ.0) THEN
	688	write(,) 'inifis: Reading radia.def !!!'
	689	READ(99,fmt='(a)') ch1
	690	READ(99,*) callemis
	691	WRITE(*,8000) ch1,callemis
	692
	693	READ(99,fmt='(a)') ch1
	694	READ(99,*) iradia
	695	WRITE(*,8001) ch1,iradia
	696
	697	READ(99,fmt='(a)') ch1
	698	READ(99,*) ilwd
	699	WRITE(*,8001) ch1,ilwd
	700
	701	READ(99,fmt='(a)') ch1
	702	READ(99,*) ilwn
	703	WRITE(*,8001) ch1,ilwn
	704
	705	READ(99,fmt='(a)') ch1
	706	READ(99,*) linear
	707	WRITE(*,8000) ch1,linear
	708
	709	READ(99,fmt='(a)') ch1
	710	READ(99,*) ncouche
	711	WRITE(*,8001) ch1,ncouche
	712
	713	READ(99,fmt='(a)') ch1
	714	READ(99,*) alphan
	715	WRITE(,) ch1,alphan
	716
	717	READ(99,fmt='(a)') ch1
	718	READ(99,*) ilwb
	719	WRITE(*,8001) ch1,ilwb
	720
	721
	722	READ(99,fmt='(a)') ch1
	723	READ(99,'(l1)') callg2d
	724	WRITE(*,8000) ch1,callg2d
	725
	726	READ(99,fmt='(a)') ch1
	727	READ(99,*) semi
	728	WRITE(,) ch1,semi
	729	end if
	730	CLOSE(99)
	731
	732	!-----------------------------------------------------------------------
	733	! Some more initialization:
	734	! ------------------------
	735
	736	! in 'comgeomfi.h'
	737	CALL SCOPY(ngrid,plon,1,long,1)
	738	CALL SCOPY(ngrid,plat,1,lati,1)
	739	CALL SCOPY(ngrid,parea,1,area,1)
	740	totarea=SSUM(ngridmx,area,1)
	741
	742	! in 'comdiurn.h'
	743	DO ig=1,ngrid
	744	sinlat(ig)=sin(plat(ig))
	745	coslat(ig)=cos(plat(ig))
	746	sinlon(ig)=sin(plon(ig))
	747	coslon(ig)=cos(plon(ig))
	748	ENDDO
	749
	750	pi=2.*asin(1.) ! NB: pi is a common in comcstfi.h
	751
	752	! managing the tracers, and tests:
	753	! -------------------------------
	754	! Ehouarn: removed; as these tests are now done in initracer.F
	755	! if(tracer) then
	756	!
	757	!! when photochem is used, nqchem_min is the rank
	758	!! of the first chemical species
	759	!
	760	!! Ehouarn: nqchem_min is now meaningless and no longer used
	761	!! nqchem_min = 1
	762	! if (photochem .or. callthermos) then
	763	! chem = .true.
	764	! end if
	765	!
	766	! if (water .or. thermoswater) h2o = .true.
	767	!
	768	!! TESTS
	769	!
	770	! print*,'inifis: TRACERS:'
	771	! write(,) " chem=",chem," h2o=",h2o
	772	!! write(,) " doubleq=",doubleq
	773	!! write(,) " dustbin=",dustbin
	774	!
	775	! if ((doubleq).and.(h2o).and.
	776	! $ (chem)) then
	777	! print*,' 2 dust tracers (doubleq)'
	778	! print*,' 1 water vapour tracer'
	779	! print*,' 1 water ice tracer'
	780	! print*,nqmx-4,' chemistry tracers'
	781	! endif
	782	!
	783	! if ((doubleq).and.(h2o).and.
	784	! $ .not.(chem)) then
	785	! print*,' 2 dust tracers (doubleq)'
	786	! print*,' 1 water vapour tracer'
	787	! print*,' 1 water ice tracer'
	788	! if (nqmx.LT.4) then
	789	! print*,'nqmx should be at least equal to'
	790	! print*,'4 with these options.'
	791	! stop
	792	! endif
	793	! endif
	794	!
	795	! if (.not.(doubleq).and.(h2o).and.
	796	! $ (chem)) then
	797	! if (dustbin.gt.0) then
	798	! print*,dustbin,' dust bins'
	799	! endif
	800	! print*,nqmx-2-dustbin,' chemistry tracers'
	801	! print*,' 1 water vapour tracer'
	802	! print*,' 1 water ice tracer'
	803	! endif
	804	!
	805	! if (.not.(doubleq).and.(h2o).and.
	806	! $ .not.(chem)) then
	807	! if (dustbin.gt.0) then
	808	! print*,dustbin,' dust bins'
	809	! endif
	810	! print*,' 1 water vapour tracer'
	811	! print*,' 1 water ice tracer'
	812	! if (nqmx.gt.(dustbin+2)) then
	813	! print*,'nqmx should be ',(dustbin+2),
	814	! $ ' with these options...'
	815	! print*,'(or check callphys.def)'
	816	! endif
	817	! if (nqmx.lt.(dustbin+2)) then
	818	! write(,) "inifis: nqmx.lt.(dustbin+2)"
	819	! stop
	820	! endif
	821	! endif
	822	!
	823	! endif ! of if (tracer)
	824	!
	825	! RETURN
[86]	826
	827	#endif
	828
[42]	829	END

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

Download in other formats:

Original Format