Context Navigation

physdem.F @ 1047

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

Mars GCM:

IMPORTANT CHANGE: Removed all reference/use of ngridmx (dimphys.h) in routines (necessary prerequisite to using parallel dynamics); in most cases this just means adding 'ngrid' as routine argument, and making local saved variables allocatable (and allocated at first call). In the process, had to convert many *.h files to equivalent modules: yomaer.h => yomaer_h.F90 , surfdat.h => surfdat_h.F90 , comsaison.h => comsaison_h.F90 , yomlw.h => yomlw_h.F90 , comdiurn.h => comdiurn_h.F90 , dimradmars.h => dimradmars_mod.F90 , comgeomfi.h => comgeomfi_h.F90, comsoil.h => comsoil_h.F90 , slope.h => slope_mod.F90
Also updated EOF routines, everything is now in eofdump_mod.F90
Removed unused routine lectfux.F (in dyn3d)

File size: 30.6 KB

Rev	Line
[1047]	1	subroutine physdem0(filename,lonfi,latfi,nsoil,ngrid,nlay,nq,
[999]	2	. phystep,day_ini,time,airefi,
	3	. alb,ith,pzmea,pzstd,pzsig,pzgam,pzthe)
	4
[1036]	5	use infotrac, only: nqtot, tnom
[1047]	6	use comsoil_h, only: inertiedat, volcapa, mlayer
	7	use comgeomfi_h, only: area
	8	use surfdat_h, only: albedodat, zmea, zstd, zsig, zgam, zthe,
	9	& z0_default, albedice, emisice, emissiv,
	10	& iceradius, dtemisice, phisfi, z0
	11	use yomaer_h, only: tauvis
[999]	12	implicit none
	13	c
	14	c Variables qui NE possèdent PAS un axe des temps :
	15	c airefi
	16	c alb
	17	c ith
	18	c pzmea
	19	c pzstd
	20	c pzsig
	21	c pzgam
	22	c pzthe
	23	c
	24	c-------------------------------------------------------------
	25	c
	26	c create physics (re-)start data file "restartfi.nc"
	27	c
	28	c
	29	c
	30	#include "dimensions.h"
	31	#include "paramet.h"
	32	c-----------------------------------------------------------------------
	33	#include "comvert.h"
	34	#include "comgeom2.h"
	35	#include "control.h"
	36	#include "comdissnew.h"
	37	#include "logic.h"
	38	#include "ener.h"
	39	#include "netcdf.inc"
[1047]	40	!#include "dimphys.h"
[1036]	41	!#include "advtrac.h"
[999]	42	#include "callkeys.h"
	43	c
	44	INTEGER nid,iq
	45	INTEGER, parameter :: ivap=1
	46	REAL, parameter :: qsolmax= 150.0
	47	character (len=*) :: filename
	48	character (len=7) :: str7
	49
	50	REAL day_ini
[1047]	51	INTEGER,INTENT(IN) :: nsoil,ngrid,nlay,nq
[999]	52	integer ierr,idim1,idim2,idim3,idim4,idim5,idim6,nvarid
	53
	54	c
	55	REAL phystep,time
[1047]	56	REAL latfi(ngrid), lonfi(ngrid)
	57	! REAL champhys(ngrid)
[999]	58	INTEGER length
	59	PARAMETER (length=100)
	60	REAL tab_cntrl(length)
	61
	62	c
	63
	64	! EXTERNAL defrun_new,iniconst,geopot,inigeom,massdair,pression
	65	! EXTERNAL exner_hyb , SSUM
	66	c
	67	#include "serre.h"
	68	#include "clesph0.h"
	69	#include "fxyprim.h"
[1047]	70	!#include "comgeomfi.h"
	71	!#include "surfdat.h"
	72	!#include "comsoil.h"
[999]	73	#include "planete.h"
[1047]	74	!#include "dimradmars.h"
	75	!#include "yomaer.h"
[999]	76	#include "comcstfi.h"
	77
[1047]	78	real airefi(ngrid)
	79	real alb(ngrid),ith(ngrid,nsoil)
	80	real pzmea(ngrid),pzstd(ngrid)
	81	real pzsig(ngrid),pzgam(ngrid),pzthe(ngrid)
[999]	82	integer ig
	83
	84	! flag which identifies if we are using old tracer names (qsurf01,...)
	85	logical :: oldtracernames=.false.
	86	integer :: count
	87	character(len=30) :: txt ! to store some text
	88	! indexes of water vapour & water ice tracers (if any):
	89	integer :: i_h2o_vap=0
	90	integer :: i_h2o_ice=0
	91	c-----------------------------------------------------------------------
	92
	93	! copy airefi(:) to area(:)
[1047]	94	CALL SCOPY(ngrid,airefi,1,area,1)
	95	! note: area() is defined in comgeomfi_h
[999]	96
	97
[1047]	98	DO ig=1,ngrid
	99	albedodat(ig)=alb(ig) ! note: albedodat() is defined in surfdat_h
	100	zmea(ig)=pzmea(ig) ! note: zmea() is defined in surfdat_h
	101	zstd(ig)=pzstd(ig) ! note: zstd() is defined in surfdat_h
	102	zsig(ig)=pzsig(ig) ! note: zsig() is defined in surfdat_h
	103	zgam(ig)=pzgam(ig) ! note: zgam() is defined in surfdat_h
	104	zthe(ig)=pzthe(ig) ! note: zthe() is defined in surfdat_h
[999]	105	ENDDO
	106
