Context Navigation

physdem.F @ 1009

Last change on this file since 1009 was 999, checked in by tnavarro, 11 years ago
Possibility to store multiple initial states in one start/startfi. This is RETROCOMPATIBLE. New option ecrithist in run.def to write data in start/startfi every ecrithist dynamical timestep. New option timestart in run.def to initialize the GCM with the time timestart stored in start
File size: 30.1 KB

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

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

Download in other formats:

Original Format