Context Navigation

phyetat0_mod.F90 @ 2392

Last change on this file since 2392 was 2378, checked in by lrossi, 5 years ago

MARS GCM

HDO
Correction of an error in newstart for inihdo.
Other minor corrections for HDO cycle.
Transition from fractionation coefficients from Merlivat et al. 1967 to Lamb et al. 2017

File size: 19.5 KB

Rev	Line
[1944]	1	module phyetat0_mod
	2
	3	implicit none
	4
	5	contains
	6
[1130]	7	subroutine phyetat0 (fichnom,tab0,Lmodif,nsoil,ngrid,nlay,nq, &
[1944]	8	day_ini,time0,tsurf,tsoil,albedo,emis,q2,qsurf,co2ice, &
	9	tauscaling,totcloudfrac,wstar,mem_Mccn_co2,mem_Nccn_co2,&
[2260]	10	mem_Mh2o_co2,watercap)
[1922]	11
[1621]	12	use tracer_mod, only: noms ! tracer names
[1130]	13	use surfdat_h, only: phisfi, albedodat, z0, z0_default,&
[2079]	14	zmea, zstd, zsig, zgam, zthe, hmons, summit, base
[1130]	15	use iostart, only: nid_start, open_startphy, close_startphy, &
	16	get_field, get_var, inquire_field, &
	17	inquire_dimension, inquire_dimension_length
[2220]	18	use nonoro_gwd_ran_mod, only: du_nonoro_gwd, dv_nonoro_gwd
[2265]	19	use compute_dtau_mod, only: dtau
	20
[2281]	21	USE ioipsl_getin_p_mod, ONLY : getin_p
	22
[1130]	23	implicit none
[1974]	24
	25	include "callkeys.h"
[1130]	26	!======================================================================
	27	! Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
	28	! Adaptation à Mars : Yann Wanherdrick
	29	! Objet: Lecture de l etat initial pour la physique
	30	! Modifs: Aug.2010 EM : use NetCDF90 to load variables (enables using
	31	! r4 or r8 restarts independently of having compiled
	32	! the GCM in r4 or r8)
	33	! June 2013 TN : Possibility to read files with a time axis
	34	! November 2013 EM : Enabeling parallel, using iostart module
[2281]	35	! March 2020 AD: Enabling initialization of physics without startfi
	36	! flag: startphy_file
[1130]	37	!======================================================================
	38	INTEGER nbsrf !Mars nbsrf a 1 au lieu de 4
	39	PARAMETER (nbsrf=1) ! nombre de sous-fractions pour une maille
	40	!======================================================================
	41	! Arguments:
	42	! ---------
	43	! inputs:
[2281]	44	! logical,intent(in) :: startphy_file ! .true. if reading start file
[1130]	45	character(),intent(in) :: fichnom ! "startfi.nc" file
	46	integer,intent(in) :: tab0
	47	integer,intent(in) :: Lmodif
	48	integer,intent(in) :: nsoil ! # of soil layers
	49	integer,intent(in) :: ngrid ! # of atmospheric columns
	50	integer,intent(in) :: nlay ! # of atmospheric layers
	51	integer,intent(in) :: nq
	52	integer :: day_ini
	53	real :: time0
	54
	55	! outputs:
	56	real,intent(out) :: tsurf(ngrid) ! surface temperature
	57	real,intent(out) :: tsoil(ngrid,nsoil) ! soil temperature
[1944]	58	real,intent(out) :: albedo(ngrid,2) ! surface albedo
[1130]	59	real,intent(out) :: emis(ngrid) ! surface emissivity
	60	real,intent(out) :: q2(ngrid,nlay+1) !
	61	real,intent(out) :: qsurf(ngrid,nq) ! tracers on surface
	62	real,intent(out) :: co2ice(ngrid) ! co2 ice cover