[1047]	107	inertiedat(:,:)=ith(:,:) ! note inertiedat() is defined in comsoil_h
[999]	108	c
	109	c things to store in the physics start file:
	110	c
	111	ierr = NF_CREATE(adjustl(filename),
	112	& IOR(NF_CLOBBER,NF_64BIT_OFFSET), nid)
	113	IF (ierr.NE.NF_NOERR) THEN
	114	WRITE(6,*)'physdem0: Problem creating file ',adjustl(filename)
	115	write(6,*) NF_STRERROR(ierr)
	116	CALL ABORT
	117	ENDIF
	118	c
	119	ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 18,
	120	. "Physics start file")
	121	c
	122	ierr = NF_DEF_DIM (nid,"index",length,idim1)
	123	if (ierr.ne.NF_NOERR) then
	124	WRITE(6,*)'physdem0: Problem defining index dimension'
	125	write(6,*) NF_STRERROR(ierr)
	126	call abort
	127	endif
	128	c
[1047]	129	ierr = NF_DEF_DIM (nid,"physical_points",ngrid,idim2)
[999]	130	if (ierr.ne.NF_NOERR) then
	131	WRITE(6,*)'physdem0: Problem defining physical_points dimension'
	132	write(6,*) NF_STRERROR(ierr)
	133	call abort
	134	endif
	135	c
	136	ierr = NF_DEF_DIM (nid,"subsurface_layers",nsoil,idim3)
	137	if (ierr.ne.NF_NOERR) then
	138	WRITE(6,*)'physdem0: Problem defining subsurface_layers dimension'
	139	write(6,*) NF_STRERROR(ierr)
	140	call abort
	141	endif
	142	c
	143	! ierr = NF_DEF_DIM (nid,"nlayer+1",llm+1,idim4)
	144	ierr = NF_DEF_DIM (nid,"nlayer_plus_1",llm+1,idim4)
	145	if (ierr.ne.NF_NOERR) then
	146	WRITE(6,*)'physdem0: Problem defining nlayer+1 dimension'
	147	write(6,*) NF_STRERROR(ierr)
	148	call abort
	149	endif
	150	c
	151	ierr = NF_DEF_DIM (nid,"number_of_advected_fields",nq,idim5)
	152	if (ierr.ne.NF_NOERR) then
	153	WRITE(6,*)'physdem0: Problem defining advected fields dimension'
	154	WRITE(6,*)' nq = ',nq,' and ierr = ', ierr
	155	write(6,*) NF_STRERROR(ierr)
	156	endif
	157	c
	158	ierr = NF_DEF_DIM (nid,"Time",NF_UNLIMITED,idim6)
	159	if (ierr.ne.NF_NOERR) then
	160	WRITE(6,*)'physdem0: Problem defining time dimension'
	161	WRITE(6,*)' ierr = ', ierr
	162	write(6,*) NF_STRERROR(ierr)
	163	endif
	164
	165	ierr = NF_ENDDEF(nid) ! exit NetCDF define mode
	166
	167	c clear tab_cntrl(:) array
	168	DO ierr = 1, length
	169	tab_cntrl(ierr) = 0.0
	170	ENDDO
	171
[1047]	172	write(,) "physdem0: ngrid: ",ngrid
[999]	173
	174	ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	175	c Fill control array tab_cntrl(:) with parameters for this run
	176	ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	177	c Informations on the physics grid
[1047]	178	tab_cntrl(1) = float(ngrid) ! number of nodes on physics grid
	179	tab_cntrl(2) = float(nlay) ! number of atmospheric layers
