source: trunk/LMDZ.MARS/libf/dynlonlat_phylonlat/phymars/start2archive.F @ 1403

Last change on this file since 1403 was 1403, checked in by emillour, 10 years ago

All models: Reorganizing the physics/dynamics interface.

  • makelmdz and makelmdz_fcm scripts adapted to handle the new directory settings
  • misc: (replaces what was the "bibio" directory)
  • Should only contain extremely generic (and non physics or dynamics-specific) routines
  • Therefore moved initdynav.F90, initfluxsto.F, inithist.F, writedynav.F90, write_field.F90, writehist.F to "dyn3d_common"
  • dynlonlat_phylonlat: (new interface directory)
  • This directory contains routines relevent to physics/dynamics grid interactions, e.g. routines gr_dyn_fi or gr_fi_dyn and calfis
  • Moreover the dynlonlat_phylonlat contains directories "phy*" corresponding to each physics package "phy*" to be used. These subdirectories should only contain specific interfaces (e.g. iniphysiq) or main programs (e.g. newstart)
  • phy*/dyn1d: this subdirectory contains the 1D model using physics from phy*

EM
File size: 14.8 KB
Line 
1c=======================================================================
2      PROGRAM start2archive
3c=======================================================================
4c
5c
6c   Date:    01/1997
7c   ----
8c
9c
10c   Objet:   Passage des  fichiers netcdf d'etat initial "start" et
11c   -----    "startfi" a un fichier netcdf unique "start_archive"
12c
13c  "start_archive" est une banque d'etats initiaux:
14c  On peut stocker plusieurs etats initiaux dans un meme fichier "start_archive"
15c    (Veiller dans ce cas avoir un day_ini different pour chacun des start)
16c
17c
18c
19c=======================================================================
20
21      use infotrac, only: iniadvtrac, nqtot, tname
22      use comsoil_h, only: nsoilmx, inertiedat
23      use surfdat_h, only: ini_surfdat_h, qsurf
24      use comsoil_h, only: ini_comsoil_h
25      use comgeomphy, only: initcomgeomphy
26      use filtreg_mod, only: inifilr
27      implicit none
28
29#include "dimensions.h"
30      integer, parameter :: ngridmx = (2+(jjm-1)*iim - 1/jjm)
31#include "paramet.h"
32#include "comconst.h"
33#include "comdissip.h"
34#include "comvert.h"
35#include "comgeom.h"
36#include "logic.h"
37#include "temps.h"
38#include "ener.h"
39#include "description.h"
40#include "netcdf.inc"
41
42c-----------------------------------------------------------------------
43c   Declarations
44c-----------------------------------------------------------------------
45
46c variables dynamiques du GCM
47c -----------------------------
48      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
49      REAL teta(ip1jmp1,llm)                    ! temperature potentielle
50      REAL,ALLOCATABLE :: q(:,:,:)   ! champs advectes
51      REAL pks(ip1jmp1)                      ! exner (f pour filtre)
52      REAL pk(ip1jmp1,llm)
53      REAL pkf(ip1jmp1,llm)
54      REAL beta(iip1,jjp1,llm)
55      REAL phis(ip1jmp1)                     ! geopotentiel au sol
56      REAL masse(ip1jmp1,llm)                ! masse de l'atmosphere
57      REAL ps(ip1jmp1)                       ! pression au sol
58      REAL p3d(iip1, jjp1, llm+1)            ! pression aux interfaces
59     
60c Variable Physiques (grille physique)
61c ------------------------------------
62      REAL tsurf(ngridmx)       ! Surface temperature
63      REAL tsoil(ngridmx,nsoilmx) ! Soil temperature
64      REAL co2ice(ngridmx)      ! CO2 ice layer
65      REAL tauscaling(ngridmx) ! dust conversion factor
66      REAL q2(ngridmx,llm+1)
67      REAL emis(ngridmx)
68      INTEGER start,length
69      PARAMETER (length = 100)
70      REAL tab_cntrl_fi(length) ! tableau des parametres de startfi
71      INTEGER*4 day_ini_fi
72
73c Variable naturelle / grille scalaire
74c ------------------------------------
75      REAL T(ip1jmp1,llm),us(ip1jmp1,llm),vs(ip1jmp1,llm)
76      REAL tsurfS(ip1jmp1)
77      REAL tsoilS(ip1jmp1,nsoilmx)
78      REAL ithS(ip1jmp1,nsoilmx) ! Soil Thermal Inertia
79      REAL co2iceS(ip1jmp1)
80      REAL tauscalingS(ip1jmp1)
81      REAL q2S(ip1jmp1,llm+1)
82      REAL,ALLOCATABLE :: qsurfS(:,:)
83      REAL emisS(ip1jmp1)
84
85c Variables intermediaires : vent naturel, mais pas coord scalaire
86c----------------------------------------------------------------
87      REAL vn(ip1jm,llm),un(ip1jmp1,llm)
88
89c Autres  variables
90c -----------------
91      LOGICAL startdrs
92      INTEGER Lmodif
93
94      REAL ptotal, co2icetotal
95      REAL timedyn,timefi !fraction du jour dans start, startfi
96      REAL date
97
98      CHARACTER*2 str2
99      CHARACTER*80 fichier
100      data  fichier /'startfi'/
101
102      INTEGER ij, l,i,j,isoil,iq
103      character*80      fichnom
104      integer :: ierr,ntime
105      integer :: nq,numvanle
106      character(len=30) :: txt ! to store some text
107
108c Netcdf
109c-------
110      integer varid,dimid,timelen
111      INTEGER nid,nid1
112
113c-----------------------------------------------------------------------
114c   Initialisations
115c-----------------------------------------------------------------------
116
117      CALL defrun_new(99, .TRUE. )
118      grireg   = .TRUE.
119! initialize "serial/parallel" related stuff
120      CALL init_phys_lmdz(iim,jjp1,llm,1,(/(jjm-1)*iim+2/))
121      call initcomgeomphy
122
123c=======================================================================
124c Lecture des donnees
125c=======================================================================
126! Load tracer number and names:
127      call iniadvtrac(nqtot,numvanle)
128
129! allocate arrays:
130      allocate(q(ip1jmp1,llm,nqtot))
131      allocate(qsurfS(ip1jmp1,nqtot))
132      call ini_surfdat_h(ngridmx,nqtot)
133      call ini_comsoil_h(ngridmx)
134     
135
136      fichnom = 'start.nc'
137      CALL dynetat0(fichnom,nqtot,vcov,ucov,teta,q,masse,
138     .       ps,phis,timedyn)
139
140
141      fichnom = 'startfi.nc'
142      Lmodif=0
143
144      CALL phyetat0 (fichnom,0,Lmodif,nsoilmx,ngridmx,llm,nqtot,
145     &      day_ini_fi,timefi,tsurf,tsoil,emis,q2,qsurf,co2ice,
146     &      tauscaling)
147
148       ierr = NF_OPEN (fichnom, NF_NOWRITE,nid1)
149       IF (ierr.NE.NF_NOERR) THEN
150         write(6,*)' Pb d''ouverture du fichier'//fichnom
151        CALL ABORT
152       ENDIF
153                                               
154      ierr = NF_INQ_VARID (nid1, "controle", varid)
155      IF (ierr .NE. NF_NOERR) THEN
156       PRINT*, "start2archive: Le champ <controle> est absent"
157       CALL abort
158      ENDIF
159#ifdef NC_DOUBLE
160       ierr = NF_GET_VAR_DOUBLE(nid1, varid, tab_cntrl_fi)
161#else
162      ierr = NF_GET_VAR_REAL(nid1, varid, tab_cntrl_fi)
163#endif
164       IF (ierr .NE. NF_NOERR) THEN
165          PRINT*, "start2archive: Lecture echoue pour <controle>"
166          CALL abort
167       ENDIF
168
169      ierr = NF_CLOSE(nid1)
170
171c-----------------------------------------------------------------------
172c Controle de la synchro
173c-----------------------------------------------------------------------
174!mars a voir      if ((day_ini_fi.ne.day_ini).or.(abs(timefi-timedyn).gt.1.e-10))
175      if ((day_ini_fi.ne.day_ini))
176     &  stop ' Probleme de Synchro entre start et startfi !!!'
177
178
179c *****************************************************************
180c    Option : Reinitialisation des dates dans la premieres annees :
181       do while (day_ini.ge.669)
182          day_ini=day_ini-669
183       enddo
184c *****************************************************************
185
186c-----------------------------------------------------------------------
187c   Initialisations
188c-----------------------------------------------------------------------
189
190      CALL defrun_new(99, .FALSE. )
191      call iniconst
192      call inigeom
193      call inifilr
194      CALL pression(ip1jmp1, ap, bp, ps, p3d)
195      call exner_hyb(ip1jmp1, ps, p3d, beta, pks, pk, pkf)
196
197c=======================================================================
198c Transformation EN VARIABLE NATURELLE / GRILLE SCALAIRE si necessaire
199c=======================================================================
200c  Les variables modeles dependent de la resolution. Il faut donc
201c  eliminer les facteurs responsables de cette dependance
202c  (pour utiliser newstart)
203c=======================================================================
204
205c-----------------------------------------------------------------------
206c Vent   (depend de la resolution horizontale)
207c-----------------------------------------------------------------------
208c
209c ucov --> un  et  vcov --> vn
210c un --> us  et   vn --> vs
211c
212c-----------------------------------------------------------------------
213
214      call covnat(llm,ucov, vcov, un, vn)
215      call wind_scal(un,vn,us,vs)
216
217c-----------------------------------------------------------------------
218c Temperature  (depend de la resolution verticale => de "sigma.def")
219c-----------------------------------------------------------------------
220c
221c h --> T
222c
223c-----------------------------------------------------------------------
224
225      DO l=1,llm
226         DO ij=1,ip1jmp1
227            T(ij,l)=teta(ij,l)*pk(ij,l)/cpp !mars deduit de l'equation dans newstart
228         ENDDO
229      ENDDO
230
231c-----------------------------------------------------------------------
232c Variable physique
233c-----------------------------------------------------------------------
234c
235c tsurf --> tsurfS
236c co2ice --> co2iceS
237c tsoil --> tsoilS
238c emis --> emisS
239c q2 --> q2S
240c qsurf --> qsurfS
241c tauscaling --> tauscalingS
242c
243c-----------------------------------------------------------------------
244
245      call gr_fi_dyn(1,ngridmx,iip1,jjp1,tsurf,tsurfS)
246      call gr_fi_dyn(1,ngridmx,iip1,jjp1,co2ice,co2iceS)
247      call gr_fi_dyn(nsoilmx,ngridmx,iip1,jjp1,tsoil,tsoilS)
248      ! Note: thermal inertia "inertiedat" is in comsoil.h
249      call gr_fi_dyn(nsoilmx,ngridmx,iip1,jjp1,inertiedat,ithS)
250      call gr_fi_dyn(1,ngridmx,iip1,jjp1,emis,emisS)
251      call gr_fi_dyn(llm+1,ngridmx,iip1,jjp1,q2,q2S)
252      call gr_fi_dyn(nqtot,ngridmx,iip1,jjp1,qsurf,qsurfS)
253      call gr_fi_dyn(1,ngridmx,iip1,jjp1,tauscaling,tauscalingS)
254
255c=======================================================================
256c Info pour controler
257c=======================================================================
258
259      ptotal =  0.
260      co2icetotal = 0.
261      DO j=1,jjp1
262         DO i=1,iim
263           ptotal=ptotal+aire(i+(iim+1)*(j-1))*ps(i+(iim+1)*(j-1))/g
264           co2icetotal = co2icetotal +
265     &            co2iceS(i+(iim+1)*(j-1))*aire(i+(iim+1)*(j-1))
266         ENDDO
267      ENDDO
268      write(*,*)'Ancienne grille : masse de l''atm :',ptotal
269      write(*,*)'Ancienne grille : masse de la glace CO2 :',co2icetotal
270
271c-----------------------------------------------------------------------
272c Passage de "ptotal" et "co2icetotal" par tab_cntrl_fi
273c-----------------------------------------------------------------------
274
275      tab_cntrl_fi(49) = ptotal
276      tab_cntrl_fi(50) = co2icetotal
277
278c=======================================================================
279c Ecriture dans le fichier  "start_archive"
280c=======================================================================
281
282c-----------------------------------------------------------------------
283c Ouverture de "start_archive"
284c-----------------------------------------------------------------------
285
286      ierr = NF_OPEN ('start_archive.nc', NF_WRITE,nid)
287 
288c-----------------------------------------------------------------------
289c  si "start_archive" n'existe pas:
290c    1_ ouverture
291c    2_ creation de l'entete dynamique ("ini_archive")
292c-----------------------------------------------------------------------
293c ini_archive:
294c On met dans l'entete le tab_cntrl dynamique (1 a 16)
295c  On y ajoute les valeurs du tab_cntrl_fi (a partir de 51)
296c  En plus les deux valeurs ptotal et co2icetotal (99 et 100)
297c-----------------------------------------------------------------------
298
299      if (ierr.ne.NF_NOERR) then
300         write(*,*)'OK, Could not open file "start_archive.nc"'
301         write(*,*)'So let s create a new "start_archive"'
302         ierr = NF_CREATE('start_archive.nc',
303     &  IOR(NF_CLOBBER,NF_64BIT_OFFSET), nid)
304         call ini_archive(nid,day_ini,phis,ithS,tab_cntrl_fi)
305      endif
306
307c-----------------------------------------------------------------------
308c Ecriture de la coordonnee temps (date en jours)
309c-----------------------------------------------------------------------
310
311      date = day_ini + hour_ini
312      ierr= NF_INQ_VARID(nid,"Time",varid)
313      ierr= NF_INQ_DIMID(nid,"Time",dimid)
314      ierr= NF_INQ_DIMLEN(nid,dimid,timelen)
315      ntime=timelen+1
316
317      write(*,*) "******************"
318      write(*,*) "ntime",ntime
319      write(*,*) "******************"
320#ifdef NC_DOUBLE
321      ierr= NF_PUT_VARA_DOUBLE(nid,varid,ntime,1,date)
322#else
323      ierr= NF_PUT_VARA_REAL(nid,varid,ntime,1,date)
324#endif
325      if (ierr.ne.NF_NOERR) then
326         write(*,*) "time matter ",NF_STRERROR(ierr)
327         stop
328      endif
329
330c-----------------------------------------------------------------------
331c Ecriture des champs  (co2ice,emis,ps,Tsurf,T,u,v,q2,q,qsurf)
332c-----------------------------------------------------------------------
333c ATTENTION: q2 a une couche de plus!!!!
334c    Pour creer un fichier netcdf lisible par grads,
335c    On passe donc une des couches de q2 a part
336c    comme une variable 2D (la couche au sol: "q2surf")
337c    Les lmm autres couches sont nommees "q2atm" (3D)
338c-----------------------------------------------------------------------
339
340      call write_archive(nid,ntime,'co2ice','couche de glace co2',
341     &  'kg/m2',2,co2iceS)
342      call write_archive(nid,ntime,'tauscaling',
343     &  'dust conversion factor',' ',2,tauscalingS)
344      call write_archive(nid,ntime,'emis','grd emis',' ',2,emisS)
345      call write_archive(nid,ntime,'ps','Psurf','Pa',2,ps)
346      call write_archive(nid,ntime,'tsurf','surf T','K',2,tsurfS)
347      call write_archive(nid,ntime,'temp','temperature','K',3,t)
348      call write_archive(nid,ntime,'u','Vent zonal','m.s-1',3,us)
349      call write_archive(nid,ntime,'v','Vent merid','m.s-1',3,vs)
350      call write_archive(nid,ntime,'q2surf','wind variance','m2.s-2',2,
351     .              q2S)
352      call write_archive(nid,ntime,'q2atm','wind variance','m2.s-2',3,
353     .              q2S(1,2))
354
355c-----------------------------------------------------------------------
356c Ecriture du champs  q  ( q[1,nqtot] )
357c-----------------------------------------------------------------------
358      do iq=1,nqtot
359c       write(str2,'(i2.2)') iq
360c        call write_archive(nid,ntime,'q'//str2,'tracer','kg/kg',
361c     .         3,q(1,1,iq))
362        call write_archive(nid,ntime,tname(iq),'tracer','kg/kg',
363     &         3,q(1,1,iq))
364      end do
365c-----------------------------------------------------------------------
366c Ecriture du champs  qsurf  ( qsurf[1,nqtot] )
367c-----------------------------------------------------------------------
368      do iq=1,nqtot
369c       write(str2,'(i2.2)') iq
370c       call write_archive(nid,ntime,'qsurf'//str2,'Tracer on surface',
371c     $  'kg.m-2',2,qsurfS(1,iq))
372        txt=trim(tname(iq))//"_surf"
373        call write_archive(nid,ntime,txt,'Tracer on surface',
374     &  'kg.m-2',2,qsurfS(1,iq))
375      enddo
376
377
378c-----------------------------------------------------------------------
379c Ecriture du champs  tsoil  ( Tg[1,10] )
380c-----------------------------------------------------------------------
381c "tsoil" Temperature au sol definie dans 10 couches dans le sol
382c   Les 10 couches sont lues comme 10 champs
383c  nommees Tg[1,10]
384
385c      do isoil=1,nsoilmx
386c       write(str2,'(i2.2)') isoil
387c       call write_archive(nid,ntime,'Tg'//str2,'Ground Temperature ',
388c     .   'K',2,tsoilS(1,isoil))
389c      enddo
390
391! Write soil temperatures tsoil
392      call write_archive(nid,ntime,'tsoil','Soil temperature',
393     &     'K',-3,tsoilS)
394
395! Write soil thermal inertia
396      call write_archive(nid,ntime,'inertiedat',
397     &     'Soil thermal inertia',
398     &     'J.s-1/2.m-2.K-1',-3,ithS)
399
400! Write (0D) volumetric heat capacity (stored in comsoil.h)
401!      call write_archive(nid,ntime,'volcapa',
402!     &     'Soil volumetric heat capacity',
403!     &     'J.m-3.K-1',0,volcapa)
404! Note: no need to write volcapa, it is stored in "controle" table
405
406      ierr=NF_CLOSE(nid)
407c-----------------------------------------------------------------------
408c Fin
409c-----------------------------------------------------------------------
410
411      write(*,*) "startarchive: all is well that ends well"
412     
413      end
Note: See TracBrowser for help on using the repository browser.