source: trunk/NEMO/OPA_SRC/restart_dimg.h90 @ 392

Last change on this file since 392 was 392, checked in by opalod, 15 years ago

RB:nemo_v1_update_038: first integration of Agrif :

  • add agrif to dynspg_flt_jki.F90
  • cosmetic change of key_AGRIF in key_agrif
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 17.1 KB
Line 
1   !!---------------------------------------------------------------------
2   !!                     ***  restart_dimg.h90  ***
3   !!---------------------------------------------------------------------
4   !!  OPA 9.0 , LOCEAN-IPSL (2005)
5   !! $Header$
6   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
7   !!----------------------------------------------------------------------
8
9   SUBROUTINE rst_write(kt)
10     !!---------------------------------------------------------------------
11     !!                  ***  ROUTINE rst_write  ***
12     !!
13     !! ** Purpose :   Write restart fields in direct access format in mpp.
14     !!      one per process
15     !!
16     !! ** Method  :   each nstock time step , save fields which are necessary
17     !!      for restart
18     !!      Record #1 hold general information on the state of the run
19     !!      Data fields (either 3D or 2D ) starts ar record #2
20     !!
21     !! History :
22     !!        !  91-03  ()  original code
23     !!        !  91-11  (G. Madec)
24     !!        !  92-06  (M. Imbard)  correction restart file
25     !!        !  92-07  (M. Imbard)  split into diawri and rstwri
26     !!        !  98-02  (M. Guyon)  FETI method
27     !!        !  98-05  (G. Roullet)  free surface
28     !!        !  99-11  (M. Imbard)  NetCDF FORMAT with ioipsl
29     !!   8.5  !  03-06  (J.M. Molines)  F90: Free form, mpp support
30     !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization
31     !!----------------------------------------------------------------------
32     !! * Arguments
33     INTEGER, INTENT( in ) ::   kt         ! ocean time-step
34
35     !! * Local declarations
36     INTEGER :: ino0, it0, ipcg0, isor0, itke0
37     INTEGER :: irecl8, irec
38     INTEGER :: jk               ! dummy loop indices
39     INTEGER :: inum = 11        ! temporary logical unit
40     INTEGER :: ios1 , ios2      ! flag for ice and bulk in the current run
41     INTEGER :: ios3             ! flag for free surface.  0 = none 1 = yes.  0 = none 1 = yes
42     INTEGER :: ios4             ! flag for coupled (1) or not (0)
43
44     CHARACTER(LEN=80)  :: clres
45
46     REAL(wp), DIMENSION( 1) ::   zfice, zfblk   ! used only in case of ice & bulk
47     !!----------------------------------------------------------------------
48
49    IF( MOD( kt, nstock ) == 0 .OR. kt == nitend ) THEN
50
51       ! 0. Initializations
52       ! ------------------
53
54       IF(lwp) THEN
55          WRITE(numout,*) ' '
56          WRITE(numout,*) ' rst_write: output done in inum = ',   &
57               inum,' at it= ',kt,' date= ',ndastp
58          WRITE(numout,*) ' -------'
59       ENDIF
60
61       ! Open direct access file, with reclength for 2D wp fields
62       irecl8= jpi * jpj * wp
63       WRITE(clres,'(a,i3.3)') 'restart.output.',narea
64#if defined key_agrif
65       inum = Agrif_Get_Unit()     
66       If(.NOT. Agrif_root() ) clres = TRIM(Agrif_CFixed())//'_'//TRIM(clres)
67#endif
68       OPEN(inum,FILE=clres,FORM='UNFORMATTED', ACCESS='DIRECT', RECL=irecl8 )
69
70
71       ino0  = no
72       it0   = kt
73       ipcg0 = 0
74       isor0 = 0
75       itke0 = 0
76       isor0 = nsolv - 1
77       ipcg0 = 2 - nsolv
78       IF ( lk_zdftke ) itke0=1
79       ! FETI method
80       IF (nsolv == 3) THEN
81          isor0 = 2
82          ipcg0 = 2
83       ENDIF
84
85       ! 1. Write in inum
86       ! ------------------
87
88       ! first record
89       ios1 = 0
90       ios2 = 0
91       ios3 = 0
92       ios4 = 0
93       IF ( lk_ice_lim                           )   ios1 = 1
94       IF ( l_bulk                               )   ios2 = 1
95       IF ( lk_dynspg_flt                        )   ios3 = 1
96       IF ( lk_cpl                               )   ios4 = 1
97
98       WRITE(inum,REC=1) irecl8, ino0, it0, isor0, ipcg0, itke0, &
99          &              nfice, nfbulk , ios1, ios2, ios3, ios4, &
100          &              ndastp, adatrj, jpi, jpj, jpk,  &
101          &              jpni, jpnj, jpnij, narea, jpiglo, jpjglo, &
102          &              nlcit, nlcjt, nldit, nldjt, nleit, nlejt, nimppt, njmppt
103
104       ! prognostic variables
105
106       irec=2
107
108       ! 'before' fields
109       DO jk = 1, jpk
110          WRITE(inum,REC=irec) ub(:,:,jk)   ;    irec = irec +1
111       END DO
112       DO jk = 1, jpk
113          WRITE(inum,REC=irec) vb(:,:,jk)   ;    irec = irec +1
114       END DO
115       DO jk = 1, jpk
116          WRITE(inum,REC=irec) tb(:,:,jk)   ;    irec = irec +1
117       END DO
118       DO jk = 1, jpk
119          WRITE(inum,REC=irec) sb(:,:,jk)   ;    irec = irec +1
120       END DO
121       DO jk = 1, jpk
122          WRITE(inum,REC=irec) rotb(:,:,jk)   ;    irec = irec +1
123       END DO
124       DO jk = 1, jpk
125          WRITE(inum,REC=irec) hdivb(:,:,jk)   ;    irec = irec +1
126       END DO
127
128       ! 'now' fields
129       DO jk = 1, jpk
130          WRITE(inum,REC=irec) un(:,:,jk)   ;   irec = irec +1
131       END DO
132       DO jk = 1, jpk
133          WRITE(inum,REC=irec) vn(:,:,jk)   ;   irec = irec +1
134       END DO
135       DO jk = 1, jpk
136          WRITE(inum,REC=irec) tn(:,:,jk)   ;   irec = irec +1
137       END DO
138       DO jk = 1, jpk
139          WRITE(inum,REC=irec) sn(:,:,jk)   ;   irec = irec +1
140       END DO
141       DO jk = 1, jpk
142          WRITE(inum,REC=irec) rotn(:,:,jk)   ;   irec = irec +1
143       END DO
144       DO jk = 1, jpk
145          WRITE(inum,REC=irec) hdivn(:,:,jk)   ;   irec = irec +1
146       END DO
147
148       ! elliptic solver arrays
149       WRITE(inum,REC=irec ) gcx(1:jpi,1:jpj)   ;   irec = irec +1
150       WRITE(inum,REC=irec ) gcxb(1:jpi,1:jpj)   ;   irec = irec +1
151#if defined key_dynspg_rl
152       ! Rigid-lid formulation (bsf)
153       WRITE(inum,REC=irec ) bsfb(:,:)   ;   irec = irec +1
154       WRITE(inum,REC=irec ) bsfn(:,:)   ;   irec = irec +1
155       WRITE(inum,REC=irec ) bsfd(:,:)   ;   irec = irec +1
156# else
157       ! free surface formulation (ssh)
158       WRITE(inum,REC=irec ) sshb(:,:)   ;   irec = irec +1
159       WRITE(inum,REC=irec ) sshn(:,:)   ;   irec = irec +1
160# if defined key_dynspg_ts
161       ! free surface formulation issued from barotropic loop
162       WRITE(inum,REC=irec ) sshb_b(:,:)   ;   irec = irec +1
163       WRITE(inum,REC=irec ) sshn_b(:,:)   ;   irec = irec +1
164
165       ! horizontal transports issued from barotropic loop
166       WRITE(inum,REC=irec) un_b(:,:)   ;   irec = irec +1
167       WRITE(inum,REC=irec) vn_b(:,:)   ;   irec = irec +1
168# endif
169#endif
170
171       ! TKE arrays
172
173#if defined key_zdftke
174         DO jk = 1, jpk
175            WRITE(inum,REC=irec) en(:,:,jk) ; irec = irec + 1
176         END DO
177#endif
178
179#if defined key_ice_lim
180          zfice(1) = FLOAT( nfice )                                      ! Louvain La Neuve Sea Ice Model
181          WRITE(inum,REC=irec) zfice(:)     ; irec = irec + 1
182          WRITE(inum,REC=irec) sst_io(:,:)  ; irec = irec + 1
183          WRITE(inum,REC=irec) sss_io(:,:)  ; irec = irec + 1
184          WRITE(inum,REC=irec) u_io  (:,:)  ; irec = irec + 1
185          WRITE(inum,REC=irec) v_io  (:,:)  ; irec = irec + 1
186#    if defined key_coupled
187          WRITE(inum,REC=irec) alb_ice(:,:)  ; irec = irec + 1
188#    endif
189#endif
190# if defined key_flx_bulk_monthly || defined key_flx_bulk_daily
191          zfblk(1) = FLOAT( nfbulk )                                 ! Bulk
192          WRITE(inum,REC=irec) zfblk(:)   ; irec = irec + 1
193          WRITE(inum,REC=irec) gsst(:,:)  ; irec = irec + 1
194# endif
195
196       CLOSE(inum)
197    ENDIF
198
199  END SUBROUTINE rst_write
200
201
202  SUBROUTINE rst_read
203    !!---------------------------------------------------------------------
204    !!                  ***  ROUTINE rst_read  ***
205    !! ** Purpose :
206    !!        Read restart fields in direct access format, one per process
207    !!
208    !! ** Method :   Just does the opposit than rst_wri
209    !!
210    !! History :
211    !!        !  91-03  ()  original code
212    !!        !  91-11  (G. Madec)
213    !!        !  92-06  (M. Imbard)  correction restart file
214    !!        !  92-07  (M. Imbard)  split into diawri and rstwri
215    !!        !  98-02  (M. Guyon)  FETI method
216    !!        !  98-05  (G. Roullet)  free surface
217    !!        !  99-11  (M. Imbard)  NetCDF FORMAT with ioipsl
218    !!   8.5  !  03-06  (J.M. Molines)  F90: Free form, mpp support
219    !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization
220    !!----------------------------------------------------------------------
221    USE lib_mpp
222    !! * Local declarations
223    INTEGER :: ino0, it0, ipcg0, isor0, itke0
224    INTEGER :: ino1, it1, isor1, ipcg1, itke1, idast1
225    INTEGER :: iice1, ibulk1
226    INTEGER :: ipi,ipj,ipk, ipni,ipnj,ipnij,iarea
227    INTEGER :: irecl8, irec
228    INTEGER :: ji,jj,jk
229    INTEGER :: ick, inum
230    INTEGER :: ios1, ios2, ios3, ios4
231
232    CHARACTER(LEN=80)  :: clres
233
234    LOGICAL   :: lstop
235
236      REAL(wp), DIMENSION( 1) ::   zfice, zfblk   ! used only in case of ice & bulk
237    !!----------------------------------------------------------------------
238
239
240    ! 0. Initialisations
241    ! ------------------
242
243    inum = 11
244    ino0  = no
245    it0   = nit000
246    ipcg0 = 0
247    isor0 = 0
248    itke0 = 0
249    isor0 = nsolv-1
250    ipcg0 = 2-nsolv
251    IF (lk_zdftke ) itke0 = 1
252    ! FETI method
253    IF( nsolv == 3 ) THEN
254       isor0=2
255       ipcg0=2
256    ENDIF
257
258    IF(lwp) THEN
259       WRITE(numout,*)
260       WRITE(numout,*) ' *** rst_read:  beginning of restart'
261       WRITE(numout,*) ' '
262       WRITE(numout,*) ' the present run :'
263       WRITE(numout,*) '   job number : ', no
264       WRITE(numout,*) '   with nit000 : ', nit000
265       WRITE(numout,*) '   with pcg option ipcg0 : ', ipcg0
266       WRITE(numout,*) '   with sor option isor0 : ', isor0
267       WRITE(numout,*) '   with FETI solver option ipcg0 & isor0 : ', ipcg0,' & ',isor0
268       WRITE(numout,*) '   with tke option itke0 : ', itke0
269       WRITE(numout,*) '   with nfice            : ', nfice
270       WRITE(numout,*) '   with nfbulk           : ', nfbulk
271    ENDIF
272
273    ! Open direct access file, with reclength for 2D wp fields
274    WRITE(clres,'(a,i3.3)') 'restart.',narea
275#if defined key_agrif       
276    inum = Agrif_Get_Unit()
277    If(.NOT. Agrif_root() ) clres = TRIM(Agrif_CFixed())//'_'//TRIM(clres)
278#endif
279    OPEN(inum,FILE=clres,FORM='UNFORMATTED', ACCESS='DIRECT', RECL=8 )
280    READ(inum,REC=1)irecl8
281    CLOSE(inum)
282
283    OPEN(inum,FILE=clres,FORM='UNFORMATTED', ACCESS='DIRECT', RECL=irecl8 )
284
285    ! 1. Read inum
286    ! --------------
287
288
289    READ(inum,REC=1) irecl8, ino1, it1, isor1, ipcg1, itke1, &
290     &  iice1, ibulk1, ios1, ios2, ios3, ios4, &
291     &  idast1, adatrj0,  ipi,ipj,ipk,ipni,ipnj,ipnij,iarea
292
293    ! Performs checks on the file
294
295    ! processor layout changed ? check only on lwp
296    lstop =.FALSE.
297    IF ( ipni /= jpni ) THEN
298       lstop=.TRUE.
299       IF (lwp)  WRITE(numout,*) ' E R R O R : Processor splitting change along I '
300       IF (lwp)  WRITE(numout,*) ' =========='
301    END IF
302
303    IF ( ipnj /= jpnj ) THEN
304       lstop=.TRUE.
305       IF (lwp)  WRITE(numout,*) ' E R R O R : Processor splitting change along J '
306       IF (lwp)  WRITE(numout,*) ' =========='
307    END IF
308
309    IF ( ipnij /= jpnij ) THEN
310       lstop=.TRUE.
311       IF (lwp)  WRITE(numout,*) ' E R R O R : Total number of processors changed '
312       IF (lwp)  WRITE(numout,*) ' =========='
313    END IF
314
315    ick = narea -iarea
316    CALL mpp_sum( ick )
317
318    IF (ick /= 0 ) THEN
319       lstop=.TRUE.
320       IF (lwp)  WRITE(numout,*) ' E R R O R : mismatch in area numbering ...'
321       IF (lwp)  WRITE(numout,*) ' =========='
322    END IF
323
324    IF(lwp) THEN
325       WRITE(numout,*)
326       WRITE(numout,*) ' READ inum with '
327       WRITE(numout,*) '   job number : ', ino1
328       WRITE(numout,*) '   with time step it : ', it1
329       WRITE(numout,*) '   with pcg option ipcg1 : ', ipcg1
330       WRITE(numout,*) '   with sor option isor1 : ', isor1
331       WRITE(numout,*) '   with tke option itke1 : ', itke1
332       WRITE(numout,*) '   with FETI solver option ipcg1 + isor1 : ', ipcg1 + isor1
333       WRITE(numout,*)
334    ENDIF
335
336    ! Control of date
337
338    IF( (it0-it1) /= 1 .AND. nrstdt /= 0 ) THEN
339       lstop=.TRUE.
340       IF(lwp) THEN
341          WRITE(numout,*) ' ===>>>> : problem with nit000 for the restart'
342          WRITE(numout,*) ' =======                               ======='
343          WRITE(numout,*) ' we stop. verify the file'
344          WRITE(numout,*) ' or rerun with the value  0 for the'
345          WRITE(numout,*) ' control of time parameter  nrstdt'
346          WRITE(numout,*)
347       ENDIF
348    ENDIF
349
350    IF ( nrstdt /= 2 ) THEN  !  Compatibility with OPA8
351                             !  the beginning of the new run is ndate0 read in the namelist
352                             !  adatrj0 is recalculated assuming constant time step.
353
354       adatrj0 =  ( FLOAT( nit000-1 ) * rdttra(1) ) / rday
355    ELSE                 !  restart option nrstdt = 2
356                         !  both adatrj0 and ndastp are read in the restart file.
357       ndastp = idast1
358    ENDIF
359
360    IF (lstop ) STOP 'rst_read'
361
362    irec=2
363
364    ! 'before' fields
365    DO jk = 1, jpk
366       READ(inum,REC=irec) ub(:,:,jk)   ;   irec = irec +1
367    END DO
368    DO jk = 1, jpk
369       READ(inum,REC=irec) vb(:,:,jk)   ;   irec = irec +1
370    END DO
371    DO jk = 1, jpk
372       READ(inum,REC=irec) tb(:,:,jk)   ;   irec = irec +1
373    END DO
374    DO jk = 1, jpk
375       READ(inum,REC=irec) sb(:,:,jk)   ;   irec = irec +1
376    END DO
377    DO jk = 1, jpk
378       READ(inum,REC=irec) rotb(:,:,jk)   ;   irec = irec +1
379    END DO
380    DO jk = 1, jpk
381       READ(inum,REC=irec) hdivb(:,:,jk)   ;   irec = irec +1
382    END DO
383
384    ! 'now' fields
385    DO jk = 1, jpk
386       READ(inum,REC=irec) un(:,:,jk)   ;   irec = irec +1
387    END DO
388    DO jk = 1, jpk
389       READ(inum,REC=irec) vn(:,:,jk)   ;   irec = irec +1
390    END DO
391    DO jk = 1, jpk
392       READ(inum,REC=irec) tn(:,:,jk)   ;   irec = irec +1
393    END DO
394    DO jk = 1, jpk
395       READ(inum,REC=irec) sn(:,:,jk)   ;   irec = irec +1
396    END DO
397    DO jk = 1, jpk
398       READ(inum,REC=irec) rotn(:,:,jk)   ;   irec = irec +1
399    END DO
400    DO jk = 1, jpk
401       READ(inum,REC=irec) hdivn(:,:,jk)   ;   irec = irec +1
402    END DO
403
404    ! elliptic solver arrays
405    READ(inum,REC=irec ) gcx(1:jpi,1:jpj)   ;   irec = irec +1
406    READ(inum,REC=irec ) gcxb(1:jpi,1:jpj)   ;   irec = irec +1
407#if defined key_dynspg_rl
408    ! Rigid-lid formulation (bsf)
409    READ(inum,REC=irec ) bsfb(:,:)   ;   irec = irec +1
410    READ(inum,REC=irec ) bsfn(:,:)   ;   irec = irec +1
411    READ(inum,REC=irec ) bsfd(:,:)   ;   irec = irec +1
412#else
413    ! free surface formulation (eta)
414    READ(inum,REC=irec ) sshb(:,:)   ;   irec = irec +1
415    READ(inum,REC=irec ) sshn(:,:)   ;   irec = irec +1
416# if defined key_dynspg_ts
417    ! free surface formulation issued from barotropic loop
418    READ(inum,REC=irec ) sshb_b(:,:)   ;   irec = irec +1
419    READ(inum,REC=irec ) sshn_b(:,:)   ;   irec = irec +1
420    ! horizontal transports issued from barotropic loop
421    READ(inum,REC=irec) un_b(:,:)   ;   irec = irec +1
422    READ(inum,REC=irec) vn_b(:,:)   ;   irec = irec +1
423# endif
424#endif
425
426    ! TKE arrays
427
428#if defined key_zdftke
429    IF ( itke1 == 1 ) THEN
430       DO jk = 1, jpk
431          READ(inum,REC=irec) en(:,:,jk)   ;   irec = irec +1
432       END DO
433    ELSE
434       IF(lwp) THEN
435          WRITE(numout,*) ' ===>>>> : the previous restart file didnot used  tke scheme'
436          WRITE(numout,*) ' =======                ======='
437       ENDIF
438       nrstdt = 2
439    ENDIF
440#endif
441#if defined key_ice_lim
442    ! Louvain La Neuve Sea Ice Model
443    ! check if it was in the previous run
444    IF ( ios1 == 1 ) THEN
445       READ(inum,REC=irec) zfice(:)    ; irec = irec + 1
446       READ(inum,REC=irec) sst_io(:,:) ; irec = irec + 1
447       READ(inum,REC=irec) sss_io(:,:) ; irec = irec + 1
448       READ(inum,REC=irec) u_io  (:,:) ; irec = irec + 1
449       READ(inum,REC=irec) v_io  (:,:) ; irec = irec + 1
450#  if defined key_coupled
451       READ(inum,REC=irec) alb_ice(:,:) ; irec = irec + 1
452#  endif
453    ENDIF
454    IF ( zfice(1) /= FLOAT(nfice) .OR. ios1 == 0 ) THEN
455         IF(lwp) WRITE(numout,*)
456         IF(lwp) WRITE(numout,*) 'rst_read :  LLN sea Ice Model => Ice initialization'
457         IF(lwp) WRITE(numout,*)
458         sst_io(:,:) = sst_io(:,:) + ( nfice-1 )*( tn(:,:,1) + rt0 )
459         sss_io(:,:) = sss_io(:,:) + ( nfice-1 )*  sn(:,:,1)
460         DO jj = 2, jpj
461            DO ji = 2, jpi
462               u_io(ji,jj) = u_io(ji,jj) + (nfice-1)*0.5*( un(ji-1,jj,1)+un(ji-1,jj-1,1) )
463               v_io(ji,jj) = v_io(ji,jj) + (nfice-1)*0.5*( vn(ji,jj-1,1)+vn(ji-1,jj-1,1) )
464            END DO
465         END DO
466#    if defined key_coupled
467         alb_ice(:,:) = 0.8 * tmask(:,:,1)
468#    endif
469    ENDIF
470#endif
471#if defined key_flx_bulk_monthly || defined key_flx_bulk_daily
472      ! bulk forcing
473      IF( ios2 == 1 ) THEN
474         READ(inum,REC=irec) zfblk(:)   ; irec = irec + 1
475         READ(inum,REC=irec) gsst (:,:) ; irec = irec + 1
476      ENDIF
477      IF( zfblk(1) /= FLOAT(nfbulk)  .OR. ios2 == 0 ) THEN
478         IF(lwp) WRITE(numout,*)
479         IF(lwp) WRITE(numout,*) 'rst_read :  Bulk forcing ==> Initialization '
480         IF(lwp) WRITE(numout,*)
481         gsst(:,:) = 0.e0
482         gsst(:,:) = gsst(:,:) + ( nfbulk-1 )*( tn(:,:,1) + rt0 )
483      ENDIF
484#endif
485    CLOSE(inum)
486  ! In case of restart with neuler = 0 then put all before fields = to now fields
487    IF ( neuler == 0 ) THEN
488      tb(:,:,:)=tn(:,:,:)
489      sb(:,:,:)=sn(:,:,:)
490      ub(:,:,:)=un(:,:,:)
491      vb(:,:,:)=vn(:,:,:)
492      rotb(:,:,:)=rotn(:,:,:)
493      hdivb(:,:,:)=hdivn(:,:,:)
494#if defined key_dynspg_rl
495      bsfb(:,:)=bsfn(:,:)      ! rigid lid
496#else
497      sshb(:,:)=sshn(:,:)      ! free surface formulation (eta)
498#endif
499    ENDIF
500
501
502  END SUBROUTINE rst_read
Note: See TracBrowser for help on using the repository browser.