[999]	180	tab_cntrl(3) = day_ini + int(time) ! initial day
	181	tab_cntrl(4) = time -int(time) ! initial time of day
	182
	183	c Informations about Mars, used by dynamics and physics
	184	tab_cntrl(5) = rad ! radius of Mars (m) ~3397200
	185	tab_cntrl(6) = omeg ! rotation rate (rad.s-1)
	186	tab_cntrl(7) = g ! gravity (m.s-2) ~3.72
	187	tab_cntrl(8) = mugaz ! Molar mass of the atmosphere (g.mol-1) ~43.49
	188	tab_cntrl(9) = rcp ! = r/cp ~0.256793 (=kappa dans dynamique)
	189	tab_cntrl(10) = daysec ! length of a sol (s) ~88775
	190
	191	tab_cntrl(11) = phystep ! time step in the physics
	192	tab_cntrl(12) = 0.
	193	tab_cntrl(13) = 0.
	194
	195	c Informations about Mars, only for physics
	196	tab_cntrl(14) = year_day ! length of year (sols) ~668.6
	197	tab_cntrl(15) = periheli ! min. Sun-Mars distance (Mkm) ~206.66
	198	tab_cntrl(16) = aphelie ! max. SUn-Mars distance (Mkm) ~249.22
	199	tab_cntrl(17) = peri_day ! date of perihelion (sols since N. spring)
	200	tab_cntrl(18) = obliquit ! Obliquity of the planet (deg) ~23.98
	201
	202	c Boundary layer and turbulence
	203	tab_cntrl(19) = z0_default ! default surface roughness (m) ~0.01
	204	tab_cntrl(20) = lmixmin ! mixing length ~100
	205	tab_cntrl(21) = emin_turb ! minimal energy ~1.e-8
	206
	207	c Optical properties of polar caps and ground emissivity
	208	tab_cntrl(22) = albedice(1) ! Albedo of northern cap ~0.5
	209	tab_cntrl(23) = albedice(2) ! Albedo of southern cap ~0.5
	210	tab_cntrl(24) = emisice(1) ! Emissivity of northern cap ~0.95
	211	tab_cntrl(25) = emisice(2) ! Emissivity of southern cap ~0.95
	212	tab_cntrl(26) = emissiv ! Emissivity of martian soil ~.95
	213	tab_cntrl(31) = iceradius(1) ! mean scat radius of CO2 snow (north)
	214	tab_cntrl(32) = iceradius(2) ! mean scat radius of CO2 snow (south)
	215	tab_cntrl(33) = dtemisice(1) ! time scale for snow metamorphism (north)
	216	tab_cntrl(34) = dtemisice(2) ! time scale for snow metamorphism (south)
	217
	218	c dust aerosol properties
	219	tab_cntrl(27) = tauvis ! mean visible optical depth
	220
	221	tab_cntrl(28) = 0.
	222	tab_cntrl(29) = 0.
	223	tab_cntrl(30) = 0.
	224
	225	! Soil properties:
	226	tab_cntrl(35) = volcapa ! soil volumetric heat capacity
	227
	228	c
	229
	230	ierr = NF_REDEF (nid) ! Enter NetCDF (re-)define mode
	231	IF (ierr.NE.NF_NOERR) THEN
	232	PRINT*, 'physdem0: Failed to swich to NetCDF define mode'
	233	CALL abort
	234	ENDIF
	235	! define variable
	236	#ifdef NC_DOUBLE
	237	ierr = NF_DEF_VAR (nid, "controle", NF_DOUBLE, 1, idim1,nvarid)
	238	#else
	239	ierr = NF_DEF_VAR (nid, "controle", NF_FLOAT, 1, idim1,nvarid)
	240	#endif
	241	IF (ierr.NE.NF_NOERR) THEN
	242	PRINT*, 'physdem0: Failed to define controle'
	243	CALL abort
	244	ENDIF
	245	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 18,
	246	. "Control parameters")
	247	IF (ierr.NE.NF_NOERR) THEN
	248	PRINT*, 'physdem0: Failed to define controle title attribute'
	249	CALL abort
	250	ENDIF
	251	ierr = NF_ENDDEF(nid) ! Leave NetCDF define mode
	252	IF (ierr.NE.NF_NOERR) THEN
	253	PRINT*, 'physdem0: Failed to swich out of NetCDF define mode'
	254	CALL abort
	255	ENDIF
	256	! write variable
	257	#ifdef NC_DOUBLE
	258	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl)
	259	#else
	260	ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl)
	261	#endif
	262	IF (ierr.NE.NF_NOERR) THEN
	263	PRINT*, 'physdem0: Failed to write controle data'
	264	CALL abort
	265	ENDIF
	266
	267