[1208]	63	real,intent(out) :: tauscaling(ngrid) ! dust conversion factor
[1711]	64	real,intent(out) :: totcloudfrac(ngrid) ! total cloud fraction
[1944]	65	real,intent(out) :: wstar(ngrid) ! Max vertical velocity in thermals (m/s)
[1922]	66	real,intent(out) :: mem_Mccn_co2(ngrid,nlay) ! Memory of CCN mass of H2O and dust used by CO2
	67	real,intent(out) :: mem_Nccn_co2(ngrid,nlay) ! Memory of CCN number of H2O and dust used by CO2
	68	real,intent(out) :: mem_Mh2o_co2(ngrid,nlay) ! Memory of H2O mass integred into CO2 crystal
[2260]	69	real,intent(out) :: watercap(ngrid) ! h2o_ice_cover
[1130]	70	!======================================================================
	71	! Local variables:
	72
	73	real surffield(ngrid) ! to temporarily store a surface field
	74	real xmin,xmax ! to display min and max of a field
	75	!
	76	INTEGER ig,iq,lmax
	77	INTEGER nid, nvarid
	78	INTEGER ierr, i, nsrf
	79	! integer isoil
	80	! INTEGER length
	81	! PARAMETER (length=100)
	82	CHARACTER*7 str7
	83	CHARACTER*2 str2
	84	CHARACTER*1 yes
	85	!
	86	REAL p_rad,p_omeg,p_g,p_mugaz,p_daysec
	87	INTEGER nqold
	88
	89	! flag which identifies if 'startfi.nc' file is using old names (qsurf01,...)
	90	logical :: oldtracernames=.false.
	91	integer :: count
	92	character(len=30) :: txt ! to store some text
	93
	94	! specific for time
	95	REAL,ALLOCATABLE :: time(:) ! times stored in start
	96	INTEGER timelen ! number of times stored in the file
	97	INTEGER indextime ! index of selected time
	98
	99	INTEGER :: edges(3),corner(3)
	100	LOGICAL :: found
	101
[1525]	102	REAL :: timestart ! to pick which initial state to start from
[2281]	103	REAL :: surfemis ! constant emissivity when no startfi
	104	REAL :: surfalbedo ! constant emissivity when no startfi
	105	CHARACTER(len=5) :: modname="phyetat0"
[1525]	106
[2281]	107	write(,) "phyetat0: startphy_file", startphy_file
[1130]	108
[2281]	109	if (startphy_file) then
	110	! open physics initial state file:
	111	call open_startphy(fichnom)
	112	! possibility to modify tab_cntrl in tabfi
	113	write(,)
	114	write(,) 'TABFI in phyeta0: Lmodif=',Lmodif," tab0=",tab0
	115	call tabfi (nid_start,Lmodif,tab0,day_ini,lmax,p_rad, &
	116	p_omeg,p_g,p_mugaz,p_daysec,time0)
	117	else ! "academic" initialization of planetary parameters via tabfi
	118	call tabfi (0,0,0,day_ini,lmax,p_rad, &
	119	p_omeg,p_g,p_mugaz,p_daysec,time0)
	120	endif ! of if (startphy_file)
[1130]	121
[2281]	122	if (startphy_file) then
	123	! Load surface geopotential:
	124	call get_field("phisfi",phisfi,found)
	125	if (.not.found) then
	126	call abort_physic(modname, &
	127	"phyetat0: Failed loading <phisfi>",1)
	128	endif
[1130]	129	else
[2281]	130	phisfi(:)=0.
	131	endif ! of if (startphy_file)
	132	write(,) "phyetat0: surface geopotential <phisfi> range:", &
	133	minval(phisfi), maxval(phisfi)
[1130]	134
	135
[2281]	136	if (startphy_file) then
	137	! Load bare ground albedo:
	138	call get_field("albedodat",albedodat,found)
	139	if (.not.found) then
	140	call abort_physic(modname, &
	141	"phyetat0: Failed loading <albedodat>",1)
	142	endif
	143	else ! If no startfi file, use parameter surfalbedo in def file
[2294]	144	surfalbedo=0.1
	145	call getin_p("surfalbedo_without_startfi",surfalbedo)
	146	print*,"surfalbedo_without_startfi",surfalbedo
[2281]	147	albedodat(:)=surfalbedo
	148	endif ! of if (startphy_file)
	149	write(,) "phyetat0: Bare ground albedo <albedodat> range:", &
[1130]	150	minval(albedodat), maxval(albedodat)
	151
	152	! ZMEA
