Context Navigation

newstart.F @ 3556

Last change on this file since 3556 was 1523, checked in by emillour, 9 years ago

All models: More updates to make planetary codes (+Earth) setups converge.

in dyn3d_common:
convmas.F => convmas.F90
enercin.F => enercin.F90
flumass.F => flumass.F90
massbar.F => massbar.F90
tourpot.F => tourpot.F90
vitvert.F => vitvert.F90

in misc:
move "q_sat" from "dyn3d_common" to "misc" (in Earth model, it is also called by the physics)
move "write_field" from "dyn3d_common" to "misc"(may be called from physics or dynamics and depends on neither).

in phy_common:
move "write_field_phy" here since it may be called from any physics package)
add module "regular_lonlat_mod" to store global information on lon-lat grid

in dynlonlat_phylonlat/phy*:
turn "iniphysiq.F90" into module "iniphysiq_mod.F90" (and of course adapt gcm.F[90] and 1D models accordingly)

File size: 38.2 KB

Line
1	C======================================================================
2	PROGRAM newstart
3	c=======================================================================
4	c
5	c
6	c Auteur: S. Lebonnois,
7	c a partir des newstart/start_archive/lect_start_archive martiens
8	c
9	c Derniere modif : 02/09 (ecriture des q*)
10	c 01/12 (inclusion dans svn dyn3d)
11	c
12	c Objet: Modify the grid for the initial state (LMD GCM VENUS/TITAN)
13	c ----- (from file NetCDF start_archive.nc)
14	c
15	c
16	c=======================================================================
17
18	use IOIPSL
19	USE filtreg_mod
20	USE startvar
21	USE control_mod
22	USE infotrac
23	use cpdet_mod, only: ini_cpdet,t2tpot
24	use exner_hyb_m, only: exner_hyb
25	use exner_milieu_m, only: exner_milieu
26	USE comconst_mod
27	USE comvert_mod, ONLY: ap,bp,presnivs,pa,preff,nivsigs,nivsig,
28	. aps,bps,scaleheight,pseudoalt,
29	. disvert_type,pressure_exner
30	USE serre_mod, ONLY: clon,clat,grossismx,grossismy,
31	& dzoomx,dzoomy,taux,tauy
32	USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
33	USE logic_mod, ONLY: iflag_trac,fxyhypb,ysinus
34	USE temps_mod, ONLY: day_ref,annee_ref
35	implicit none
36
37	#include "dimensions.h"
38	#include "paramet.h"
39	#include "comdissnew.h"
40	#include "comgeom2.h"
41	#include "description.h"
42	#include "dimsoil.h"
43	#include "netcdf.inc"
44
45	c-----------------------------------------------------------------------
46	c Declarations
47	c-----------------------------------------------------------------------
48
49	c Variables pour fichier "ini"
50	c------------------------------------
51	INTEGER imold,jmold,lmold,nqold,ip1jmp1old
52	INTEGER length
53	parameter (length = 100)
54	real tab_cntrl(2*length)
55	INTEGER isoil,iq,iqmax
56	CHARACTER*2 str2
57
58	c Variable histoire
59	c------------------
60	REAL vcov(iip1,jjm,llm),ucov(iip1,jjp1,llm) ! vents covariants
61	REAL teta(iip1,jjp1,llm),ps(iip1,jjp1)
62	REAL phis(iip1,jjp1) ! geopotentiel au sol
63	REAL masse(ip1jmp1,llm) ! masse de l'atmosphere
64	REAL, ALLOCATABLE, DIMENSION(:,:,:,:):: q! champs advectes
65	REAL tab_cntrl_dyn(length) ! tableau des parametres de start
66
67	c variable physique
68	c------------------
69	integer ngridmx
70	parameter (ngridmx=(2+(jjm-1)*iim - 1/jjm))
71	REAL tab_cntrl_fi(length) ! tableau des parametres de startfi
72	real rlat(ngridmx),rlon(ngridmx)
73	REAL tsurf(ngridmx),tsoil(ngridmx,nsoilmx)
74	REAL albe(ngridmx),radsol(ngridmx),sollw(ngridmx)
75	real solsw(ngridmx),dlw(ngridmx)
76	REAL zmea(ngridmx), zstd(ngridmx)
77	REAL zsig(ngridmx), zgam(ngridmx), zthe(ngridmx)
78	REAL zpic(ngridmx), zval(ngridmx)
79	real t_fi(ngridmx,llm)
80
81	c Variable nouvelle grille naturelle au point scalaire
82	c------------------------------------------------------
83	real us(iip1,jjp1,llm),vs(iip1,jjp1,llm)
84	REAL p3d(iip1,jjp1,llm+1) ! pression aux interfaces
85	REAL phisold_newgrid(iip1,jjp1)
86	REAL T(iip1,jjp1,llm)
87	real rlonS(iip1,jjp1),rlatS(iip1,jjp1)
88	real tsurfS(iip1,jjp1),tsoilS(iip1,jjp1,nsoilmx)
89	real albeS(ip1jmp1),radsolS(ip1jmp1),sollwS(ip1jmp1)
90	real solswS(ip1jmp1),dlwS(ip1jmp1)
91	real zmeaS(ip1jmp1),zstdS(ip1jmp1),zsigS(ip1jmp1)
92	real zgamS(ip1jmp1),ztheS(ip1jmp1),zpicS(ip1jmp1)
93	real zvalS(ip1jmp1)
94
95	real ptotal
96
97	c Var intermediaires : vent naturel, mais pas coord scalaire
98	c-----------------------------------------------------------
99	real vnat(iip1,jjm,llm),unat(iip1,jjp1,llm)
100
101	REAL pks(iip1,jjp1) ! exner (f pour filtre)
102	REAL pk(iip1,jjp1,llm)
103	REAL pkf(iip1,jjp1,llm)
104	REAL alpha(iip1,jjp1,llm),beta(iip1,jjp1,llm)
105
106
107	c Variable de l'ancienne grille
108	c---------------------------------------------------------
109
110	real, dimension(:), allocatable :: rlonuold, rlatvold
111	real, dimension(:), allocatable :: rlonvold, rlatuold
112	real, dimension(:), allocatable :: nivsigsold,nivsigold
113	real, dimension(:), allocatable :: apold,bpold
114	real, dimension(:), allocatable :: presnivsold
115	real, dimension(:,:,:), allocatable :: uold,vold,told
116	real, dimension(:,:,:,:), allocatable :: qold
117	real, dimension(:,:,:), allocatable :: tsoilold
118	real, dimension(:,:), allocatable :: psold,phisold
119	real, dimension(:,:), allocatable :: tsurfold
120	real, dimension(:,:), allocatable :: albeold,radsolold
121	real, dimension(:,:), allocatable :: sollwold,solswold
122	real, dimension(:,:), allocatable :: dlwold
123	real, dimension(:,:), allocatable :: zmeaold,zstdold,zsigold
124	real, dimension(:,:), allocatable :: zgamold,ztheold,zpicold
125	real, dimension(:,:), allocatable :: zvalold
126
127	real ptotalold
128
129	c Variable intermediaires iutilise pour l'extrapolation verticale
130	c----------------------------------------------------------------
131	real, dimension(:,:,:), allocatable :: var,varp1
132
133	c divers local
134	c-----------------
135
136	integer ierr,nid,nvarid
137	INTEGER ij, l,i,j
138	character*80 fichnom
139	integer, dimension(4) :: start,counter
140	REAL phisinverse(iip1,jjp1) ! geopotentiel au sol avant inversion
141	logical topoflag,albedoflag,razvitu,razvitv
142	real albedo
143
144	c=======================================================================
145	c INITIALISATIONS DIVERSES
146	c=======================================================================
147
148	c VENUS/TITAN
149
150	iflag_trac = 1
151	c-----------------------------------------------------------------------
152	c Initialisation des traceurs
153	c ---------------------------
154	c Choix du nombre de traceurs et du schema pour l'advection
155	c dans fichier traceur.def, par default ou via INCA
156	call infotrac_init
157
158	c Allocation de la tableau q : champs advectes
159	allocate(q(iip1,jjp1,llm,nqtot))
160
161	c-----------------------------------------------------------------------
162	c Ouverture du fichier a modifier (start_archive.nc)
163	c-----------------------------------------------------------------------
164
165	write(,) 'Creation d un etat initial a partir de'
166	write(,) './start_archive.nc'
167	write(,)
168	fichnom = 'start_archive.nc'
169	ierr = NF_OPEN (fichnom, NF_NOWRITE,nid)
170	IF (ierr.NE.NF_NOERR) THEN
171	write(6,*)' Pb d''ouverture du fichier ',fichnom
172	write(6,*)' ierr = ', ierr
173	CALL ABORT
174	ENDIF
175
176	c-----------------------------------------------------------------------
177	c Lecture du tableau des parametres du run (pour la dynamique)
178	c-----------------------------------------------------------------------
179
180	write(,) 'lecture tab_cntrl START_ARCHIVE'
181	c
182	ierr = NF_INQ_VARID (nid, "controle", nvarid)
183	#ifdef NC_DOUBLE
184	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tab_cntrl)
185	#else
186	ierr = NF_GET_VAR_REAL(nid, nvarid, tab_cntrl)
187	#endif
188	c
189	write(,) 'Impression de tab_cntrl'
190	do i=1,200
191	write(,) i,tab_cntrl(i)
192	enddo
193
194	c-----------------------------------------------------------------------
195	c Initialisation des constantes
196	c-----------------------------------------------------------------------
197
198	imold = tab_cntrl(1)
199	jmold = tab_cntrl(2)
200	lmold = tab_cntrl(3)
201	day_ref = tab_cntrl(4)
202	annee_ref = tab_cntrl(5)
203	rad = tab_cntrl(6)
204	omeg = tab_cntrl(7)
205	g = tab_cntrl(8)
206	cpp = tab_cntrl(9)
207	kappa = tab_cntrl(10)
208	daysec = tab_cntrl(11)
209	dtvr = tab_cntrl(12)
210	etot0 = tab_cntrl(13)
211	ptot0 = tab_cntrl(14)
212	ztot0 = tab_cntrl(15)
213	stot0 = tab_cntrl(16)
214	ang0 = tab_cntrl(17)
215	pa = tab_cntrl(18)
216	preff = tab_cntrl(19)
217	c
218	clon = tab_cntrl(20)
219	clat = tab_cntrl(21)
220	grossismx = tab_cntrl(22)
221	grossismy = tab_cntrl(23)
222
223	IF ( tab_cntrl(24).EQ.1. ) THEN
224	fxyhypb = . TRUE .
225	dzoomx = tab_cntrl(25)
226	dzoomy = tab_cntrl(26)
227	taux = tab_cntrl(28)
228	tauy = tab_cntrl(29)
229	ELSE
230	fxyhypb = . FALSE .
231	ysinus = . FALSE .
232	IF( tab_cntrl(27).EQ.1. ) ysinus = . TRUE.
233	ENDIF
234
235	ptotalold = tab_cntrl(2*length)
236
237	write(,) "Old dimensions:"
238	write(,) "longitude: ",imold
239	write(,) "latitude: ",jmold
240	write(,) "altitude: ",lmold
241	write(,)
242
243	ip1jmp1old = (imold+1-1/imold)*(jmold+1-1/jmold)
244
245	c dans run.def
246	CALL getin('planet_type',planet_type)
247	call ini_cpdet
248
249	c=======================================================================
250	c CHANGEMENT DE CONSTANTES CONTENUES DANS tab_cntrl
251	c=======================================================================
252	c changement de la resolution dans le fichier de controle
253	tab_cntrl(1) = REAL(iim)
254	tab_cntrl(2) = REAL(jjm)
255	tab_cntrl(3) = REAL(llm)
256
257	c--------------------------------
258	c We use a specific run.def to read new parameters that need to be changed
259	c--------------------------------
260
261	c Changement de cpp:
262	CALL getin('cpp',cpp)
263	kappa = (8.314511/43.44E-3)/cpp
264	tab_cntrl(9) = cpp
265	tab_cntrl(10) = kappa
266
267	c CHANGING THE ZOOM PARAMETERS TO CHANGE THE GRID
268	CALL getin('clon',clon)
269	CALL getin('clat',clat)
270	tab_cntrl(20) = clon
271	tab_cntrl(21) = clat
272	CALL getin('grossismx',grossismx)
273	CALL getin('grossismy',grossismy)
274	tab_cntrl(22) = grossismx
275	tab_cntrl(23) = grossismy
276	CALL getin('fxyhypb',fxyhypb)
277	IF ( fxyhypb ) THEN
278	CALL getin('dzoomx',dzoomx)
279	CALL getin('dzoomy',dzoomy)
280	tab_cntrl(25) = dzoomx
281	tab_cntrl(26) = dzoomy
282	CALL getin('taux',taux)
283	CALL getin('tauy',tauy)
284	tab_cntrl(28) = taux
285	tab_cntrl(29) = tauy
286	ELSE
287	CALL getin('ysinus',ysinus)
288	IF (ysinus) THEN
289	tab_cntrl(27) = 1
290	ELSE
291	tab_cntrl(27) = 0
292	ENDIF
293	ENDIF
294
295	c changes must be done BEFORE these lines...
296	DO l=1,length
297	tab_cntrl_dyn(l)= tab_cntrl(l)
298	tab_cntrl_fi(l) = tab_cntrl(l+length)
299	ENDDO
300
301	c-----------------------------------------------------------------------
302	c Autres initialisations
303	c-----------------------------------------------------------------------
304
305	CALL iniconst
306	CALL inigeom
307	call inifilr
308
309	c-----------------------------------------------------------------------
310	c Allocations des anciennes variables
311	c-----------------------------------------------------------------------
312
313	allocate(rlonuold(imold+1), rlatvold(jmold))
314	allocate(rlonvold(imold+1), rlatuold(jmold+1))
315	allocate(nivsigsold(lmold+1),nivsigold(lmold))
316	allocate(apold(lmold),bpold(lmold))
317	allocate(presnivsold(lmold))
318	allocate(uold(imold+1,jmold+1,lmold))
319	allocate(vold(imold+1,jmold+1,lmold))
320	allocate(told(imold+1,jmold+1,lmold))
321	allocate(qold(imold+1,jmold+1,lmold,nqtot))
322	allocate(psold(imold+1,jmold+1))
323	allocate(phisold(imold+1,jmold+1))
324	allocate(tsurfold(imold+1,jmold+1))
325	allocate(tsoilold(imold+1,jmold+1,nsoilmx))
326	allocate(albeold(imold+1,jmold+1),radsolold(imold+1,jmold+1))
327	allocate(sollwold(imold+1,jmold+1),solswold(imold+1,jmold+1))
328	allocate(dlwold(imold+1,jmold+1))
329	allocate(zmeaold(imold+1,jmold+1),zstdold(imold+1,jmold+1))
330	allocate(zsigold(imold+1,jmold+1),zgamold(imold+1,jmold+1))
331	allocate(ztheold(imold+1,jmold+1),zpicold(imold+1,jmold+1))
332	allocate(zvalold(imold+1,jmold+1))
333
334	allocate(var (imold+1,jmold+1,llm))
335	allocate(varp1 (imold+1,jmold+1,llm+1))
336
337	print*,"Initialisations OK"
338
339	c-----------------------------------------------------------------------
340	c Lecture des longitudes et latitudes
341	c-----------------------------------------------------------------------
342	c
343	ierr = NF_INQ_VARID (nid, "rlonv", nvarid)
344	IF (ierr .NE. NF_NOERR) THEN
345	PRINT*, "new_grid: Le champ <rlonv> est absent de start.nc"
346	CALL abort
347	ENDIF
348	#ifdef NC_DOUBLE
349	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonvold)
350	#else
351	ierr = NF_GET_VAR_REAL(nid, nvarid, rlonvold)
352	#endif
353	IF (ierr .NE. NF_NOERR) THEN
354	PRINT*, "new_grid: Lecture echouee pour <rlonv>"
355	CALL abort
356	ENDIF
357	c
358	ierr = NF_INQ_VARID (nid, "rlatu", nvarid)
359	IF (ierr .NE. NF_NOERR) THEN
360	PRINT*, "new_grid: Le champ <rlatu> est absent de start.nc"
361	CALL abort
362	ENDIF
363	#ifdef NC_DOUBLE
364	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatuold)
365	#else
366	ierr = NF_GET_VAR_REAL(nid, nvarid, rlatuold)
367	#endif
368	IF (ierr .NE. NF_NOERR) THEN
369	PRINT*, "new_grid: Lecture echouee pour <rlatu>"
370	CALL abort
371	ENDIF
372	c
373	ierr = NF_INQ_VARID (nid, "rlonu", nvarid)
374	IF (ierr .NE. NF_NOERR) THEN
375	PRINT*, "new_grid: Le champ <rlonu> est absent de start.nc"
376	CALL abort
377	ENDIF
378	#ifdef NC_DOUBLE
379	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonuold)
380	#else
381	ierr = NF_GET_VAR_REAL(nid, nvarid, rlonuold)
382	#endif
383	IF (ierr .NE. NF_NOERR) THEN
384	PRINT*, "new_grid: Lecture echouee pour <rlonu>"
385	CALL abort
386	ENDIF
387	c
388	ierr = NF_INQ_VARID (nid, "rlatv", nvarid)
389	IF (ierr .NE. NF_NOERR) THEN
390	PRINT*, "new_grid: Le champ <rlatv> est absent de start.nc"
391	CALL abort
392	ENDIF
393	#ifdef NC_DOUBLE
394	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatvold)
395	#else
396	ierr = NF_GET_VAR_REAL(nid, nvarid, rlatvold)
397	#endif
398	IF (ierr .NE. NF_NOERR) THEN
399	PRINT*, "new_grid: Lecture echouee pour <rlatv>"
400	CALL abort
401	ENDIF
402	c
403
404	print*,"Lecture lon/lat OK"
405
406	c-----------------------------------------------------------------------
407	c Lecture des niveaux verticaux
408	c-----------------------------------------------------------------------
409	c
410	ierr = NF_INQ_VARID (nid, "nivsigs", nvarid)
411	IF (ierr .NE. NF_NOERR) THEN
412	PRINT*, "dynetat0: Le champ <nivsigs> est absent"
413	CALL abort
414	ENDIF
415	#ifdef NC_DOUBLE
416	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, nivsigsold)
417	#else
418	ierr = NF_GET_VAR_REAL(nid, nvarid, nivsigsold)
419	#endif
420	IF (ierr .NE. NF_NOERR) THEN
421	PRINT*, "dynetat0: Lecture echouee pour <nivsigs>"
422	CALL abort
423	ENDIF
424
425	ierr = NF_INQ_VARID (nid, "nivsig", nvarid)
426	IF (ierr .NE. NF_NOERR) THEN
427	PRINT*, "dynetat0: Le champ <nivsig> est absent"
428	CALL abort
429	ENDIF
430	#ifdef NC_DOUBLE
431	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, nivsigold)
432	#else
433	ierr = NF_GET_VAR_REAL(nid, nvarid, nivsigold)
434	#endif
435	IF (ierr .NE. NF_NOERR) THEN
436	PRINT*, "dynetat0: Lecture echouee pour <nivsig>"
437	CALL abort
438	ENDIF
439
440	ierr = NF_INQ_VARID (nid, "ap", nvarid)
441	IF (ierr .NE. NF_NOERR) THEN
442	PRINT*, "new_grid: Le champ <ap> est absent de start.nc"
443	CALL abort
444	ELSE
445	#ifdef NC_DOUBLE
446	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, apold)
447	#else
448	ierr = NF_GET_VAR_REAL(nid, nvarid, apold)
449	#endif
450	IF (ierr .NE. NF_NOERR) THEN
451	PRINT*, "new_grid: Lecture echouee pour <ap>"
452	ENDIF
453	ENDIF
454	c
455	ierr = NF_INQ_VARID (nid, "bp", nvarid)
456	IF (ierr .NE. NF_NOERR) THEN
457	PRINT*, "new_grid: Le champ <bp> est absent de start.nc"
458	CALL abort
459	ENDIF
460	#ifdef NC_DOUBLE
461	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, bpold)
462	#else
463	ierr = NF_GET_VAR_REAL(nid, nvarid, bpold)
464	#endif
465	IF (ierr .NE. NF_NOERR) THEN
466	PRINT*, "new_grid: Lecture echouee pour <bp>"
467	CALL abort
468	END IF
469
470	ierr = NF_INQ_VARID (nid, "presnivs", nvarid)
471	IF (ierr .NE. NF_NOERR) THEN
472	PRINT*, "dynetat0: Le champ <presnivs> est absent"
473	CALL abort
474	ENDIF
475	#ifdef NC_DOUBLE
476	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, presnivsold)
477	#else
478	ierr = NF_GET_VAR_REAL(nid, nvarid, presnivsold)
479	#endif
480	IF (ierr .NE. NF_NOERR) THEN
481	PRINT*, "dynetat0: Lecture echouee pour <presnivs>"
482	CALL abort
483	ENDIF
484
485	c-----------------------------------------------------------------------
486	c Lecture geopotentiel au sol
487	c-----------------------------------------------------------------------
488	c
489	ierr = NF_INQ_VARID (nid, "phisinit", nvarid)
490	IF (ierr .NE. NF_NOERR) THEN
491	PRINT*, "new_grid: Le champ <phisinit> est absent de start.nc"
492	CALL abort
493	ENDIF
494	#ifdef NC_DOUBLE
495	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, phisold)
496	#else
497	ierr = NF_GET_VAR_REAL(nid, nvarid, phisold)
498	#endif
499	IF (ierr .NE. NF_NOERR) THEN
500	PRINT*, "new_grid: Lecture echouee pour <phisinit>"
501	CALL abort
502	ENDIF
503
504	print*,"Lecture niveaux et geopot OK"
505
506	c-----------------------------------------------------------------------
507	c Lecture des champs 2D ()
508	c-----------------------------------------------------------------------
509
510	start=(/1,1,1,0/)
511	counter=(/imold+1,jmold+1,1,0/)
512
513	ierr = NF_INQ_VARID (nid, "ps", nvarid)
514	IF (ierr .NE. NF_NOERR) THEN
515	PRINT*, "lect_start_archive: Le champ <ps> est absent"
516	CALL abort
517	ENDIF
518	#ifdef NC_DOUBLE
519	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,psold)
520	#else
521	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,psold)
522	#endif
523	IF (ierr .NE. NF_NOERR) THEN
524	PRINT*, "lect_start_archive: Lecture echouee pour <ps>"
525	CALL abort
526	ENDIF
527	c
528	ierr = NF_INQ_VARID (nid, "tsurf", nvarid)
529	IF (ierr .NE. NF_NOERR) THEN
530	PRINT*, "lect_start_archive: Le champ <tsurf> est absent"
531	CALL abort
532	ENDIF
533	#ifdef NC_DOUBLE
534	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,tsurfold)
535	#else
536	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,tsurfold)
537	#endif
538	IF (ierr .NE. NF_NOERR) THEN
539	PRINT*, "lect_start_archive: Lecture echouee pour <tsurf>"
540	CALL abort
541	ENDIF
542	c
543	do isoil=1,nsoilmx
544	write(str2,'(i2.2)') isoil
545	c
546	ierr = NF_INQ_VARID (nid, "Tsoil"//str2, nvarid)
547	IF (ierr .NE. NF_NOERR) THEN
548	PRINT*, "lect_start_archive:
549	. Le champ <","Tsoil"//str2,"> est absent"
550	CALL abort
551	ENDIF
552	#ifdef NC_DOUBLE
553	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,
554	. tsoilold(1,1,isoil))
555	#else
556	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,
557	. tsoilold(1,1,isoil))
558	#endif
559	IF (ierr .NE. NF_NOERR) THEN
560	PRINT*, "lect_start_archive:
561	. Lecture echouee pour <","Tsoil"//str2,">"
562	CALL abort
563	ENDIF
564	end do
565	c
566	ierr = NF_INQ_VARID (nid, "albe", nvarid)
567	IF (ierr .NE. NF_NOERR) THEN
568	PRINT*, "lect_start_archive: Le champ <albe> est absent"
569	CALL abort
570	ENDIF
571	#ifdef NC_DOUBLE
572	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,albeold)
573	#else
574	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,albeold)
575	#endif
576	IF (ierr .NE. NF_NOERR) THEN
577	PRINT*, "lect_start_archive: Lecture echouee pour <albe>"
578	CALL abort
579	ENDIF
580	c
581	ierr = NF_INQ_VARID (nid, "radsol", nvarid)
582	IF (ierr .NE. NF_NOERR) THEN
583	PRINT*, "lect_start_archive: Le champ <radsol> est absent"
584	CALL abort
585	ENDIF
586	#ifdef NC_DOUBLE
587	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,radsolold)
588	#else
589	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,radsolold)
590	#endif
591	IF (ierr .NE. NF_NOERR) THEN
592	PRINT*, "lect_start_archive: Lecture echouee pour <radsol>"
593	CALL abort
594	ENDIF
595	c
596	ierr = NF_INQ_VARID (nid, "sollw", nvarid)
597	IF (ierr .NE. NF_NOERR) THEN
598	PRINT*, "lect_start_archive: Le champ <sollw> est absent"
599	CALL abort
600	ENDIF
601	#ifdef NC_DOUBLE
602	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,sollwold)
603	#else
604	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,sollwold)
605	#endif
606	IF (ierr .NE. NF_NOERR) THEN
607	PRINT*, "lect_start_archive: Lecture echouee pour <sollw>"
608	CALL abort
609	ENDIF
610	c
611	ierr = NF_INQ_VARID (nid, "solsw", nvarid)
612	IF (ierr .NE. NF_NOERR) THEN
613	PRINT*, "lect_start_archive: Le champ <solsw> est absent"
614	CALL abort
615	ENDIF
616	#ifdef NC_DOUBLE
617	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,solswold)
618	#else
619	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,solswold)
620	#endif
621	IF (ierr .NE. NF_NOERR) THEN
622	PRINT*, "lect_start_archive: Lecture echouee pour <solsw>"
623	CALL abort
624	ENDIF
625	c
626	ierr = NF_INQ_VARID (nid, "dlw", nvarid)
627	IF (ierr .NE. NF_NOERR) THEN
628	PRINT*, "lect_start_archive: Le champ <dlw> est absent"
629	CALL abort
630	ENDIF
631	#ifdef NC_DOUBLE
632	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,dlwold)
633	#else
634	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,dlwold)
635	#endif
636	IF (ierr .NE. NF_NOERR) THEN
637	PRINT*, "lect_start_archive: Lecture echouee pour <dlw>"
638	CALL abort
639	ENDIF
640	c
641	ierr = NF_INQ_VARID (nid, "zmea", nvarid)
642	IF (ierr .NE. NF_NOERR) THEN
643	PRINT*, "lect_start_archive: Le champ <zmea> est absent"
644	PRINT*, " Il est donc initialise a zero"
645	zmeaold=0.
646	ENDIF
647	#ifdef NC_DOUBLE
648	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zmeaold)
649	#else
650	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zmeaold)
651	#endif
652	IF (ierr .NE. NF_NOERR) THEN
653	PRINT*, "lect_start_archive: Lecture echouee pour <zmea>"
654	CALL abort
655	ENDIF
656	c
657	ierr = NF_INQ_VARID (nid, "zstd", nvarid)
658	IF (ierr .NE. NF_NOERR) THEN
659	PRINT*, "lect_start_archive: Le champ <zstd> est absent"
660	PRINT*, " Il est donc initialise a zero"
661	zstdold=0.
662	ENDIF
663	#ifdef NC_DOUBLE
664	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zstdold)
665	#else
666	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zstdold)
667	#endif
668	IF (ierr .NE. NF_NOERR) THEN
669	PRINT*, "lect_start_archive: Lecture echouee pour <zstd>"
670	CALL abort
671	ENDIF
672	c
673	ierr = NF_INQ_VARID (nid, "zsig", nvarid)
674	IF (ierr .NE. NF_NOERR) THEN
675	PRINT*, "lect_start_archive: Le champ <zsig> est absent"
676	PRINT*, " Il est donc initialise a zero"
677	zsigold=0.
678	ENDIF
679	#ifdef NC_DOUBLE
680	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zsigold)
681	#else
682	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zsigold)
683	#endif
684	IF (ierr .NE. NF_NOERR) THEN
685	PRINT*, "lect_start_archive: Lecture echouee pour <zsig>"
686	CALL abort
687	ENDIF
688	c
689	ierr = NF_INQ_VARID (nid, "zgam", nvarid)
690	IF (ierr .NE. NF_NOERR) THEN
691	PRINT*, "lect_start_archive: Le champ <zgam> est absent"
692	PRINT*, " Il est donc initialise a zero"
693	zgamold=0.
694	ENDIF
695	#ifdef NC_DOUBLE
696	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zgamold)
697	#else
698	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zgamold)
699	#endif
700	IF (ierr .NE. NF_NOERR) THEN
701	PRINT*, "lect_start_archive: Lecture echouee pour <zgam>"
702	CALL abort
703	ENDIF
704	c
705	ierr = NF_INQ_VARID (nid, "zthe", nvarid)
706	IF (ierr .NE. NF_NOERR) THEN
707	PRINT*, "lect_start_archive: Le champ <zthe> est absent"
708	PRINT*, " Il est donc initialise a zero"
709	ztheold=0.
710	ENDIF
711	#ifdef NC_DOUBLE
712	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,ztheold)
713	#else
714	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,ztheold)
715	#endif
716	IF (ierr .NE. NF_NOERR) THEN
717	PRINT*, "lect_start_archive: Lecture echouee pour <zthe>"
718	CALL abort
719	ENDIF
720	c
721	ierr = NF_INQ_VARID (nid, "zpic", nvarid)
722	IF (ierr .NE. NF_NOERR) THEN
723	PRINT*, "lect_start_archive: Le champ <zpic> est absent"
724	PRINT*, " Il est donc initialise a zero"
725	zpicold=0.
726	ENDIF
727	#ifdef NC_DOUBLE
728	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zpicold)
729	#else
730	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zpicold)
731	#endif
732	IF (ierr .NE. NF_NOERR) THEN
733	PRINT*, "lect_start_archive: Lecture echouee pour <zpic>"
734	CALL abort
735	ENDIF
736	c
737	ierr = NF_INQ_VARID (nid, "zval", nvarid)
738	IF (ierr .NE. NF_NOERR) THEN
739	PRINT*, "lect_start_archive: Le champ <zval> est absent"
740	PRINT*, " Il est donc initialise a zero"
741	zvalold=0.
742	ENDIF
743	#ifdef NC_DOUBLE
744	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zvalold)
745	#else
746	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zvalold)
747	#endif
748	IF (ierr .NE. NF_NOERR) THEN
749	PRINT*, "lect_start_archive: Lecture echouee pour <zval>"
750	CALL abort
751	ENDIF
752	c
753
754	print*,"Lecture champs 2D OK"
755
756	c-----------------------------------------------------------------------
757	c Lecture des champs 3D ()
758	c-----------------------------------------------------------------------
759
760	start=(/1,1,1,1/)
761	counter=(/imold+1,jmold+1,lmold,1/)
762	c
763	ierr = NF_INQ_VARID (nid,"temp", nvarid)
764	IF (ierr .NE. NF_NOERR) THEN
765	PRINT*, "lect_start_archive: Le champ <temp> est absent"
766	CALL abort
767	ENDIF
768	#ifdef NC_DOUBLE
769	ierr = NF_GET_VARA_DOUBLE(nid, nvarid, start, counter, told)
770	#else
771	ierr = NF_GET_VARA_REAL(nid, nvarid, start, counter, told)
772	#endif
773	IF (ierr .NE. NF_NOERR) THEN
774	PRINT*, "lect_start_archive: Lecture echouee pour <temp>"
775	CALL abort
776	ENDIF
777	c
778	ierr = NF_INQ_VARID (nid,"u", nvarid)
779	IF (ierr .NE. NF_NOERR) THEN
780	PRINT*, "lect_start_archive: Le champ <u> est absent"
781	CALL abort
782	ENDIF
783	#ifdef NC_DOUBLE
784	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,uold)
785	#else
786	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,uold)
787	#endif
788	IF (ierr .NE. NF_NOERR) THEN
789	PRINT*, "lect_start_archive: Lecture echouee pour <u>"
790	CALL abort
791	ENDIF
792	c
793	ierr = NF_INQ_VARID (nid,"v", nvarid)
794	IF (ierr .NE. NF_NOERR) THEN
795	PRINT*, "lect_start_archive: Le champ <v> est absent"
796	CALL abort
797	ENDIF
798	#ifdef NC_DOUBLE
799	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,vold)
800	#else
801	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,vold)
802	#endif
803	IF (ierr .NE. NF_NOERR) THEN
804	PRINT*, "lect_start_archive: Lecture echouee pour <v>"
805	CALL abort
806	ENDIF
807	c
808
809	c TNAME: IL EST LU A PARTIR DE traceur.def (mettre l'ancien si
810	c changement du nombre de traceurs)
811
812	IF(iflag_trac.eq.1) THEN
813	DO iq=1,nqtot
814	ierr = NF_INQ_VARID (nid, tname(iq), nvarid)
815	IF (ierr .NE. NF_NOERR) THEN
816	PRINT*, "dynetat0: Le champ <"//tname(iq)//"> est absent"
817	PRINT*, " Il est donc initialise a zero"
818	qold(:,:,:,iq)=0.
819	ELSE
820	#ifdef NC_DOUBLE
821	ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,counter,qold(1,1,1,iq))
822	#else
823	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,qold(1,1,1,iq))
824	#endif
825	IF (ierr .NE. NF_NOERR) THEN
826	PRINT*, "dynetat0: Lecture echouee pour "//tname(iq)
827	CALL abort
828	ENDIF
829	ENDIF
830	ENDDO
831	ENDIF
832
833
834	print*,"Lecture champs 3D OK"
835
836	c=======================================================================
837	c INTERPOLATION DANS LA NOUVELLE GRILLE et initialisation des variables
838	c=======================================================================
839	c Interpolation horizontale puis passage dans la grille physique pour
840	c les variables physique
841	c Interpolation verticale puis horizontale pour chaque variable 3D
842	c=======================================================================
843
844	c-----------------------------------------------------------------------
845	c Variable 2d :
846	c-----------------------------------------------------------------------
847
848	c Relief
849	! topoflag = F: we keep the existing topography
850	! = T: we read it from the Relief.nc file
851	! topoflag need to be in the specific run.def for newstart
852	topoflag = . FALSE .
853	CALL getin('topoflag',topoflag)
854
855	IF ( topoflag ) then
856	print*,"Topography (phis) read in file Relief.nc"
857	print*,"For Venus, map directly inverted in Relief.nc"
858	CALL startget_phys2d('surfgeo',iip1,jjp1,rlonv,rlatu,phis,0.0,
859	. jjm ,rlonu,rlatv,.true.)
860	c needed ?
861	phis(iip1,:) = phis(1,:)
862
863	CALL startget_phys1d('zmea',iip1,jjp1,rlonv,rlatu,ngridmx,zmea,
864	. 0.0,jjm,rlonu,rlatv,.true.)
865	CALL startget_phys1d('zstd',iip1,jjp1,rlonv,rlatu,ngridmx,zstd,
866	. 0.0,jjm,rlonu,rlatv,.true.)
867	CALL startget_phys1d('zsig',iip1,jjp1,rlonv,rlatu,ngridmx,zsig,
868	. 0.0,jjm,rlonu,rlatv,.true.)
869	CALL startget_phys1d('zgam',iip1,jjp1,rlonv,rlatu,ngridmx,zgam,
870	. 0.0,jjm,rlonu,rlatv,.true.)
871	CALL startget_phys1d('zthe',iip1,jjp1,rlonv,rlatu,ngridmx,zthe,
872	. 0.0,jjm,rlonu,rlatv,.true.)
873	CALL startget_phys1d('zpic',iip1,jjp1,rlonv,rlatu,ngridmx,zpic,
874	. 0.0,jjm,rlonu,rlatv,.true.)
875	CALL startget_phys1d('zval',iip1,jjp1,rlonv,rlatu,ngridmx,zval,
876	. 0.0,jjm,rlonu,rlatv,.true.)
877
878	ELSE
879	print*,'Using existing topography'
880	call interp_horiz (phisold,phis,imold,jmold,iim,jjm,1,
881	& rlonuold,rlatvold,rlonu,rlatv)
882
883	call interp_horiz (zmeaold,zmeaS,imold,jmold,iim,jjm,1,
884	& rlonuold,rlatvold,rlonu,rlatv)
885	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zmeaS,zmea)
886	call interp_horiz (zstdold,zstdS,imold,jmold,iim,jjm,1,
887	& rlonuold,rlatvold,rlonu,rlatv)
888	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zstdS,zstd)
889	call interp_horiz (zsigold,zsigS,imold,jmold,iim,jjm,1,
890	& rlonuold,rlatvold,rlonu,rlatv)
891	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zsigS,zsig)
892	call interp_horiz (zgamold,zgamS,imold,jmold,iim,jjm,1,
893	& rlonuold,rlatvold,rlonu,rlatv)
894	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zgamS,zgam)
895	call interp_horiz (ztheold,ztheS,imold,jmold,iim,jjm,1,
896	& rlonuold,rlatvold,rlonu,rlatv)
897	call gr_dyn_fi (1,iip1,jjp1,ngridmx,ztheS,zthe)
898	call interp_horiz (zpicold,zpicS,imold,jmold,iim,jjm,1,
899	& rlonuold,rlatvold,rlonu,rlatv)
900	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zpicS,zpic)
901	call interp_horiz (zvalold,zvalS,imold,jmold,iim,jjm,1,
902	& rlonuold,rlatvold,rlonu,rlatv)
903	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zvalS,zval)
904	ENDIF
905
906	print*,"New surface geopotential OK"
907
908	c Lat et lon pour physique
909	do i=1,iip1
910	rlatS(i,:)=rlatu(:)*180./pi
911	enddo
912	call gr_dyn_fi (1,iip1,jjp1,ngridmx,rlatS,rlat)
913
914	do j=2,jjm
915	rlonS(:,j)=rlonv(:)*180./pi
916	enddo
917	rlonS(:,1)=0.
918	rlonS(:,jjp1)=0.
919	call gr_dyn_fi (1,iip1,jjp1,ngridmx,rlonS,rlon)
920
921	c Temperature de surface
922	call interp_horiz (tsurfold,tsurfS,imold,jmold,iim,jjm,1,
923	& rlonuold,rlatvold,rlonu,rlatv)
924	call gr_dyn_fi (1,iip1,jjp1,ngridmx,tsurfS,tsurf)
925	c write(44,*) 'tsurf', tsurf
926
927	c Temperature du sous-sol
928	call interp_horiz(tsoilold,tsoilS,
929	& imold,jmold,iim,jjm,nsoilmx,
930	& rlonuold,rlatvold,rlonu,rlatv)
931	call gr_dyn_fi (nsoilmx,iip1,jjp1,ngridmx,tsoilS,tsoil)
932	c write(45,*) 'tsoil',tsoil
933
934	! CHANGING ALBEDO: may be done through run.def
935	albedoflag = . FALSE .
936	CALL getin('albedoflag',albedoflag)
937
938	IF ( albedoflag ) then
939	print*,"Albedo is reinitialized to the albedo value in run.def"
940	print*,"We may want to consider a map later on..."
941	albedo=0.1
942	CALL getin('albedo',albedo)
943	albe=albedo
944	ELSE
945	call interp_horiz (albeold,albeS,imold,jmold,iim,jjm,1,
946	& rlonuold,rlatvold,rlonu,rlatv)
947	call gr_dyn_fi (1,iip1,jjp1,ngridmx,albeS,albe)
948	ENDIF
949
950	call interp_horiz (radsolold,radsolS,imold,jmold,iim,jjm,1,
951	& rlonuold,rlatvold,rlonu,rlatv)
952	call gr_dyn_fi (1,iip1,jjp1,ngridmx,radsolS,radsol)
953
954	call interp_horiz (sollwold,sollwS,imold,jmold,iim,jjm,1,
955	& rlonuold,rlatvold,rlonu,rlatv)
956	call gr_dyn_fi (1,iip1,jjp1,ngridmx,sollwS,sollw)
957
958	call interp_horiz (solswold,solswS,imold,jmold,iim,jjm,1,
959	& rlonuold,rlatvold,rlonu,rlatv)
960	call gr_dyn_fi (1,iip1,jjp1,ngridmx,solswS,solsw)
961
962	call interp_horiz (dlwold,dlwS,imold,jmold,iim,jjm,1,
963	& rlonuold,rlatvold,rlonu,rlatv)
964	call gr_dyn_fi (1,iip1,jjp1,ngridmx,dlwS,dlw)
965
966	print*,"Nouvelles var physiques OK"
967
968	c-----------------------------------------------------------------------
969	c Traitement special de la pression au sol :
970	c-----------------------------------------------------------------------
971
972	c Extrapolation la pression dans la nouvelle grille
973	call interp_horiz(psold,ps,imold,jmold,iim,jjm,1,
974	& rlonuold,rlatvold,rlonu,rlatv)
975
976	c On assure la conservation de la masse de l'atmosphere
977	c--------------------------------------------------------------
978
979	!!! ATTENTION TEMPORAIRE
980	c ps(:,:)=146700.
981
982	ptotal = 0.
983	DO j=1,jjp1
984	DO i=1,iim
985	ptotal=ptotal+ps(i,j)*aire(i,j)/g
986	END DO
987	END DO
988
989	write(,)
990	write(,)'Ancienne grille: masse de l atm :',ptotalold
991	write(,)'Nouvelle grille: masse de l atm :',ptotal
992	write (,) 'Ratio new atm./ old atm =', ptotal/ptotalold
993	write(,)
994
995
996	DO j=1,jjp1
997	DO i=1,iip1
998	ps(i,j)=ps(i,j) * ptotalold/ptotal
999	END DO
1000	END DO
1001
1002	c la pression de surface et les temperatures ne sont pas reequilibrees en fonction
1003	c de la nouvelle topographie...
1004	c Si l'ajustement inevitable du debut pose des problemes, voir le newstart martien.
1005
1006	print*,"Nouvelle ps OK"
1007	print*,"If changes done on topography, beware !"
1008	print*,"Some time may be needed for adjustments at the beginning"
1009	print*,"so if unstable, relax the timestep and/or dissipation."
1010
1011	c-----------------------------------------------------------------------
1012	c Variable 3d :
1013	c-----------------------------------------------------------------------
1014
1015	CALL pression(ip1jmp1, ap, bp, ps, p3d)
1016	if (disvert_type==1) then
1017	CALL exner_hyb( ip1jmp1, ps, p3d, pks, pk, pkf )
1018	else ! we assume that we are in the disvert_type==2 case
1019	CALL exner_milieu( ip1jmp1, ps, p3d, pks, pk, pkf )
1020	endif
1021
1022	c temperatures atmospheriques
1023
1024	c ATTENTION: peut servir, mais bon...
1025	c do l=1,lmold
1026	c do j=1,jmold+1
1027	c do i=1,imold+1
1028	c modif: profil uniforme
1029	c told(i,j,l)=told(1,jmold/2,l)
1030	c mean T profile:
1031	c told(i,j,l) = 142.1exp(-((p3d(i,j,l)/100.+21.45)/40.11)*2.)
1032	c . + 106.3exp(-((p3d(i,j,l)/100.-3183.)/4737.)*2.)
1033	c enddo
1034	c enddo
1035	c enddo
1036
1037	write (,) 'told ', told (1,jmold+1,1) ! INFO
1038	call interp_vert
1039	& (told,var,lmold,llm,apold,bpold,ap,bp,
1040	& psold,(imold+1)*(jmold+1))
1041	write (,) 'var ', var (1,jmold+1,1) ! INFO
1042	call interp_horiz(var,t,imold,jmold,iim,jjm,llm,
1043	& rlonuold,rlatvold,rlonu,rlatv)
1044	write (,) 'T ', T(1,jjp1,1) ! INFO
1045	! pour info:
1046	! Si extension verticale, la T est extrapolee constante au-dessus de lmold
1047
1048	c passage grille physique pour restartphy.nc
1049	call gr_dyn_fi (llm,iip1,jjp1,ngridmx,T,t_fi)
1050
1051	! ADAPTATION GCM POUR CP(T)
1052	c passage en temperature potentielle
1053	call t2tpot(ip1jmp1*llm,T,teta,pk)
1054	c on assure la periodicite
1055	teta(iip1,:,:) = teta(1,:,:)
1056
1057	! RESETING U TO 0: may be done through run.def
1058	razvitu = . FALSE .
1059	CALL getin('razvitu',razvitu)
1060	razvitv = . FALSE .
1061	CALL getin('razvitv',razvitv)
1062
1063	c calcul des champ de vent; passage en vent covariant
1064	write (,) 'uold ', uold (1,2,1) ! INFO
1065	call interp_vert
1066	& (uold,var,lmold,llm,apold,bpold,ap,bp,
1067	& psold,(imold+1)*(jmold+1))
1068	write (,) 'var ', var (1,2,1) ! INFO
1069	call interp_horiz(var,us,imold,jmold,iim,jjm,llm,
1070	& rlonuold,rlatvold,rlonu,rlatv)
1071	write (,) 'us ', us (1,2,1) ! INFO
1072
1073	call interp_vert
1074	& (vold,var,lmold,llm,
1075	& apold,bpold,ap,bp,psold,(imold+1)*(jmold+1))
1076	call interp_horiz(var,vs,imold,jmold,iim,jjm,llm,
1077	& rlonuold,rlatvold,rlonu,rlatv)
1078	call scal_wind(us,vs,unat,vnat)
1079	! Reseting u=0
1080	if (razvitu) then
1081	unat(:,:,:) = 0.
1082	endif
1083	write (,) 'unat ', unat (1,2,1) ! INFO
1084	do l=1,llm
1085	do j = 1, jjp1
1086	do i=1,iip1
1087	ucov( i,j,l ) = unat( i,j,l ) * cu(i,j)
1088	! pour info:
1089	! Si extension verticale, on impose u=0 au-dessus de lmold
1090	if (l.gt.lmold) ucov( i,j,l ) = 0
1091	end do
1092	end do
1093	end do
1094	write (,) 'ucov ', ucov (1,2,1) ! INFO
1095	c write(48,*) 'ucov',ucov
1096	! Reseting v=0
1097	if (razvitv) then
1098	vnat(:,:,:) = 0.
1099	endif
1100	write (,) 'vnat ', vnat (1,2,1) ! INFO
1101	do l=1,llm
1102	do j = 1, jjm
1103	do i=1,iim
1104	vcov( i,j,l ) = vnat( i,j,l ) * cv(i,j)
1105	! pour info:
1106	! Si extension verticale, on impose v=0 au-dessus de lmold
1107	if (l.gt.lmold) vcov( i,j,l ) = 0
1108	end do
1109	vcov( iip1,j,l ) = vcov( 1,j,l )
1110	end do
1111	end do
1112	c write(49,*) 'ucov',vcov
1113
1114	c masse: on la recalcule (ps a été ajustée pour conserver la masse totale)
1115	call massdair(p3d,masse)
1116
1117	c traceurs 3D
1118	do iq = 1, nqtot
1119	call interp_vert(qold(1,1,1,iq),var,lmold,llm,
1120	& apold,bpold,ap,bp,psold,(imold+1)*(jmold+1))
1121	call interp_horiz(var,q(1,1,1,iq),imold,jmold,iim,jjm,llm,
1122	& rlonuold,rlatvold,rlonu,rlatv)
1123	enddo
1124
1125	print*,"Nouvelles var dynamiques OK"
1126
1127	c=======================================================================
1128	c Ecriture des restart.nc et restartphy.nc :
1129	c=======================================================================
1130
1131	call writerestart('restart.nc',tab_cntrl_dyn,
1132	. phis,vcov,ucov,teta,masse,q,ps)
1133
1134	print*,"restart OK"
1135
1136	call writerestartphy('restartphy.nc',tab_cntrl_fi,ngridmx,llm,
1137	. rlat,rlon,tsurf,tsoil,albe,solsw, sollw,dlw,
1138	. radsol,
1139	. zmea,zstd,zsig,zgam,zthe,zpic,zval,
1140	. t_fi)
1141
1142	print*,"restartphy OK"
1143
1144	end

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

Download in other formats:

Original Format