[1047]	268	! write mid-layer depths mlayer() !known from comsoil_h
[999]	269
	270	ierr = NF_REDEF (nid) ! Enter NetCDF (re-)define mode
	271	! define variable
	272	#ifdef NC_DOUBLE
	273	ierr = NF_DEF_VAR (nid,"soildepth",NF_DOUBLE,1,idim3,nvarid)
	274	#else
	275	ierr = NF_DEF_VAR (nid,"soildepth",NF_FLOAT,1,idim3,nvarid)
	276	#endif
	277	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 20,
	278	. "Soil mid-layer depth")
	279	ierr = NF_ENDDEF(nid) ! Leave NetCDF define mode
	280	! write variable
	281	#ifdef NC_DOUBLE
	282	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,mlayer)
	283	#else
	284	ierr = NF_PUT_VAR_REAL (nid,nvarid,mlayer)
	285	#endif
	286
	287	c
	288
	289	ierr = NF_REDEF (nid)
	290	#ifdef NC_DOUBLE
	291	ierr = NF_DEF_VAR (nid, "longitude", NF_DOUBLE, 1, idim2,nvarid)
	292	#else
	293	ierr = NF_DEF_VAR (nid, "longitude", NF_FLOAT, 1, idim2,nvarid)
	294	#endif
	295	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 26,
	296	. "Longitudes of physics grid")
	297	ierr = NF_ENDDEF(nid)
	298
	299	#ifdef NC_DOUBLE
	300	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,lonfi)
	301	#else
	302	ierr = NF_PUT_VAR_REAL (nid,nvarid,lonfi)
	303	#endif
	304
	305	c
	306
	307	ierr = NF_REDEF (nid)
	308	#ifdef NC_DOUBLE
	309	ierr = NF_DEF_VAR (nid, "latitude", NF_DOUBLE, 1, idim2,nvarid)
	310	#else
	311	ierr = NF_DEF_VAR (nid, "latitude", NF_FLOAT, 1, idim2,nvarid)
	312	#endif
	313	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 25,
	314	. "Latitudes of physics grid")
	315	ierr = NF_ENDDEF(nid)
	316	#ifdef NC_DOUBLE
	317	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,latfi)
	318	#else
	319	ierr = NF_PUT_VAR_REAL (nid,nvarid,latfi)
	320	#endif
	321
	322	c
	323
	324	ierr = NF_REDEF (nid)
	325	#ifdef NC_DOUBLE
	326	ierr = NF_DEF_VAR (nid, "area", NF_DOUBLE, 1, idim2,nvarid)
	327	#else
	328	ierr = NF_DEF_VAR (nid, "area", NF_FLOAT, 1, idim2,nvarid)
	329	#endif
	330	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 9,
	331	. "Mesh area")
	332	ierr = NF_ENDDEF(nid)
	333	#ifdef NC_DOUBLE
	334	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,area)
	335	#else
	336	ierr = NF_PUT_VAR_REAL (nid,nvarid,area)
	337	#endif
	338
	339	c
	340
	341	ierr = NF_REDEF (nid)
	342	#ifdef NC_DOUBLE
	343	ierr = NF_DEF_VAR (nid, "phisfi", NF_DOUBLE, 1, idim2,nvarid)
	344	#else
	345	ierr = NF_DEF_VAR (nid, "phisfi", NF_FLOAT, 1, idim2,nvarid)
	346	#endif
	347	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 27,
	348	. "Geopotential at the surface")
	349	ierr = NF_ENDDEF(nid)
	350	#ifdef NC_DOUBLE
	351	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,phisfi)
	352	#else
	353	ierr = NF_PUT_VAR_REAL (nid,nvarid,phisfi)
	354	#endif
	355
	356	c
	357
	358	ierr = NF_REDEF (nid)
	359	#ifdef NC_DOUBLE
	360	ierr = NF_DEF_VAR (nid, "albedodat", NF_DOUBLE, 1, idim2,nvarid)
	361	#else
	362	ierr = NF_DEF_VAR (nid, "albedodat", NF_FLOAT, 1, idim2,nvarid)
	363	#endif
	364	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 21,
	365	. "Albedo of bare ground")
	366	ierr = NF_ENDDEF(nid)
	367	#ifdef NC_DOUBLE
	368	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,albedodat)
	369	#else
	370	ierr = NF_PUT_VAR_REAL (nid,nvarid,albedodat)
	371	#endif
	372
	373	c
	374	c some data for Francois Lott's programs
	375	c
	376
	377	ierr = NF_REDEF (nid)
	378	#ifdef NC_DOUBLE
	379	ierr = NF_DEF_VAR (nid, "ZMEA", NF_DOUBLE, 1, idim2,nvarid)
	380	#else
	381	ierr = NF_DEF_VAR (nid, "ZMEA", NF_FLOAT, 1, idim2,nvarid)
	382	#endif
	383	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
	384	. "Relief: mean relief")
	385	ierr = NF_ENDDEF(nid)
	386	#ifdef NC_DOUBLE
	387	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zmea)
	388	#else
	389	ierr = NF_PUT_VAR_REAL (nid,nvarid,zmea)
	390	#endif
	391	c
	392	ierr = NF_REDEF (nid)
	393	#ifdef NC_DOUBLE
	394	ierr = NF_DEF_VAR (nid, "ZSTD", NF_DOUBLE, 1, idim2,nvarid)
	395	#else
	396	ierr = NF_DEF_VAR (nid, "ZSTD", NF_FLOAT, 1, idim2,nvarid)
	397	#endif
	398	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 26,
	399	. "Relief: standard deviation")
	400	ierr = NF_ENDDEF(nid)
	401	#ifdef NC_DOUBLE
	402	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zstd)
	403	#else
	404	ierr = NF_PUT_VAR_REAL (nid,nvarid,zstd)
	405	#endif
	406	c
	407	ierr = NF_REDEF (nid)
	408	#ifdef NC_DOUBLE
	409	ierr = NF_DEF_VAR (nid, "ZSIG", NF_DOUBLE, 1, idim2,nvarid)
	410	#else
	411	ierr = NF_DEF_VAR (nid, "ZSIG", NF_FLOAT, 1, idim2,nvarid)
	412	#endif
	413	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
	414	. "Relief: sigma parameter")
	415	ierr = NF_ENDDEF(nid)
	416	#ifdef NC_DOUBLE
	417	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zsig)
	418	#else
	419	ierr = NF_PUT_VAR_REAL (nid,nvarid,zsig)
	420	#endif
	421	c
	422	ierr = NF_REDEF (nid)
	423	#ifdef NC_DOUBLE
	424	ierr = NF_DEF_VAR (nid, "ZGAM", NF_DOUBLE, 1, idim2,nvarid)
	425	#else
	426	ierr = NF_DEF_VAR (nid, "ZGAM", NF_FLOAT, 1, idim2,nvarid)
	427	#endif
	428	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
	429	. "Relief: gamma parameter")
	430	ierr = NF_ENDDEF(nid)
	431	#ifdef NC_DOUBLE
	432	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zgam)
	433	#else
	434	ierr = NF_PUT_VAR_REAL (nid,nvarid,zgam)
	435	#endif
	436	c
	437	ierr = NF_REDEF (nid)
	438	#ifdef NC_DOUBLE
	439	ierr = NF_DEF_VAR (nid, "ZTHE", NF_DOUBLE, 1, idim2,nvarid)
	440	#else
	441	ierr = NF_DEF_VAR (nid, "ZTHE", NF_FLOAT, 1, idim2,nvarid)
	442	#endif
	443	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
	444	. "Relief: theta parameter")
	445	ierr = NF_ENDDEF(nid)
	446	#ifdef NC_DOUBLE
	447	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zthe)
	448	#else
	449	ierr = NF_PUT_VAR_REAL (nid,nvarid,zthe)
	450	#endif
	451
	452
	453	! Soil Thermal inertia
	454
	455	ierr = NF_REDEF (nid)
	456	#ifdef NC_DOUBLE
	457	ierr = NF_DEF_VAR (nid,"inertiedat",NF_DOUBLE,
	458	& 2,(/idim2,idim3/),nvarid)
	459	#else
	460	ierr = NF_DEF_VAR (nid,"inertiedat",NF_FLOAT,
	461	& 2,(/idim2,idim3/),nvarid)
	462	#endif
	463	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 20,
	464	. "Soil thermal inertia")
	465	ierr = NF_ENDDEF(nid)
	466	#ifdef NC_DOUBLE
	467	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,inertiedat)
	468	#else
	469	ierr = NF_PUT_VAR_REAL (nid,nvarid,inertiedat)
	470	#endif
	471