[2281]	153	if (startphy_file) then
	154	call get_field("ZMEA",zmea,found)
	155	if (.not.found) then
	156	call abort_physic(modname, &
	157	"phyetat0: Failed loading <ZMEA>",1)
	158	endif
[1130]	159	else
[2281]	160	zmea(:)=0.
	161	endif ! of if (startphy_file)
	162	write(,) "phyetat0: <ZMEA> range:", &
	163	minval(zmea), maxval(zmea)
[1130]	164
	165
	166	! ZSTD
[2281]	167	if (startphy_file) then
	168	call get_field("ZSTD",zstd,found)
	169	if (.not.found) then
	170	call abort_physic(modname, &
	171	"phyetat0: Failed loading <ZSTD>",1)
	172	endif
[1130]	173	else
[2281]	174	zstd(:)=0.
	175	endif ! of if (startphy_file)
	176	write(,) "phyetat0: <ZSTD> range:", &
	177	minval(zstd), maxval(zstd)
[1130]	178
	179
	180	! ZSIG
[2281]	181	if (startphy_file) then
	182	call get_field("ZSIG",zsig,found)
	183	if (.not.found) then
	184	call abort_physic(modname, &
	185	"phyetat0: Failed loading <ZSIG>",1)
	186	endif
[1130]	187	else
[2281]	188	zsig(:)=0.
	189	endif ! of if (startphy_file)
	190	write(,) "phyetat0: <ZSIG> range:", &
	191	minval(zsig), maxval(zsig)
[1130]	192
	193
	194	! ZGAM
[2281]	195	if (startphy_file) then
	196	call get_field("ZGAM",zgam,found)
	197	if (.not.found) then
	198	call abort_physic(modname, &
	199	"phyetat0: Failed loading <ZGAM>",1)
	200	endif
[1130]	201	else
[2281]	202	zgam(:)=0.
	203	endif ! of if (startphy_file)
	204	write(,) "phyetat0: <ZGAM> range:", &
	205	minval(zgam), maxval(zgam)
[1130]	206
	207
	208	! ZTHE
[2281]	209	if (startphy_file) then
	210	call get_field("ZTHE",zthe,found)
	211	if (.not.found) then
	212	call abort_physic(modname, &
	213	"phyetat0: Failed loading <ZTHE>",1)
	214	endif
[1130]	215	else
[2281]	216	zthe(:)=0.
	217	endif ! of if (startphy_file)
	218	write(,) "phyetat0: <ZTHE> range:", &
[1130]	219	minval(zthe), maxval(zthe)
	220
[2281]	221	! hmons
	222	if (startphy_file) then
	223	call get_field("hmons",hmons,found)
	224	if (.not.found) then
	225	write(,) "WARNING: phyetat0: Failed loading <hmons>"
	226	if (rdstorm) then
	227	call abort_physic(modname, &
	228	"phyetat0: Failed loading <hmons>",1)
	229	else
	230	write(,) "will continue anyway..."
	231	write(,) "because you may not need it."
	232	hmons(:)=0.
	233	end if ! if (rdstorm)
	234	else
	235	do ig=1,ngrid
	236	if (hmons(ig) .eq. -999999.) hmons(ig)=0.
	237	enddo
	238	endif ! (.not.found)
[1974]	239	else
[2281]	240	hmons(:)=0.
	241	endif ! if (startphy_file)
	242	write(,) "phyetat0: <hmons> range:", &
	243	minval(hmons), maxval(hmons)
[2079]	244
[2281]	245
	246	! summit
	247	if (startphy_file) then
	248	call get_field("summit",summit,found)
	249	if (.not.found) then
	250	write(,) "WARNING: phyetat0: Failed loading <summit>"
	251	if (rdstorm) then
	252	call abort_physic(modname, &
	253	"phyetat0: Failed loading <summit>",1)
	254	else
	255	write(,) "will continue anyway..."
	256	write(,) "because you may not need it."
	257	summit(:)=0.
	258	end if
	259	else
	260	do ig=1,ngrid
	261	if (summit(ig) .eq. -999999.) summit(ig)=0.
	262	enddo
	263	endif ! if (.not.found)
