Context Navigation

phyetat0_mod.F90 @ 2443

Last change on this file since 2443 was 2417, checked in by emillour, 5 years ago

Mars GCM:
Add a new scheme to handle the dust and its radiative impact. Triggered using
a new flag dustscaling_mode=2 (dustscaling_mod=0: no rescaling at all, and
dustscaling_mode=1: full rescaling using tauscaling, GCMv5.3 style). Rescaling
is then only done on the radiative impact (see dust_scaling_mod.F90) of dust.
Moreover the scaling factor "dust_rad_adjust" is evaluated using the target dust
scenario opacity for the next sol and left to evolve linearly until then to not
impose the diurnal evolution of dust.
In practice, main changes or additions in the code are:

renamed flag "tauscaling_mode" as "dustscaling_mode"
moved parameter "t_scenario_sol" to "dust_param_mod"
adapted "compute_dustscaling" routine in "dust_scaling_mod"
added module "dust_rad_adjust_mod"
2D fields "dust_rad_adjust_prev" and "dust_rad_adjust_next" required to compute coefficient "dust_rad_adjust" need to be stored in (re)startfi files

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

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

Download in other formats:

Original Format