[1047]	472	! surface roughness length (z0 is a common in surfdat_h)
[999]	473
	474	ierr = NF_REDEF (nid)
	475	#ifdef NC_DOUBLE
	476	ierr = NF_DEF_VAR (nid, "z0", NF_DOUBLE, 1, idim2,nvarid)
	477	#else
	478	ierr = NF_DEF_VAR (nid, "z0", NF_FLOAT, 1, idim2,nvarid)
	479	#endif
	480	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 24,
	481	. "Surface roughness length")
	482	ierr = NF_ENDDEF(nid)
	483	#ifdef NC_DOUBLE
	484	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,z0)
	485	#else
	486	ierr = NF_PUT_VAR_REAL (nid,nvarid,z0)
	487	#endif
	488
	489
	490
	491
	492	c--------------------- Time dependent variables ---------------------
	493	c-------------------------------------------------------------------
	494	c-------------------------------------------------------------------
	495
	496
	497	! Time
	498
	499	ierr = NF_REDEF (nid)
	500	#ifdef NC_DOUBLE
	501	ierr = NF_DEF_VAR (nid, "Time", NF_DOUBLE, 1,
	502	. (/idim6/), nvarid)
	503	#else
	504	ierr = NF_DEF_VAR (nid, "Time", NF_FLOAT, 1,
	505	. (/idim6/), nvarid)
	506	#endif
	507	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
	508	. "Temps de simulation")
	509	ierr = NF_ENDDEF(nid)
	510
	511
	512	! CO2 Ice Cover
	513
	514	ierr = NF_REDEF (nid)
	515	#ifdef NC_DOUBLE
	516	ierr = NF_DEF_VAR (nid, "co2ice", NF_DOUBLE, 2,
	517	. (/idim2,idim6/), nvarid)
	518	#else
	519	ierr = NF_DEF_VAR (nid, "co2ice", NF_FLOAT, 2,
	520	. (/idim2,idim6/), nvarid)
	521	#endif
	522	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 13,
	523	. "CO2 ice cover")
	524	ierr = NF_ENDDEF(nid)
	525
	526
	527	! Surface temperature
	528
	529	ierr = NF_REDEF (nid)
	530	#ifdef NC_DOUBLE
	531	ierr = NF_DEF_VAR (nid, "tsurf", NF_DOUBLE, 2,
	532	. (/idim2,idim6/), nvarid)
	533	#else
	534	ierr = NF_DEF_VAR (nid, "tsurf", NF_FLOAT, 2,
	535	. (/idim2,idim6/), nvarid)
	536	#endif
	537	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
	538	. "Surface temperature")
	539	ierr = NF_ENDDEF(nid)
	540
	541
	542	! Soil temperature
	543
	544	ierr = NF_REDEF (nid)
	545	#ifdef NC_DOUBLE
	546	ierr = NF_DEF_VAR (nid, "tsoil", NF_DOUBLE, 3,
	547	. (/idim2,idim3,idim6/), nvarid)
	548	#else
	549	ierr = NF_DEF_VAR (nid, "tsoil", NF_FLOAT, 3,
	550	. (/idim2,idim3,idim6/), nvarid)
	551	#endif
	552	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 16,
	553	. "Soil temperature")
	554	ierr = NF_ENDDEF(nid)
	555
	556
	557	c emissivity
	558
	559	ierr = NF_REDEF (nid)
	560	#ifdef NC_DOUBLE
	561	ierr = NF_DEF_VAR (nid,"emis", NF_DOUBLE, 2,
	562	. (/idim2,idim6/), nvarid)
	563	#else
	564	ierr = NF_DEF_VAR (nid, "emis", NF_FLOAT, 2,
	565	. (/idim2,idim6/), nvarid)
	566	#endif
	567	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 18,
	568	. "Surface emissivity")
	569	ierr = NF_ENDDEF(nid)
	570
	571
	572	c planetary boundary layer
	573
	574	ierr = NF_REDEF (nid)
	575	#ifdef NC_DOUBLE
	576	ierr = NF_DEF_VAR (nid, "q2", NF_DOUBLE, 3,
	577	. (/idim2,idim4,idim6/), nvarid)
	578	#else
	579	ierr = NF_DEF_VAR (nid, "q2", NF_FLOAT, 3,
	580	. (/idim2,idim4,idim6/), nvarid)
	581	#endif
	582	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 17,
	583	. "pbl wind variance")
	584	ierr = NF_ENDDEF(nid)
	585
	586
	587
	588
	589	c tracers
	590	!
	591	! Preliminary stuff: check if tracers follow old naming convention (qsurf01,
	592	! qsurf02, ...)
	593	count=0
[1036]	594	do iq=1,nqtot
[999]	595	txt= " "
	596	write(txt,'(a1,i2.2)')'q',iq
	597	if (txt.ne.tnom(iq)) then ! use tracer names stored in dynamics
	598	! did not find old tracer name
	599	exit ! might as well stop here
	600	else
	601	! found old tracer name
	602	count=count+1
	603	endif
	604	enddo