[2079]	264	else
[2281]	265	summit(:)=0.
	266	endif ! if (startphy_file)
	267	write(,) "phyetat0: <summit> range:", &
	268	minval(summit), maxval(summit)
[2079]	269
[2281]	270
	271	! base
	272	if (startphy_file) then
	273	call get_field("base",base,found)
	274	if (.not.found) then
	275	write(,) "WARNING: phyetat0: Failed loading <base>"
	276	if (rdstorm) then
	277	call abort_physic(modname, &
	278	"phyetat0: Failed loading <base>",1)
	279	else
	280	write(,) "will continue anyway..."
	281	write(,) "because you may not need it."
	282	base(:)=0.
	283	end if
	284	else
	285	do ig=1,ngrid
	286	if (base(ig) .eq. -999999.) base(ig)=0.
	287	enddo
	288	endif ! if(.not.found)
[2079]	289	else
[2281]	290	base(:)=0.
	291	endif ! if (startphy_file)
	292	write(,) "phyetat0: <base> range:", &
	293	minval(base), maxval(base)
	294
[1130]	295	! Time axis
[1525]	296	! obtain timestart from run.def
	297	timestart=-9999 ! default value
[2281]	298	call getin_p("timestart",timestart)
	299	if (startphy_file) then
	300	found=inquire_dimension("Time")
	301	if (.not.found) then
	302	indextime = 1
	303	write(,) "phyetat0: No time axis found in "//trim(fichnom)
	304	else
	305	write(,) "phyetat0: Time axis found in "//trim(fichnom)
	306	timelen=inquire_dimension_length("Time")
	307	allocate(time(timelen))
	308	! load "Time" array:
	309	call get_var("Time",time,found)
	310	if (.not.found) then
	311	call abort_physic(modname, &
	312	"phyetat0: Failed loading <Time>",1)
	313	endif
	314	! seclect the desired time index
	315	IF (timestart .lt. 0) THEN ! default: we use the last time value
	316	indextime = timelen
	317	ELSE ! else we look for the desired value in the time axis
	318	indextime = 0
	319	DO i=1,timelen
	320	IF (abs(time(i) - timestart) .lt. 0.01) THEN
	321	indextime = i
	322	EXIT
	323	ENDIF
	324	ENDDO
	325	IF (indextime .eq. 0) THEN
	326	PRINT*, "Time", timestart," is not in "//trim(fichnom)//"!!"
	327	PRINT*, "Stored times are:"
	328	DO i=1,timelen
	329	PRINT*, time(i)
	330	ENDDO
	331	call abort_physic(modname,"phyetat0: Time error",1)
	332	ENDIF
	333	ENDIF ! of IF (timestart .lt. 0)
	334	! In startfi the absolute date is day_ini + time0 + time
	335	! For now on, in the GCM physics, it is day_ini + time0
	336	time0 = time(indextime) + time0
	337	day_ini = day_ini + INT(time0)
	338	time0 = time0 - INT(time0)
	339	PRINT*, "phyetat0: Selected time ",time(indextime), &
	340	" at index ",indextime
	341	DEALLOCATE(time)
	342	endif ! of if Time not found in file
	343	else
[1130]	344	indextime = 1
[2281]	345	endif ! if (startphy_file)
[1130]	346
	347	! CO2 ice cover
[2281]	348	if (startphy_file) then
	349	call get_field("co2ice",co2ice,found,indextime)
	350	if (.not.found) then
	351	call abort_physic(modname, &
	352	"phyetat0: Failed loading <co2ice>",1)
	353	endif
[1130]	354	else
[2281]	355	co2ice(:)=0.
	356	endif !if (startphy_file)
	357	write(,) "phyetat0: CO2 ice cover <co2ice> range:", &
	358	minval(co2ice), maxval(co2ice)
[1130]	359
[1922]	360	! Memory of the origin of the co2 particles
[2281]	361	if (startphy_file) then
	362	call get_field("mem_Mccn_co2",mem_Mccn_co2,found,indextime)
	363	if (.not.found) then
	364	write(,) "phyetat0: <mem_Mccn_co2> not in file"
	365	mem_Mccn_co2(:,:)=0
	366	endif
