My Project
 All Classes Files Functions Variables Macros
dynredem_loc.F
Go to the documentation of this file.
1 !
2 ! $Id$
3 !
4 c
5  SUBROUTINE dynredem0_loc(fichnom,iday_end,phis)
6 #ifdef CPP_IOIPSL
7  USE ioipsl
8 #endif
9  USE parallel
10  USE mod_hallo
11  USE infotrac
12  IMPLICIT NONE
13 c=======================================================================
14 c Ecriture du fichier de redemarrage sous format NetCDF (initialisation)
15 c=======================================================================
16 c Declarations:
17 c -------------
18 #include "dimensions.h"
19 #include "paramet.h"
20 #include "comconst.h"
21 #include "comvert.h"
22 #include "comgeom.h"
23 #include "temps.h"
24 #include "ener.h"
25 #include "logic.h"
26 #include "netcdf.inc"
27 #include "description.h"
28 #include "serre.h"
29 #include "iniprint.h"
30 
31 c Arguments:
32 c ----------
33  INTEGER iday_end
34  REAL phis(ijb_u:ije_u)
35  CHARACTER*(*) fichnom
36 
37 c Local:
38 c ------
39  INTEGER iq,l
40  INTEGER length
41  parameter(length = 100)
42  REAL tab_cntrl(length) ! tableau des parametres du run
43  INTEGER ierr
44  character*20 modname
45  character*80 abort_message
46 
47 c Variables locales pour NetCDF:
48 c
49  INTEGER dims2(2), dims3(3), dims4(4)
50  INTEGER idim_index
51  INTEGER idim_rlonu, idim_rlonv, idim_rlatu, idim_rlatv
52  INTEGER idim_s, idim_sig
53  INTEGER idim_tim
54  INTEGER nid,nvarid
55 
56  REAL zan0,zjulian,hours
57  INTEGER yyears0,jjour0, mmois0
58  character*30 unites
59  REAL :: phis_glo(ip1jmp1)
60 
61  CALL gather_field_u(phis,phis_glo,1)
62 
63 
64 c-----------------------------------------------------------------------
65  if (mpi_rank==0) then
66 
67  modname='dynredem0_loc'
68 
69 #ifdef CPP_IOIPSL
70  call ymds2ju(annee_ref, 1, iday_end, 0.0, zjulian)
71  call ju2ymds(zjulian, yyears0, mmois0, jjour0, hours)
72 #else
73 ! set yyears0, mmois0, jjour0 to 0,1,1 (hours is not used)
74  yyears0=0
75  mmois0=1
76  jjour0=1
77 #endif
78 
79  DO l=1,length
80  tab_cntrl(l) = 0.
81  ENDDO
82  tab_cntrl(1) = REAL(iim)
83  tab_cntrl(2) = REAL(jjm)
84  tab_cntrl(3) = REAL(llm)
85  tab_cntrl(4) = REAL(day_ref)
86  tab_cntrl(5) = REAL(annee_ref)
87  tab_cntrl(6) = rad
88  tab_cntrl(7) = omeg
89  tab_cntrl(8) = g
90  tab_cntrl(9) = cpp
91  tab_cntrl(10) = kappa
92  tab_cntrl(11) = daysec
93  tab_cntrl(12) = dtvr
94  tab_cntrl(13) = etot0
95  tab_cntrl(14) = ptot0
96  tab_cntrl(15) = ztot0
97  tab_cntrl(16) = stot0
98  tab_cntrl(17) = ang0
99  tab_cntrl(18) = pa
100  tab_cntrl(19) = preff
101 c
102 c ..... parametres pour le zoom ......
103 
104  tab_cntrl(20) = clon
105  tab_cntrl(21) = clat
106  tab_cntrl(22) = grossismx
107  tab_cntrl(23) = grossismy
108 c
109  IF ( fxyhypb ) THEN
110  tab_cntrl(24) = 1.
111  tab_cntrl(25) = dzoomx
112  tab_cntrl(26) = dzoomy
113  tab_cntrl(27) = 0.
114  tab_cntrl(28) = taux
115  tab_cntrl(29) = tauy
116  ELSE
117  tab_cntrl(24) = 0.
118  tab_cntrl(25) = dzoomx
119  tab_cntrl(26) = dzoomy
120  tab_cntrl(27) = 0.
121  tab_cntrl(28) = 0.
122  tab_cntrl(29) = 0.
123  IF( ysinus ) tab_cntrl(27) = 1.
124  ENDIF
125 
126  tab_cntrl(30) = REAL(iday_end)
127  tab_cntrl(31) = REAL(itau_dyn + itaufin)
128 c start_time: start_time of simulation (not necessarily 0.)
129  tab_cntrl(32) = start_time
130 c
131 c .........................................................
132 c
133 c Creation du fichier:
134 c
135  ierr = nf_create(fichnom, nf_clobber, nid)
136  IF (ierr.NE.nf_noerr) THEN
137  write(lunout,*)"dynredem0: Pb d ouverture du fichier "
138  & //trim(fichnom)
139  write(lunout,*)' ierr = ', ierr
140  CALL abort
141  ENDIF
142 c
143 c Preciser quelques attributs globaux:
144 c
145  ierr = nf_put_att_text(nid, nf_global, "title", 27,
146  . "Fichier demmarage dynamique")
147 c
148 c Definir les dimensions du fichiers:
149 c
150  ierr = nf_def_dim(nid, "index", length, idim_index)
151  ierr = nf_def_dim(nid, "rlonu", iip1, idim_rlonu)
152  ierr = nf_def_dim(nid, "rlatu", jjp1, idim_rlatu)
153  ierr = nf_def_dim(nid, "rlonv", iip1, idim_rlonv)
154  ierr = nf_def_dim(nid, "rlatv", jjm, idim_rlatv)
155  ierr = nf_def_dim(nid, "sigs", llm, idim_s)
156  ierr = nf_def_dim(nid, "sig", llmp1, idim_sig)
157  ierr = nf_def_dim(nid, "temps", nf_unlimited, idim_tim)
158 c
159  ierr = nf_enddef(nid) ! sortir du mode de definition
160 c
161 c Definir et enregistrer certains champs invariants:
162 c
163  ierr = nf_redef(nid)
164 cIM 220306 BEG
165 #ifdef NC_DOUBLE
166  ierr = nf_def_var(nid,"controle",nf_double,1,idim_index,nvarid)
167 #else
168  ierr = nf_def_var(nid,"controle",nf_float,1,idim_index,nvarid)
169 #endif
170 cIM 220306 END
171  ierr = nf_put_att_text(nid, nvarid, "title", 22,
172  . "Parametres de controle")
173  ierr = nf_enddef(nid)
174 #ifdef NC_DOUBLE
175  ierr = nf_put_var_double(nid,nvarid,tab_cntrl)
176 #else
177  ierr = nf_put_var_real(nid,nvarid,tab_cntrl)
178 #endif
179 c
180  ierr = nf_redef(nid)
181 cIM 220306 BEG
182 #ifdef NC_DOUBLE
183  ierr = nf_def_var(nid,"rlonu",nf_double,1,idim_rlonu,nvarid)
184 #else
185  ierr = nf_def_var(nid,"rlonu",nf_float,1,idim_rlonu,nvarid)
186 #endif
187 cIM 220306 END
188  ierr = nf_put_att_text(nid, nvarid, "title", 23,
189  . "Longitudes des points U")
190  ierr = nf_enddef(nid)
191 #ifdef NC_DOUBLE
192  ierr = nf_put_var_double(nid,nvarid,rlonu)
193 #else
194  ierr = nf_put_var_real(nid,nvarid,rlonu)
195 #endif
196 c
197  ierr = nf_redef(nid)
198 cIM 220306 BEG
199 #ifdef NC_DOUBLE
200  ierr = nf_def_var(nid,"rlatu",nf_double,1,idim_rlatu,nvarid)
201 #else
202  ierr = nf_def_var(nid,"rlatu",nf_float,1,idim_rlatu,nvarid)
203 #endif
204 cIM 220306 END
205  ierr = nf_put_att_text(nid, nvarid, "title", 22,
206  . "Latitudes des points U")
207  ierr = nf_enddef(nid)
208 #ifdef NC_DOUBLE
209  ierr = nf_put_var_double(nid,nvarid,rlatu)
210 #else
211  ierr = nf_put_var_real(nid,nvarid,rlatu)
212 #endif
213 c
214  ierr = nf_redef(nid)
215 cIM 220306 BEG
216 #ifdef NC_DOUBLE
217  ierr = nf_def_var(nid,"rlonv",nf_double,1,idim_rlonv,nvarid)
218 #else
219  ierr = nf_def_var(nid,"rlonv",nf_float,1,idim_rlonv,nvarid)
220 #endif
221 cIM 220306 END
222  ierr = nf_put_att_text(nid, nvarid, "title", 23,
223  . "Longitudes des points V")
224  ierr = nf_enddef(nid)
225 #ifdef NC_DOUBLE
226  ierr = nf_put_var_double(nid,nvarid,rlonv)
227 #else
228  ierr = nf_put_var_real(nid,nvarid,rlonv)
229 #endif
230 c
231  ierr = nf_redef(nid)
232 cIM 220306 BEG
233 #ifdef NC_DOUBLE
234  ierr = nf_def_var(nid,"rlatv",nf_double,1,idim_rlatv,nvarid)
235 #else
236  ierr = nf_def_var(nid,"rlatv",nf_float,1,idim_rlatv,nvarid)
237 #endif
238 cIM 220306 END
239  ierr = nf_put_att_text(nid, nvarid, "title", 22,
240  . "Latitudes des points V")
241  ierr = nf_enddef(nid)
242 #ifdef NC_DOUBLE
243  ierr = nf_put_var_double(nid,nvarid,rlatv)
244 #else
245  ierr = nf_put_var_real(nid,nvarid,rlatv)
246 #endif
247 c
248  ierr = nf_redef(nid)
249 cIM 220306 BEG
250 #ifdef NC_DOUBLE
251  ierr = nf_def_var(nid,"nivsigs",nf_double,1,idim_s,nvarid)
252 #else
253  ierr = nf_def_var(nid,"nivsigs",nf_float,1,idim_s,nvarid)
254 #endif
255 cIM 220306 END
256  ierr = nf_put_att_text(nid, nvarid, "title", 28,
257  . "Numero naturel des couches s")
258  ierr = nf_enddef(nid)
259 #ifdef NC_DOUBLE
260  ierr = nf_put_var_double(nid,nvarid,nivsigs)
261 #else
262  ierr = nf_put_var_real(nid,nvarid,nivsigs)
263 #endif
264 c
265  ierr = nf_redef(nid)
266 cIM 220306 BEG
267 #ifdef NC_DOUBLE
268  ierr = nf_def_var(nid,"nivsig",nf_double,1,idim_sig,nvarid)
269 #else
270  ierr = nf_def_var(nid,"nivsig",nf_float,1,idim_sig,nvarid)
271 #endif
272 cIM 220306 END
273  ierr = nf_put_att_text(nid, nvarid, "title", 32,
274  . "Numero naturel des couches sigma")
275  ierr = nf_enddef(nid)
276 #ifdef NC_DOUBLE
277  ierr = nf_put_var_double(nid,nvarid,nivsig)
278 #else
279  ierr = nf_put_var_real(nid,nvarid,nivsig)
280 #endif
281 c
282  ierr = nf_redef(nid)
283 cIM 220306 BEG
284 #ifdef NC_DOUBLE
285  ierr = nf_def_var(nid,"ap",nf_double,1,idim_sig,nvarid)
286 #else
287  ierr = nf_def_var(nid,"ap",nf_float,1,idim_sig,nvarid)
288 #endif
289 cIM 220306 END
290  ierr = nf_put_att_text(nid, nvarid, "title", 26,
291  . "Coefficient A pour hybride")
292  ierr = nf_enddef(nid)
293 #ifdef NC_DOUBLE
294  ierr = nf_put_var_double(nid,nvarid,ap)
295 #else
296  ierr = nf_put_var_real(nid,nvarid,ap)
297 #endif
298 c
299  ierr = nf_redef(nid)
300 cIM 220306 BEG
301 #ifdef NC_DOUBLE
302  ierr = nf_def_var(nid,"bp",nf_double,1,idim_sig,nvarid)
303 #else
304  ierr = nf_def_var(nid,"bp",nf_float,1,idim_sig,nvarid)
305 #endif
306 cIM 220306 END
307  ierr = nf_put_att_text(nid, nvarid, "title", 26,
308  . "Coefficient B pour hybride")
309  ierr = nf_enddef(nid)
310 #ifdef NC_DOUBLE
311  ierr = nf_put_var_double(nid,nvarid,bp)
312 #else
313  ierr = nf_put_var_real(nid,nvarid,bp)
314 #endif
315 c
316  ierr = nf_redef(nid)
317 cIM 220306 BEG
318 #ifdef NC_DOUBLE
319  ierr = nf_def_var(nid,"presnivs",nf_double,1,idim_s,nvarid)
320 #else
321  ierr = nf_def_var(nid,"presnivs",nf_float,1,idim_s,nvarid)
322 #endif
323 cIM 220306 END
324  ierr = nf_enddef(nid)
325 #ifdef NC_DOUBLE
326  ierr = nf_put_var_double(nid,nvarid,presnivs)
327 #else
328  ierr = nf_put_var_real(nid,nvarid,presnivs)
329 #endif
330 c
331 c Coefficients de passage cov. <-> contra. <--> naturel
332 c
333  ierr = nf_redef(nid)
334  dims2(1) = idim_rlonu
335  dims2(2) = idim_rlatu
336 cIM 220306 BEG
337 #ifdef NC_DOUBLE
338  ierr = nf_def_var(nid,"cu",nf_double,2,dims2,nvarid)
339 #else
340  ierr = nf_def_var(nid,"cu",nf_float,2,dims2,nvarid)
341 #endif
342 cIM 220306 END
343  ierr = nf_put_att_text(nid, nvarid, "title", 29,
344  . "Coefficient de passage pour U")
345  ierr = nf_enddef(nid)
346 #ifdef NC_DOUBLE
347  ierr = nf_put_var_double(nid,nvarid,cu)
348 #else
349  ierr = nf_put_var_real(nid,nvarid,cu)
350 #endif
351 c
352  ierr = nf_redef(nid)
353  dims2(1) = idim_rlonv
354  dims2(2) = idim_rlatv
355 cIM 220306 BEG
356 #ifdef NC_DOUBLE
357  ierr = nf_def_var(nid,"cv",nf_double,2,dims2,nvarid)
358 #else
359  ierr = nf_def_var(nid,"cv",nf_float,2,dims2,nvarid)
360 #endif
361 cIM 220306 END
362  ierr = nf_put_att_text(nid, nvarid, "title", 29,
363  . "Coefficient de passage pour V")
364  ierr = nf_enddef(nid)
365 #ifdef NC_DOUBLE
366  ierr = nf_put_var_double(nid,nvarid,cv)
367 #else
368  ierr = nf_put_var_real(nid,nvarid,cv)
369 #endif
370 c
371 c Aire de chaque maille:
372 c
373  ierr = nf_redef(nid)
374  dims2(1) = idim_rlonv
375  dims2(2) = idim_rlatu
376 cIM 220306 BEG
377 #ifdef NC_DOUBLE
378  ierr = nf_def_var(nid,"aire",nf_double,2,dims2,nvarid)
379 #else
380  ierr = nf_def_var(nid,"aire",nf_float,2,dims2,nvarid)
381 #endif
382 cIM 220306 END
383  ierr = nf_put_att_text(nid, nvarid, "title", 22,
384  . "Aires de chaque maille")
385  ierr = nf_enddef(nid)
386 #ifdef NC_DOUBLE
387  ierr = nf_put_var_double(nid,nvarid,aire)
388 #else
389  ierr = nf_put_var_real(nid,nvarid,aire)
390 #endif
391 c
392 c Geopentiel au sol:
393 c
394  ierr = nf_redef(nid)
395  dims2(1) = idim_rlonv
396  dims2(2) = idim_rlatu
397 cIM 220306 BEG
398 #ifdef NC_DOUBLE
399  ierr = nf_def_var(nid,"phisinit",nf_double,2,dims2,nvarid)
400 #else
401  ierr = nf_def_var(nid,"phisinit",nf_float,2,dims2,nvarid)
402 #endif
403 cIM 220306 END
404  ierr = nf_put_att_text(nid, nvarid, "title", 19,
405  . "Geopotentiel au sol")
406  ierr = nf_enddef(nid)
407 #ifdef NC_DOUBLE
408  ierr = nf_put_var_double(nid,nvarid,phis_glo)
409 #else
410  ierr = nf_put_var_real(nid,nvarid,phis_glo)
411 #endif
412 c
413 c Definir les variables pour pouvoir les enregistrer plus tard:
414 c
415  ierr = nf_redef(nid) ! entrer dans le mode de definition
416 c
417 cIM 220306 BEG
418 #ifdef NC_DOUBLE
419  ierr = nf_def_var(nid,"temps",nf_double,1,idim_tim,nvarid)
420 #else
421  ierr = nf_def_var(nid,"temps",nf_float,1,idim_tim,nvarid)
422 #endif
423 cIM 220306 END
424  ierr = nf_put_att_text(nid, nvarid, "title", 19,
425  . "Temps de simulation")
426  write(unites,200)yyears0,mmois0,jjour0
427 200 format('days since ',i4,'-',i2.2,'-',i2.2,' 00:00:00')
428  ierr = nf_put_att_text(nid, nvarid, "units", 30,
429  . unites)
430 
431 c
432  dims4(1) = idim_rlonu
433  dims4(2) = idim_rlatu
434  dims4(3) = idim_s
435  dims4(4) = idim_tim
436 cIM 220306 BEG
437 #ifdef NC_DOUBLE
438  ierr = nf_def_var(nid,"ucov",nf_double,4,dims4,nvarid)
439 #else
440  ierr = nf_def_var(nid,"ucov",nf_float,4,dims4,nvarid)
441 #endif
442 cIM 220306 END
443  ierr = nf_put_att_text(nid, nvarid, "title", 9,
444  . "Vitesse U")
445 c
446  dims4(1) = idim_rlonv
447  dims4(2) = idim_rlatv
448  dims4(3) = idim_s
449  dims4(4) = idim_tim
450 cIM 220306 BEG
451 #ifdef NC_DOUBLE
452  ierr = nf_def_var(nid,"vcov",nf_double,4,dims4,nvarid)
453 #else
454  ierr = nf_def_var(nid,"vcov",nf_float,4,dims4,nvarid)
455 #endif
456 cIM 220306 END
457  ierr = nf_put_att_text(nid, nvarid, "title", 9,
458  . "Vitesse V")
459 c
460  dims4(1) = idim_rlonv
461  dims4(2) = idim_rlatu
462  dims4(3) = idim_s
463  dims4(4) = idim_tim
464 cIM 220306 BEG
465 #ifdef NC_DOUBLE
466  ierr = nf_def_var(nid,"teta",nf_double,4,dims4,nvarid)
467 #else
468  ierr = nf_def_var(nid,"teta",nf_float,4,dims4,nvarid)
469 #endif
470 cIM 220306 END
471  ierr = nf_put_att_text(nid, nvarid, "title", 11,
472  . "Temperature")
473 c
474  dims4(1) = idim_rlonv
475  dims4(2) = idim_rlatu
476  dims4(3) = idim_s
477  dims4(4) = idim_tim
478 
479  DO iq=1,nqtot
480 cIM 220306 BEG
481 #ifdef NC_DOUBLE
482  ierr = nf_def_var(nid,tname(iq),nf_double,4,dims4,nvarid)
483 #else
484  ierr = nf_def_var(nid,tname(iq),nf_float,4,dims4,nvarid)
485 #endif
486 cIM 220306 END
487  ierr = nf_put_att_text(nid, nvarid, "title", 12,ttext(iq))
488  ENDDO
489 c
490  dims4(1) = idim_rlonv
491  dims4(2) = idim_rlatu
492  dims4(3) = idim_s
493  dims4(4) = idim_tim
494 cIM 220306 BEG
495 #ifdef NC_DOUBLE
496  ierr = nf_def_var(nid,"masse",nf_double,4,dims4,nvarid)
497 #else
498  ierr = nf_def_var(nid,"masse",nf_float,4,dims4,nvarid)
499 #endif
500 cIM 220306 END
501  ierr = nf_put_att_text(nid, nvarid, "title", 12,
502  . "C est quoi ?")
503 c
504  dims3(1) = idim_rlonv
505  dims3(2) = idim_rlatu
506  dims3(3) = idim_tim
507 cIM 220306 BEG
508 #ifdef NC_DOUBLE
509  ierr = nf_def_var(nid,"ps",nf_double,3,dims3,nvarid)
510 #else
511  ierr = nf_def_var(nid,"ps",nf_float,3,dims3,nvarid)
512 #endif
513 cIM 220306 END
514  ierr = nf_put_att_text(nid, nvarid, "title", 15,
515  . "Pression au sol")
516 c
517  ierr = nf_enddef(nid) ! sortir du mode de definition
518  ierr = nf_close(nid) ! fermer le fichier
519 
520  write(lunout,*)'dynredem_loc: iim,jjm,llm,iday_end',
521  & iim,jjm,llm,iday_end
522  write(lunout,*)'dynredem_loc: rad,omeg,g,cpp,kappa',
523  & rad,omeg,g,cpp,kappa
524 
525  endif ! mpi_rank==0
526  RETURN
527  END
528  SUBROUTINE dynredem1_loc(fichnom,time,
529  . vcov,ucov,teta,q,masse,ps)
530  USE parallel
531  USE mod_hallo
532  USE infotrac
533  USE control_mod
534  USE dynredem_mod
535  IMPLICIT NONE
536 c=================================================================
537 c Ecriture du fichier de redemarrage sous format NetCDF
538 c=================================================================
539 #include "dimensions.h"
540 #include "paramet.h"
541 #include "description.h"
542 #include "netcdf.inc"
543 #include "comvert.h"
544 #include "comgeom.h"
545 #include "temps.h"
546 #include "iniprint.h"
547 
548  INTEGER l
549  REAL vcov(ijb_v:ije_v,llm),ucov(ijb_u:ije_u,llm)
550  REAL teta(ijb_u:ije_u,llm)
551  REAL ps(ijb_u:ije_u),masse(ijb_u:ije_u,llm)
552  REAL q(ijb_u:ije_u,llm,nqtot)
553  CHARACTER*(*) fichnom
554 
555  REAL time
556  INTEGER nid, nvarid, nid_trac, nvarid_trac
557  REAL trac_tmp(ijb_u:ije_u,llm)
558  INTEGER ierr, ierr_file
559  INTEGER iq
560  INTEGER length
561  parameter(length = 100)
562  REAL tab_cntrl(length) ! tableau des parametres du run
563  character*20 modname
564  character*80 abort_message
565 c
566  INTEGER nb
567  SAVE nb
568  DATA nb / 0 /
569  REAL,SAVE,ALLOCATABLE :: ucov_glo(:,:),vcov_glo(:,:),teta_glo(:,:)
570  REAL,SAVE,ALLOCATABLE :: masse_glo(:,:),ps_glo(:),q_glo(:,:)
571  LOGICAL,SAVE :: exist_file
572  INTEGER,SAVE :: ierr_var
573 
574 ! call Gather_Field(ucov,ip1jmp1,llm,0)
575 ! call Gather_Field(vcov,ip1jm,llm,0)
576 ! call Gather_Field(teta,ip1jmp1,llm,0)
577 ! call Gather_Field(masse,ip1jmp1,llm,0)
578 ! call Gather_Field(ps,ip1jmp1,1,0)
579 
580 ! do iq=1,nqtot
581 ! call Gather_Field(q(1,1,iq),ip1jmp1,llm,0)
582  ! enddo
583 
584 !$OMP MASTER
585  if (mpi_rank==0) then
586  modname = 'dynredem1_loc'
587  ierr = nf_open(fichnom, nf_write, nid)
588  IF (ierr .NE. nf_noerr) THEN
589  write(lunout,*)"dynredem1: Pb. d ouverture "//trim(fichnom)
590  CALL abort
591  ENDIF
592 
593 c Ecriture/extension de la coordonnee temps
594 
595  nb = nb + 1
596  ierr = nf_inq_varid(nid, "temps", nvarid)
597  IF (ierr .NE. nf_noerr) THEN
598  write(lunout,*) nf_strerror(ierr)
599  abort_message='Variable temps n est pas definie'
600  CALL abort_gcm(modname,abort_message,ierr)
601  ENDIF
602 #ifdef NC_DOUBLE
603  ierr = nf_put_var1_double(nid,nvarid,nb,time)
604 #else
605  ierr = nf_put_var1_real(nid,nvarid,nb,time)
606 #endif
607  write(lunout,*) "dynredem1_loc: Enregistrement pour ", nb, time
608 
609 c
610 c Re-ecriture du tableau de controle, itaufin n'est plus defini quand
611 c on passe dans dynredem0
612  ierr = nf_inq_varid(nid, "controle", nvarid)
613  IF (ierr .NE. nf_noerr) THEN
614  abort_message="dynredem1: Le champ <controle> est absent"
615  ierr = 1
616  CALL abort_gcm(modname,abort_message,ierr)
617  ENDIF
618 #ifdef NC_DOUBLE
619  ierr = nf_get_var_double(nid, nvarid, tab_cntrl)
620 #else
621  ierr = nf_get_var_real(nid, nvarid, tab_cntrl)
622 #endif
623  tab_cntrl(31) = REAL(itau_dyn + itaufin)
624 #ifdef NC_DOUBLE
625  ierr = nf_put_var_double(nid,nvarid,tab_cntrl)
626 #else
627  ierr = nf_put_var_real(nid,nvarid,tab_cntrl)
628 #endif
629  endif
630 !$OMP END MASTER
631 
632 !
633  CALL dynredem_write_u(nid,"ucov",ucov,llm)
634  CALL dynredem_write_v(nid,"vcov",vcov,llm)
635  CALL dynredem_write_u(nid,"teta",teta,llm)
636  CALL dynredem_write_u(nid,"masse",masse,llm)
637  CALL dynredem_write_u(nid,"ps",ps,1)
638 
639  IF (type_trac /= 'inca') THEN
640  DO iq=1,nqtot
641  CALL dynredem_write_u(nid,tname(iq),q(:,:,iq),llm)
642  ENDDO
643  ELSE
644 
645 !$OMP MASTER
646  INQUIRE(file="start_trac.nc", exist=exist_file)
647  print *, "EXIST", exist_file
648 !$OMP END MASTER
649 !$OMP BARRIER
650 
651  IF (exist_file) THEN
652 !$OMP MASTER
653  ierr_file = nf_open("start_trac.nc", nf_nowrite,nid_trac)
654  IF (ierr_file .NE.nf_noerr) THEN
655  WRITE(6,*)' Pb d''ouverture du fichier start_trac.nc'
656  WRITE(6,*)' ierr = ', ierr_file
657  ENDIF
658 !$OMP END MASTER
659 
660  DO iq=1,nqtot
661 
662 !$OMP MASTER
663  ierr_var = nf_inq_varid(nid_trac, tname(iq), nvarid_trac)
664 !$OMP END MASTER
665 !$OMP BARRIER
666  IF (ierr == nf_noerr) THEN
667  CALL dynredem_read_u(nid_trac,tname(iq),q(:,:,iq),llm)
668  ENDIF
669  CALL dynredem_write_u(nid,tname(iq),q(:,:,iq),llm)
670  ENDDO
671 
672  ELSE ! pas de fichier start_tract
673  DO iq=1,nqtot
674  CALL dynredem_write_u(nid,tname(iq),q(:,:,iq),llm)
675  ENDDO
676  ENDIF
677  ENDIF
678 
679 
680 !$OMP MASTER
681  IF (mpi_rank==0) THEN
682  ierr = nf_close(nid)
683  ENDIF ! mpi_rank==0
684 !$OMP END MASTER
685 
686  RETURN
687  END
688