[1036]	605	if (count.eq.nqtot) then
[999]	606	write(,) "physdem0:tracers seem to follow old naming ",
	607	& "convention (qsurf01,qsurf02,...)"
	608	write(,) " => will work for now ... "
	609	write(,) " but you should run newstart to rename them"
	610	oldtracernames=.true.
[1036]	611	endif ! of if (count.eq.nqtot)
[999]	612
	613	IF(nq.GE.1) THEN
	614	! preliminary stuff: look for water vapour & water ice tracers (if any)
	615	if (.not.oldtracernames) then
	616	do iq=1,nq
	617	if (tnom(iq).eq."h2o_vap") then
	618	i_h2o_vap=iq
	619	endif
	620	if (tnom(iq).eq."h2o_ice") then
	621	i_h2o_ice=iq
	622	endif
	623	enddo ! of iq=1,nq
	624	! handle special case of only water vapour tracer (no ice)
	625	if ((i_h2o_vap.ne.0).and.(i_h2o_ice.eq.0)) then
	626	! then the index of (surface) ice is i_h2o_vap
	627	i_h2o_ice=i_h2o_vap
	628	endif
	629	endif ! of if (.not.oldtracernames)
	630
	631	DO iq=1,nq
	632	IF (oldtracernames) THEN
	633	txt=" "
	634	write(txt,'(a5,i2.2)')'qsurf',iq
	635	ELSE
	636	txt=tnom(iq)
	637	! Exception: there is no water vapour surface tracer
	638	if (txt.eq."h2o_vap") then
	639	write(,)"physdem0: skipping water vapour tracer"
	640	if (i_h2o_ice.eq.i_h2o_vap) then
	641	! then there is no "water ice" tracer; but still
	642	! there is some water ice on the surface
	643	write(,)" writing water ice instead"
	644	txt="h2o_ice"
	645	else
	646	! there is a "water ice" tracer which has been / will be
	647	! delt with in due time
	648	cycle
	649	endif ! of if (igcm_h2o_ice.eq.igcm_h2o_vap)
	650	endif ! of if (txt.eq."h2o_vap")
	651	ENDIF ! of IF (oldtracernames)
	652
	653	ierr=NF_REDEF(nid)
	654	IF (ierr.NE.NF_NOERR) THEN
	655	PRINT*, 'physdem0: Failed to swich to NetCDF define mode'
	656	CALL abort
	657	ENDIF
	658	#ifdef NC_DOUBLE
	659	ierr=NF_DEF_VAR(nid,txt,NF_DOUBLE,2,(/idim2,idim6/),nvarid)
	660	#else
	661	ierr=NF_DEF_VAR(nid,txt,NF_FLOAT,2,(/idim2,idim6/),nvarid)
	662	#endif
	663	IF (ierr.NE.NF_NOERR) THEN
	664	PRINT*, 'physdem0: Failed to define ',trim(txt)
	665	CALL abort
	666	ENDIF
	667	ierr=NF_PUT_ATT_TEXT (nid, nvarid, "title", 17,
	668	& "tracer on surface")
	669	IF (ierr.NE.NF_NOERR) THEN
	670	PRINT*, 'physdem0: Failed to define ',trim(txt),
	671	& ' title attribute'
	672	CALL abort
	673	ENDIF
	674	ierr=NF_ENDDEF(nid)
	675	IF (ierr.NE.NF_NOERR) THEN
	676	PRINT*, 'physdem0: Failed to swich out of define mode'
	677	CALL abort
	678	ENDIF
	679	ENDDO
	680	ENDIF
	681
	682	c close file
	683	ierr = NF_CLOSE(nid)
	684
	685
	686
	687	END
	688
	689
	690
	691	cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	692	cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	693	cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	694	cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	695
	696
	697
	698
[1047]	699	subroutine physdem1(filename,nsoil,ngrid,nlay,nq,
[999]	700	. phystep,time,
	701	. tsurf,tsoil,co2ice,emis,q2,qsurf)
[1036]	702	use infotrac, only: nqtot, tnom
[1047]	703	use surfdat_h, only: emis
[999]	704	implicit none
	705	c
	706	c Variables qui possèdent un axe des temps :
	707	c tsurf
	708	c tsoil
	709	c co2ice
	710	c q2
	711	c qsurf
	712	c emis (modifié dans newcondens --> co2snow)
	713	c
	714	c-------------------------------------------------------------
	715	c
	716	c write physics (re-)start data file "restartfi.nc"
	717	c
	718	c
	719	c
	720	#include "dimensions.h"
	721	#include "paramet.h"
	722	c-----------------------------------------------------------------------
	723	#include "comvert.h"
	724	#include "comgeom2.h"
	725	#include "control.h"
	726	#include "comdissnew.h"
	727	#include "logic.h"
	728	#include "ener.h"
	729	#include "netcdf.inc"
[1047]	730	!#include "dimphys.h"
[1036]	731	!#include "advtrac.h"
[999]	732	#include "callkeys.h"
	733	c
	734	INTEGER nid,iq
	735	INTEGER, parameter :: ivap=1
	736	REAL, parameter :: qsolmax= 150.0
	737	character (len=*) :: filename
	738	character (len=7) :: str7
	739
	740	REAL day_ini
[1047]	741	INTEGER,INTENT(IN) :: nsoil,ngrid,nlay,nq
[999]	742	integer ierr,nvarid
	743
	744	c
	745	REAL phystep,time
[1047]	746	! REAL champhys(ngrid)
	747	REAL tsurf(ngrid)