[1922]	367	else
[2281]	368	mem_Mccn_co2(:,:)=0
	369	endif !if (startphy_file)
	370	write(,) "phyetat0: Memory of CCN mass of H2O and dust used by CO2"
	371	write(,) " <mem_Mccn_co2> range:", &
[1922]	372	minval(mem_Mccn_co2), maxval(mem_Mccn_co2)
[1130]	373
[2281]	374	if (startphy_file) then
	375	call get_field("mem_Nccn_co2",mem_Nccn_co2,found,indextime)
	376	if (.not.found) then
	377	write(,) "phyetat0: <mem_Nccn_co2> not in file"
	378	mem_Nccn_co2(:,:)=0
	379	endif
	380	else
[1922]	381	mem_Nccn_co2(:,:)=0
[2281]	382	endif ! if (startphy_file)
	383	write(,) "phyetat0: Memory of CCN number of H2O and dust used by CO2"
	384	write(,) " <mem_Nccn_co2> range:", &
[1922]	385	minval(mem_Nccn_co2), maxval(mem_Nccn_co2)
	386
[2281]	387	if (startphy_file) then
	388	call get_field("mem_Mh2o_co2",mem_Mh2o_co2,found,indextime)
	389	if (.not.found) then
	390	write(,) "phyetat0: <mem_Mh2o_co2> not in file"
	391	mem_Mh2o_co2(:,:)=0
	392	endif
[1922]	393	else
[2281]	394	mem_Mh2o_co2(:,:)=0
	395	endif ! if (startphy_file)
	396	write(,) "phyetat0: Memory of H2O mass integred into CO2 crystal"
	397	write(,) " <mem_Mh2o_co2> range:", &
[1922]	398	minval(mem_Mh2o_co2), maxval(mem_Mh2o_co2)
	399
[1208]	400	! Dust conversion factor
[2281]	401	if (startphy_file) then
	402	call get_field("tauscaling",tauscaling,found,indextime)
	403	if (.not.found) then
	404	write(,) "phyetat0: <tauscaling> not in file"
	405	tauscaling(:) = 1
	406	endif
[1208]	407	else
[2281]	408	tauscaling(:) = 1
	409	endif ! if (startphy_file)
	410	write(,) "phyetat0: dust conversion factor <tauscaling> range:", &
	411	minval(tauscaling), maxval(tauscaling)
[1208]	412
[2281]	413
[2265]	414	! dtau: opacity difference between GCM and dust scenario
[2281]	415	if (startphy_file) then
	416	call get_field("dtau",dtau,found,indextime)
	417	if (.not.found) then
	418	write(,) "phyetat0: <dtau> not in file; set to zero"
	419	dtau(:) = 0
	420	endif
[2265]	421	else
[2281]	422	dtau(:)= 0
	423	endif ! if (startphy_file)
	424	write(,) "phyetat0: opacity diff wrt scenario <dtau> range:", &
	425	minval(dtau), maxval(dtau)
[2265]	426
	427
[1711]	428	! Sub-grid cloud fraction
[2281]	429	if (startphy_file) then
	430	call get_field("totcloudfrac",totcloudfrac,found,indextime)
	431	if (.not.found) then
	432	write(,) "phyetat0: <totcloudfrac> not in file WARNING put to 1"
	433	totcloudfrac(:) = 1.0 !valeur par defaut (CLFfixval par defaut)
	434	endif
[1711]	435	else
[2281]	436	totcloudfrac(:)=1.0
	437	endif ! if (startphy_file)
	438	write(,) "phyetat0: total cloud fraction <totcloudfrac> range:", &
	439	minval(totcloudfrac), maxval(totcloudfrac)
[1208]	440
[2281]	441
[1944]	442	! Max vertical velocity in thermals
[2281]	443	if (startphy_file) then
	444	call get_field("wstar",wstar,found,indextime)
	445	if (.not.found) then
	446	write(,) "phyetat0: <wstar> not in file! Set to zero"
	447	wstar(:)=0
	448	endif
[1944]	449	else
[2281]	450	wstar(:)=0
	451	endif ! if (startphy_file)
	452	write(,) "phyetat0: Max vertical velocity in thermals <wstar> range:", &
	453	minval(wstar),maxval(wstar)
