Context Navigation

meso_inifis.F @ 51

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

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

Download in other formats:

Original Format