[999]	748
	749	INTEGER nb
	750	INTEGER edges(3),corner(3)
	751
	752
	753
	754	c
	755
	756	! EXTERNAL defrun_new,iniconst,geopot,inigeom,massdair,pression
	757	! EXTERNAL exner_hyb , SSUM
	758	c
	759	#include "serre.h"
	760	#include "clesph0.h"
	761	#include "fxyprim.h"
[1047]	762	!#include "comgeomfi.h"
	763	!#include "surfdat.h"
	764	!#include "comsoil.h"
[999]	765	#include "planete.h"
[1047]	766	!#include "dimradmars.h"
	767	!#include "yomaer.h"
[999]	768	#include "comcstfi.h"
	769
[1047]	770	real co2ice(ngrid),tsoil(ngrid,nsoil),emis(ngrid)
	771	real q2(ngrid,nlay+1),qsurf(ngrid,nq)
[999]	772	integer ig
	773
	774	! flag which identifies if we are using old tracer names (qsurf01,...)
	775	logical,save :: oldtracernames=.false.
	776	logical,save :: firstcall = .true.
	777	integer :: count
	778
	779	character(len=30) :: txt ! to store some text
	780	! indexes of water vapour & water ice tracers (if any):
	781	integer,save :: i_h2o_vap=0
	782	integer,save :: i_h2o_ice=0
	783	c-----------------------------------------------------------------------
	784
	785
	786	ierr = NF_OPEN(filename, NF_WRITE, nid)
	787	IF (ierr .NE. NF_NOERR) THEN
	788	PRINT*, "Pb. d ouverture "//filename
	789	CALL abort
	790	ENDIF
	791
	792	c On a single run, different files can be written with physdem1.
	793	c Therefore, get the last time index from the file itself:
	794	ierr = NF_INQ_DIMID(nid,"Time",nvarid)
	795	ierr = NF_INQ_DIMLEN(nid,nvarid,nb)
	796
	797	c Ecriture/extension de la coordonnee temps
	798
	799	nb = nb + 1
	800	ierr = NF_INQ_VARID(nid, "Time", nvarid)
	801	IF (ierr .NE. NF_NOERR) THEN
	802	print *, NF_STRERROR(ierr)
	803	print*,'Variable Time n est pas definie'
	804	CALL abort
	805	ENDIF
	806	#ifdef NC_DOUBLE
	807	ierr = NF_PUT_VAR1_DOUBLE (nid,nvarid,nb,time)
	808	#else
	809	ierr = NF_PUT_VAR1_REAL (nid,nvarid,nb,time)
	810	#endif
	811	IF (ierr .NE. NF_NOERR) THEN
	812	print*, "Erreur ecriture temps!!"
	813	print*, NF_STRERROR(ierr)
	814	ENDIF
	815	!PRINT*, "Enregistrement pour ", nb, time
	816
	817
	818	c Write the physical fields
	819
	820	! Soil temperature
	821	corner(1)=1
	822	corner(2)=1
	823	corner(3)=nb
[1047]	824	edges(1)=ngrid
[999]	825	edges(2)=nsoil
	826	edges(3)=1
	827	ierr = NF_INQ_VARID(nid, "tsoil", nvarid)
	828	IF (ierr .NE. NF_NOERR) THEN
	829	PRINT*, "Variable tsoil n est pas definie"
	830	CALL abort
	831	ENDIF
	832	#ifdef NC_DOUBLE
	833	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,tsoil)
	834	#else
	835	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,tsoil)
	836	#endif
	837	IF (ierr .NE. NF_NOERR) THEN
	838	print*, "Erreur ecriture tsoil!!"
	839	print*, NF_STRERROR(ierr)
	840	ENDIF
	841
	842
	843	c planetary boundary layer
	844	corner(1)=1
	845	corner(2)=1
	846	corner(3)=nb
[1047]	847	edges(1)=ngrid
[999]	848	edges(2)=llm+1
	849	edges(3)=1
	850	ierr = NF_INQ_VARID(nid, "q2", nvarid)
	851	IF (ierr .NE. NF_NOERR) THEN
	852	PRINT*, "Variable q2 n est pas definie"
	853	CALL abort
	854	ENDIF
	855	#ifdef NC_DOUBLE
	856	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,q2)
	857	#else
	858	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,q2)
	859	#endif
	860	IF (ierr .NE. NF_NOERR) THEN
	861	print*, "Erreur ecriture q2!!"
	862	print*, NF_STRERROR(ierr)
	863	ENDIF
	864
	865
	866
	867	c Following corner and egdes are the same for co2ice,tsurf,qusrf,emis and tracers
	868	corner(1)=1
	869	corner(2)=nb