[1944]	454
[2281]	455
[1130]	456	! Surface temperature :
[2281]	457	if (startphy_file) then !tsurf
	458	call get_field("tsurf",tsurf,found,indextime)
	459	if (.not.found) then
	460	call abort_physic(modname, &
	461	"phyetat0: Failed loading <tsurf>",1)
	462	endif
[1130]	463	else
[2281]	464	tsurf(:)=0. ! will be updated afterwards in physiq !
	465	endif ! of if (startphy_file)
	466	write(,) "phyetat0: Surface temperature <tsurf> range:", &
	467	minval(tsurf), maxval(tsurf)
[1130]	468
[1944]	469	! Surface albedo
[2281]	470	if (startphy_file) then
	471	call get_field("albedo",albedo(:,1),found,indextime)
	472	if (.not.found) then
	473	write(,) "phyetat0: Failed loading <albedo>"
	474	albedo(:,1)=albedodat(:)
	475	endif
[1944]	476	else
[2281]	477	albedo(:,1)=albedodat(:)
	478	endif ! of if (startphy_file)
	479	write(,) "phyetat0: Surface albedo <albedo> range:", &
	480	minval(albedo(:,1)), maxval(albedo(:,1))
[1944]	481	albedo(:,2)=albedo(:,1)
	482
[1130]	483	! Surface emissivity
[2281]	484	if (startphy_file) then
	485	call get_field("emis",emis,found,indextime)
	486	if (.not.found) then
	487	call abort_physic(modname, &
	488	"phyetat0: Failed loading <emis>",1)
	489	endif
[1130]	490	else
[2281]	491	! If no startfi file, use parameter surfemis in def file
[2294]	492	surfemis=0.95
	493	call getin_p("surfemis_without_startfi",surfemis)
	494	print*,"surfemis_without_startfi",surfemis
[2281]	495	emis(:)=surfemis
	496	endif ! of if (startphy_file)
	497	write(,) "phyetat0: Surface emissivity <emis> range:", &
	498	minval(emis), maxval(emis)
[1130]	499
	500
	501	! surface roughness length (NB: z0 is a common in surfdat_h)
[2281]	502	if (startphy_file) then
	503	call get_field("z0",z0,found)
	504	if (.not.found) then
	505	write(,) "phyetat0: Failed loading <z0>"
	506	write(,) 'will use constant value of z0_default:',z0_default
	507	z0(:)=z0_default
	508	endif
[1130]	509	else
[2281]	510	z0(:)=z0_default
	511	endif ! of if (startphy_file)
	512	write(,) "phyetat0: Surface roughness <z0> range:", &
	513	minval(z0), maxval(z0)
[1130]	514
	515
	516	! pbl wind variance
[2281]	517	if (startphy_file) then
	518	call get_field("q2",q2,found,indextime)
	519	if (.not.found) then
	520	call abort_physic(modname, &
	521	"phyetat0: Failed loading <q2>",1)
	522	endif
[1130]	523	else
[2281]	524	q2(:,:)=0.
	525	endif ! of if (startphy_file)
	526	write(,) "phyetat0: PBL wind variance <q2> range:", &
	527	minval(q2), maxval(q2)
[1130]	528
[2281]	529	! Non-orographic gravity waves
	530	if (startphy_file) then
	531	call get_field("du_nonoro_gwd",du_nonoro_gwd,found,indextime)
	532	if (.not.found) then
	533	write(,) "phyetat0: <du_nonoro_gwd> not in file"
	534	du_nonoro_gwd(:,:)=0.
	535	endif
[2220]	536	else
[2281]	537	du_nonoro_gwd(:,:)=0.
	538	endif ! if (startphy_file)
	539	write(,) "phyetat0: Memory of zonal wind tendency due to non-orographic GW"
	540	write(,) " <du_nonoro_gwd> range:", &
	541	minval(du_nonoro_gwd), maxval(du_nonoro_gwd)