[1047]	870	edges(1)=ngrid
[999]	871	edges(2)=1
	872
	873	! CO2 Ice Cover
	874	ierr = NF_INQ_VARID(nid, "co2ice", nvarid)
	875	IF (ierr .NE. NF_NOERR) THEN
	876	PRINT*, "Variable co2ice n est pas definie"
	877	CALL abort
	878	ENDIF
	879	#ifdef NC_DOUBLE
	880	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,co2ice)
	881	#else
	882	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,co2ice)
	883	#endif
	884	IF (ierr .NE. NF_NOERR) THEN
	885	print*, "Erreur ecriture co2ice!!"
	886	print*, NF_STRERROR(ierr)
	887	ENDIF
	888
	889	! Surface temperature
	890	ierr = NF_INQ_VARID(nid, "tsurf", nvarid)
	891	IF (ierr .NE. NF_NOERR) THEN
	892	PRINT*, "Variable tsurf n est pas definie"
	893	CALL abort
	894	ENDIF
	895	#ifdef NC_DOUBLE
	896	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,tsurf)
	897	#else
	898	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,tsurf)
	899	#endif
	900	IF (ierr .NE. NF_NOERR) THEN
	901	print*, "Erreur ecriture tsurf!!"
	902	print*, NF_STRERROR(ierr)
	903	ENDIF
	904
	905
	906	c emissivity
	907	ierr = NF_INQ_VARID(nid, "emis", nvarid)
	908	IF (ierr .NE. NF_NOERR) THEN
	909	PRINT*, "Variable emis n est pas definie"
	910	CALL abort
	911	ENDIF
	912	#ifdef NC_DOUBLE
	913	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,emis)
	914	#else
	915	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,emis)
	916	#endif
	917	IF (ierr .NE. NF_NOERR) THEN
	918	print*, "Erreur ecriture emis!!"
	919	print*, NF_STRERROR(ierr)
	920	ENDIF
	921
	922
	923	c tracers
	924
	925	! Preliminary stuff: check if tracers follow old naming convention (qsurf01,
	926	! qsurf02, ...)
	927	IF (firstcall) THEN
	928	count=0
[1036]	929	do iq=1,nqtot
[999]	930	txt= " "
	931	write(txt,'(a1,i2.2)')'q',iq
	932	if (txt.ne.tnom(iq)) then ! use tracer names stored in dynamics
	933	! did not find old tracer name
	934	exit ! might as well stop here
	935	else
	936	! found old tracer name
	937	count=count+1
	938	endif
	939	enddo
[1036]	940	if (count.eq.nqtot) then
[999]	941	write(,) "physdem1:tracers seem to follow old naming ",
	942	& "convention (qsurf01,qsurf02,...)"
	943	write(,) " => will work for now ... "
	944	write(,) " but you should run newstart to rename them"
	945	oldtracernames=.true.
[1036]	946	endif ! of if (count.eq.nqtot)
[999]	947	ENDIF ! of if(firstcall)
	948
	949	! If computing water cycle with ice, move surface ice
[1036]	950	! back to qsurf(nqtot)
[999]	951	IF (oldtracernames .and. water) THEN
	952	!"loop" to avoid potential out-of-bounds on arrays
[1036]	953	write(,)'physdem1: moving surface water ice to index ',nqtot
	954	do iq=nqtot,nqtot
[1047]	955	qsurf(1:ngrid,iq)=qsurf(1:ngrid,iq-1)
	956	qsurf(1:ngrid,iq-1)=0
[999]	957	enddo
	958	ENDIF
	959
	960
	961
	962	IF(nq.GE.1) THEN
	963
	964	IF (firstcall) THEN
	965	! preliminary stuff: look for water vapour & water ice tracers (if any)
	966	if (.not.oldtracernames) then
	967	do iq=1,nq
	968	if (tnom(iq).eq."h2o_vap") then
	969	i_h2o_vap=iq
	970	endif
	971	if (tnom(iq).eq."h2o_ice") then
	972	i_h2o_ice=iq
	973	endif
	974	enddo ! of iq=1,nq
	975	! handle special case of only water vapour tracer (no ice)
	976	if ((i_h2o_vap.ne.0).and.(i_h2o_ice.eq.0)) then
	977	! then the index of (surface) ice is i_h2o_vap
	978	i_h2o_ice=i_h2o_vap
	979	endif
	980	endif ! of if (.not.oldtracernames)
	981	ENDIF ! of if(firstcall)
	982	firstcall = .false.
	983
	984
	985	DO iq=1,nq
	986	IF (oldtracernames) THEN
	987	txt=" "
	988	write(txt,'(a5,i2.2)')'qsurf',iq
	989	ELSE
	990	txt=tnom(iq)
	991	! Exception: there is no water vapour surface tracer
	992	if (txt.eq."h2o_vap") then
	993	if (i_h2o_ice.eq.i_h2o_vap) then
	994	! then there is no "water ice" tracer; but still
	995	! there is some water ice on the surface
	996	txt="h2o_ice"
	997	else
	998	! there is a "water ice" tracer which has been / will be
	999	! delt with in due time
	1000	cycle
	1001	endif ! of if (igcm_h2o_ice.eq.igcm_h2o_vap)
	1002	endif ! of if (txt.eq."h2o_vap")
	1003	ENDIF ! of IF (oldtracernames)
	1004
	1005
	1006
	1007	ierr = NF_INQ_VARID(nid, trim(txt), nvarid)
	1008	IF (ierr .NE. NF_NOERR) THEN
	1009	PRINT*, "Variable "//trim(txt)//" n est pas definie"
	1010	CALL abort
	1011	ENDIF
	1012	#ifdef NC_DOUBLE
	1013	ierr = NF_PUT_VARA_DOUBLE (nid,nvarid,corner,edges,
	1014	. qsurf(:,iq))
	1015	#else
	1016	ierr = NF_PUT_VARA_REAL (nid,nvarid,corner,edges,
	1017	. qsurf(:,iq))
	1018	#endif
	1019	IF (ierr .NE. NF_NOERR) THEN
	1020	print*, "Erreur ecriture "//trim(txt)//"!!"
	1021	print*, NF_STRERROR(ierr)
	1022	ENDIF
	1023
	1024	ENDDO ! of DO iq=1,nq
	1025	ENDIF ! of IF(nq.GE.1)
	1026
	1027	c close file
	1028	ierr = NF_CLOSE(nid)
	1029
	1030	RETURN
	1031
	1032	END

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

Download in other formats:

Original Format