[1130]	542
[2281]	543	if (startphy_file) then
	544	call get_field("dv_nonoro_gwd",dv_nonoro_gwd,found,indextime)
	545	if (.not.found) then
	546	write(,) "phyetat0: <dv_nonoro_gwd> not in file"
	547	dv_nonoro_gwd(:,:)=0.
	548	endif
	549	else ! ! if (startphy_file)
	550	dv_nonoro_gwd(:,:)=0.
	551	endif ! if (startphy_file)
	552	write(,) "phyetat0: Memory of meridional wind tendency due to non-orographic GW"
	553	write(,) " <dv_nonoro_gwd> range:", &
	554	minval(dv_nonoro_gwd), maxval(dv_nonoro_gwd)
[2220]	555
[1130]	556	! tracer on surface
	557	if (nq.ge.1) then
	558	do iq=1,nq
[1621]	559	txt=noms(iq)
[1130]	560	if (txt.eq."h2o_vap") then
	561	! There is no surface tracer for h2o_vap;
	562	! "h2o_ice" should be loaded instead
	563	txt="h2o_ice"
	564	write(,) 'phyetat0: loading surface tracer', &
	565	' h2o_ice instead of h2o_vap'
[2260]	566	write(,) 'iq = ', iq
[1130]	567	endif
[2312]	568
	569	if (hdo) then
	570	if (txt.eq."hdo_vap") then
	571	txt="hdo_ice"
	572	write(,) 'phyetat0: loading surface tracer', &
	573	' hdo_ice instead of hdo_vap'
	574	endif
	575	endif !hdo
	576
[2281]	577	if (startphy_file) then
	578	call get_field(txt,qsurf(:,iq),found,indextime)
	579	if (.not.found) then
	580	write(,) "phyetat0: Failed loading <",trim(txt),">"
	581	write(,) " ",trim(txt)," is set to zero"
	582	endif
[1130]	583	else
[2281]	584	qsurf(:,iq)=0.
	585	endif ! of if (startphy_file)
	586	write(,) "phyetat0: Surface tracer <",trim(txt),"> range:", &
[1130]	587	minval(qsurf(:,iq)), maxval(qsurf(:,iq))
[2281]	588	enddo ! of do iq=1,nq
[1130]	589	endif ! of if (nq.ge.1)
	590
[2281]	591	if (startphy_file) then
	592	call get_field("watercap",watercap,found,indextime)
	593	if (.not.found) then
	594	write(,) "phyetat0: Failed loading <watercap> : ", &
	595	"<watercap> is set to zero"
	596	watercap(:)=0
[2260]	597
[2281]	598	write(,) 'Now transfer negative surface water ice to', &
	599	' watercap !'
	600	if (nq.ge.1) then
	601	do iq=1,nq
	602	txt=noms(iq)
	603	if (txt.eq."h2o_ice") then
	604	do ig=1,ngrid
	605	if (qsurf(ig,iq).lt.0.0) then
	606	watercap(ig) = qsurf(ig,iq)
	607	qsurf(ig,iq) = 0.0
	608	end if
	609	end do
	610	endif
[2378]	611
	612	if (txt.eq."hdo_ice") then
	613	do ig=1,ngrid
	614	if (qsurf(ig,iq).lt.0.0) then
	615	qsurf(ig,iq) = 0.0
	616	end if
	617	end do
	618	endif
	619
[2281]	620	enddo
	621	endif ! of if (nq.ge.1)
	622	endif ! of if (.not.found)
[2264]	623	else
[2281]	624	watercap(:)=0
	625	endif ! of if (startphy_file)
	626	write(,) "phyetat0: Surface water ice <watercap> range:", &
	627	minval(watercap), maxval(watercap)
	628
[2260]	629
	630
[2281]	631	if (startphy_file) then
	632	! Call to soil_settings, in order to read soil temperatures,
	633	! as well as thermal inertia and volumetric heat capacity
	634	call soil_settings(nid_start,ngrid,nsoil,tsurf,tsoil,indextime)
	635	endif ! of if (startphy_file)
[1130]	636	!
	637	! close file:
	638	!
[2281]	639	if (startphy_file) call close_startphy
[1130]	640
[1944]	641	end subroutine phyetat0
	642
	643	end module phyetat0_mod

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

Download in other formats:

Original Format