source: trunk/NEMO/OPA_SRC/restart.F90 @ 783

Last change on this file since 783 was 783, checked in by smasson, 13 years ago

write multiple restarts, see ticket:44

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 19.0 KB
Line 
1MODULE restart
2   !!======================================================================
3   !!                     ***  MODULE  restart  ***
4   !! Ocean restart :  write the ocean restart file
5   !!======================================================================
6   !! History :        !  99-11  (M. Imbard)  Original code
7   !!             8.5  !  02-08  (G. Madec)  F90: Free form
8   !!             9.0  !  05-11  (V. Garnier) Surface pressure gradient organization
9   !!             9.0  !  06-07  (S. Masson)  use IOM for restart
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   rst_opn    : open the ocean restart file
14   !!   rst_write  : write the ocean restart file
15   !!   rst_read   : read the ocean restart file
16   !!----------------------------------------------------------------------
17   USE dom_oce         ! ocean space and time domain
18   USE oce             ! ocean dynamics and tracers
19   USE phycst          ! physical constants
20   USE daymod          ! calendar
21   USE ice_oce         ! ice variables
22   USE blk_oce         ! bulk variables
23   USE cpl_oce, ONLY : lk_cpl              !
24   USE in_out_manager  ! I/O manager
25   USE iom             ! I/O module
26   USE ini1d           ! re-initialization of u-v mask for the 1D configuration
27   USE zpshde          ! partial step: hor. derivative (zps_hde routine)
28   USE eosbn2          ! equation of state            (eos bn2 routine)
29   USE trdmld_oce      ! ocean active mixed layer tracers trends variables
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   rst_opn    ! routine called by step module
35   PUBLIC   rst_write  ! routine called by step module
36   PUBLIC   rst_read   ! routine called by opa  module
37
38   LOGICAL, PUBLIC ::   lrst_oce                  !: logical to control the oce restart write
39   INTEGER, PUBLIC ::   numror, numrow            !: logical unit for cean restart (read and write)
40
41   !! * Substitutions
42#  include "vectopt_loop_substitute.h90"
43   !!----------------------------------------------------------------------
44   !!  OPA 9.0 , LOCEAN-IPSL (2006)
45   !! $Header: /home/opalod/NEMOCVSROOT/NEMO/OPA_SRC/restart.F90,v 1.27 2007/06/05 10:35:19 opalod Exp $
46   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
47   !!----------------------------------------------------------------------
48
49CONTAINS
50
51   SUBROUTINE rst_opn( kt )
52      !!---------------------------------------------------------------------
53      !!                   ***  ROUTINE rst_opn  ***
54      !!                     
55      !! ** Purpose : + initialization (should be read in the namelist) of nitrst
56      !!              + open the restart when we are one time step before nitrst
57      !!                   - restart header is defined when kt = nitrst-1
58      !!                   - restart data  are written when kt = nitrst
59      !!              + define lrst_oce to .TRUE. when we need to define or write the restart
60      !!----------------------------------------------------------------------
61      INTEGER, INTENT(in) ::   kt     ! ocean time-step
62      !!
63      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step deine as a character
64      CHARACTER(LEN=50)   ::   clname   ! ice output restart file name
65      !!----------------------------------------------------------------------
66      !
67      IF( kt == nit000 ) THEN   ! default definitions
68         lrst_oce = .FALSE.   
69         nitrst = nitend
70      ENDIF
71      IF( MOD( kt - 1, nstock ) == 0 ) THEN   
72         ! we use kt - 1 and not kt - nit000 to keep the same periodicity from the beginning of the experiment
73         nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing
74         IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run
75      ENDIF
76
77      ! to get better performances with NetCDF format:
78      ! we open and define the ocean restart file one time step before writing the data (-> at nitrst - 1)
79      ! except if we write ocean restart files every time step or if an ocean restart file was writen at nitend - 1
80      IF( kt == nitrst - 1 .OR. nstock == 1 .OR. ( kt == nitend .AND. .NOT. lrst_oce ) ) THEN
81         ! beware of the format used to write kt (default is i8.8, that should be large enough...)
82         IF( nitrst > 999999999 ) THEN   ;   WRITE(clkt, *       ) nitrst
83         ELSE                            ;   WRITE(clkt, '(i8.8)') nitrst
84         ENDIF
85         ! create the file
86         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart"
87         IF(lwp) THEN
88            WRITE(numout,*)
89            SELECT CASE ( jprstlib )
90            CASE ( jprstdimg )   ;   WRITE(numout,*) '             open ocean restart binary file: '//clname
91            CASE DEFAULT         ;   WRITE(numout,*) '             open ocean restart NetCDF file: '//clname
92            END SELECT
93            IF( kt == nitrst - 1 ) THEN   ;   WRITE(numout,*) '             kt = nitrst - 1 = ', kt,' date= ', ndastp
94            ELSE                          ;   WRITE(numout,*) '             kt = '             , kt,' date= ', ndastp
95            ENDIF
96         ENDIF
97
98         CALL iom_open( clname, numrow, ldwrt = .TRUE., kiolib = jprstlib )
99         lrst_oce = .TRUE.
100      ENDIF
101      !
102   END SUBROUTINE rst_opn
103
104
105   SUBROUTINE rst_write( kt )
106      !!---------------------------------------------------------------------
107      !!                   ***  ROUTINE rstwrite  ***
108      !!                     
109      !! ** Purpose :   Write restart fields in the format corresponding to jprstlib
110      !!
111      !! ** Method  :   Write in numrow when kt == nitrst in NetCDF
112      !!      file, save fields which are necessary for restart
113      !!----------------------------------------------------------------------
114      INTEGER, INTENT(in) ::   kt   ! ocean time-step
115      !!----------------------------------------------------------------------
116
117      IF( kt == nitrst ) THEN
118         IF(lwp) WRITE(numout,*)
119         IF(lwp) WRITE(numout,*) 'rst_write : write oce restart file  kt =', kt
120         IF(lwp) WRITE(numout,*) '~~~~~~~'         
121      ENDIF
122
123      ! calendar control
124      CALL iom_rstput( kt, nitrst, numrow, 'kt'     , REAL( kt    , wp) )   ! time-step
125      CALL iom_rstput( kt, nitrst, numrow, 'ndastp' , REAL( ndastp, wp) )   ! date
126      CALL iom_rstput( kt, nitrst, numrow, 'adatrj' , adatrj            )   ! number of elapsed days since
127      !                                                                     ! the begining of the run [s]
128      CALL iom_rstput( kt, nitrst, numrow, 'rdt'    , rdt               )   ! dynamics time step
129      CALL iom_rstput( kt, nitrst, numrow, 'rdttra1', rdttra(1)         )   ! surface tracer time step
130
131      ! prognostic variables
132      CALL iom_rstput( kt, nitrst, numrow, 'ub'     , ub      )   
133      CALL iom_rstput( kt, nitrst, numrow, 'vb'     , vb      )
134      CALL iom_rstput( kt, nitrst, numrow, 'tb'     , tb      )
135      CALL iom_rstput( kt, nitrst, numrow, 'sb'     , sb      )
136      CALL iom_rstput( kt, nitrst, numrow, 'rotb'   , rotb    )
137      CALL iom_rstput( kt, nitrst, numrow, 'hdivb'  , hdivb   )
138      CALL iom_rstput( kt, nitrst, numrow, 'un'     , un      )
139      CALL iom_rstput( kt, nitrst, numrow, 'vn'     , vn      )
140      IF( lk_vvl ) CALL iom_rstput( kt, nitrst, numrow, 'wn'     , wn      )
141      CALL iom_rstput( kt, nitrst, numrow, 'tn'     , tn      )
142      CALL iom_rstput( kt, nitrst, numrow, 'sn'     , sn      )
143      CALL iom_rstput( kt, nitrst, numrow, 'rotn'   , rotn    )
144      CALL iom_rstput( kt, nitrst, numrow, 'hdivn'  , hdivn   )
145
146#if defined key_ice_lim       
147      CALL iom_rstput( kt, nitrst, numrow, 'nfice'  , REAL( nfice, wp) )   !  ice computation frequency
148      CALL iom_rstput( kt, nitrst, numrow, 'sst_io' , sst_io  )
149      CALL iom_rstput( kt, nitrst, numrow, 'sss_io' , sss_io  )
150      CALL iom_rstput( kt, nitrst, numrow, 'u_io'   , u_io    )
151      CALL iom_rstput( kt, nitrst, numrow, 'v_io'   , v_io    )
152# if defined key_coupled
153      CALL iom_rstput( kt, nitrst, numrow, 'alb_ice', alb_ice )
154# endif
155#endif
156#if defined key_flx_bulk_monthly || defined key_flx_bulk_daily || defined key_flx_core 
157      CALL iom_rstput( kt, nitrst, numrow, 'nfbulk' , REAL( nfbulk, wp) )   !  bulk computation frequency
158      CALL iom_rstput( kt, nitrst, numrow, 'gsst'   , gsst    )
159#endif
160
161      IF( nn_dynhpg_rst == 1 .OR. lk_vvl ) THEN
162         CALL iom_rstput( kt, nitrst, numrow, 'rhd' , rhd  )
163         CALL iom_rstput( kt, nitrst, numrow, 'rhop', rhop )
164         IF( ln_zps ) THEN
165            CALL iom_rstput( kt, nitrst, numrow, 'gtu' , gtu )
166            CALL iom_rstput( kt, nitrst, numrow, 'gsu' , gsu )
167            CALL iom_rstput( kt, nitrst, numrow, 'gru' , gru )
168            CALL iom_rstput( kt, nitrst, numrow, 'gtv' , gtv )
169            CALL iom_rstput( kt, nitrst, numrow, 'gsv' , gsv )
170            CALL iom_rstput( kt, nitrst, numrow, 'grv' , grv )
171         ENDIF
172      ENDIF
173
174      IF( kt == nitrst ) THEN
175         CALL iom_close( numrow )     ! close the restart file (only at last time step)
176         IF( .NOT. lk_trdmld )   lrst_oce = .FALSE.
177      ENDIF
178      !
179   END SUBROUTINE rst_write
180
181
182   SUBROUTINE rst_read
183      !!----------------------------------------------------------------------
184      !!                   ***  ROUTINE rst_read  ***
185      !!
186      !! ** Purpose :   Read files for restart (format fixed by jprstlib )
187      !!
188      !! ** Method  :   Read the previous fields on the NetCDF/binary file
189      !!      the first record indicates previous characterics
190      !!      after control with the present run, we read :
191      !!      - prognostic variables on the second record
192      !!      - elliptic solver arrays
193      !!      - barotropic stream function arrays ("key_dynspg_rl" defined)
194      !!        or free surface arrays
195      !!      - tke arrays (lk_zdftke=T)
196      !!      for this last three records,  the previous characteristics
197      !!      could be different with those used in the present run.
198      !!
199      !!   According to namelist parameter nrstdt,
200      !!       nrstdt = 0  no control on the date (nit000 is  arbitrary).
201      !!       nrstdt = 1  we verify that nit000 is equal to the last
202      !!                   time step of previous run + 1.
203      !!       In both those options, the  exact duration of the experiment
204      !!       since the beginning (cumulated duration of all previous restart runs)
205      !!       is not stored in the restart and is assumed to be (nit000-1)*rdt.
206      !!       This is valid is the time step has remained constant.
207      !!
208      !!       nrstdt = 2  the duration of the experiment in days (adatrj)
209      !!                    has been stored in the restart file.
210      !!----------------------------------------------------------------------
211      REAL(wp) ::   zcoef, zkt, zrdt, zrdttra1, zndastp, znfice, znfbulk
212#if defined key_ice_lim
213      INTEGER  ::   ji, jj
214#endif
215      !!----------------------------------------------------------------------
216
217      IF(lwp) THEN                                             ! Contol prints
218         WRITE(numout,*)
219         SELECT CASE ( jprstlib )
220         CASE ( jpnf90 )
221            WRITE(numout,*) 'rst_read : read oce NetCDF restart file'
222         CASE ( jprstdimg )
223            WRITE(numout,*) 'rst_read : read oce binary restart file'
224         END SELECT
225         WRITE(numout,*) '~~~~~~~~'
226
227         WRITE(numout,*) ' *** Info on the present job : '
228         WRITE(numout,*) '   time-step           : ', nit000
229         WRITE(numout,*) '   date ndastp         : ', ndastp
230         WRITE(numout,*)
231         WRITE(numout,*) ' *** restart option'
232         SELECT CASE ( nrstdt )
233         CASE ( 0 ) 
234            WRITE(numout,*) ' nrstdt = 0 no control of nit000'
235         CASE ( 1 ) 
236            WRITE(numout,*) ' nrstdt = 1 we control the date of nit000'
237         CASE ( 2 )
238            WRITE(numout,*) ' nrstdt = 2 the date adatrj is read in restart file'
239         CASE DEFAULT
240            WRITE(numout,*) '  ===>>>> nrstdt not equal 0, 1 or 2 : no control of the date'
241            WRITE(numout,*) '  =======                  ========='
242         END SELECT
243         WRITE(numout,*)
244      ENDIF
245
246      CALL iom_open( 'restart', numror, kiolib = jprstlib )
247
248      ! Calendar informations
249      CALL iom_get( numror, 'kt'     , zkt      )   ! time-step
250      CALL iom_get( numror, 'ndastp' , zndastp  )   ! date
251      IF(lwp) THEN
252         WRITE(numout,*)
253         WRITE(numout,*) ' *** Info on the restart file read : '
254         WRITE(numout,*) '   time-step           : ', NINT( zkt )
255         WRITE(numout,*) '   date ndastp         : ', NINT( zndastp )
256         WRITE(numout,*)
257      ENDIF
258      ! Control of date
259      IF( nit000 - NINT( zkt )  /= 1 .AND. nrstdt /= 0 ) &
260           & CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart', &
261           & ' verify the restart file or rerun with nrstdt = 0 (namelist)' )
262      ! re-initialisation of  adatrj0
263      adatrj0 = ( REAL( nit000-1, wp ) * rdttra(1) ) / rday
264      IF ( nrstdt == 2 ) THEN
265         ! by default ndatsp has been set to ndate0 in dom_nam
266         ! ndate0 has been read in the namelist (standard OPA 8)
267         ! here when nrstdt=2 we keep the  final date of previous run
268         ndastp = NINT( zndastp )
269         CALL iom_get( numror, 'adatrj', adatrj )   ! number of elapsed days since the begining of last run
270      ENDIF
271      ! Check dynamics and tracer time-step consistency and force Euler restart if changed
272      IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 )   THEN
273         CALL iom_get( numror, 'rdt', zrdt )
274         IF( zrdt /= rdt )   neuler = 0
275      ENDIF
276      IF( iom_varid( numror, 'rdttra1', ldstop = .FALSE. ) > 0 )   THEN
277         CALL iom_get( numror, 'rdttra1', zrdttra1 )
278         IF( zrdttra1 /= rdttra(1) )   neuler = 0
279      ENDIF
280      !
281      !                                                       ! Read prognostic variables
282      CALL iom_get( numror, jpdom_autoglo, 'ub'   , ub    )        ! before i-component velocity
283      CALL iom_get( numror, jpdom_autoglo, 'vb'   , vb    )        ! before j-component velocity
284      CALL iom_get( numror, jpdom_autoglo, 'tb'   , tb    )        ! before temperature
285      CALL iom_get( numror, jpdom_autoglo, 'sb'   , sb    )        ! before salinity
286      CALL iom_get( numror, jpdom_autoglo, 'rotb' , rotb  )        ! before curl
287      CALL iom_get( numror, jpdom_autoglo, 'hdivb', hdivb )        ! before horizontal divergence
288      CALL iom_get( numror, jpdom_autoglo, 'un'   , un    )        ! now    i-component velocity
289      CALL iom_get( numror, jpdom_autoglo, 'vn'   , vn    )        ! now    j-component velocity
290      IF( lk_vvl ) CALL iom_get( numror, jpdom_autoglo, 'wn'   , wn    )        ! now    k-component velocity
291      CALL iom_get( numror, jpdom_autoglo, 'tn'   , tn    )        ! now    temperature
292      CALL iom_get( numror, jpdom_autoglo, 'sn'   , sn    )        ! now    salinity
293      CALL iom_get( numror, jpdom_autoglo, 'rotn' , rotn  )        ! now    curl
294      CALL iom_get( numror, jpdom_autoglo, 'hdivn', hdivn )        ! now    horizontal divergence
295
296
297      IF( neuler == 0 ) THEN                                  ! Euler restart (neuler=0)
298         tb   (:,:,:) = tn   (:,:,:)                             ! all before fields set to now field values
299         sb   (:,:,:) = sn   (:,:,:)
300         ub   (:,:,:) = un   (:,:,:)
301         vb   (:,:,:) = vn   (:,:,:)
302         rotb (:,:,:) = rotn (:,:,:)
303         hdivb(:,:,:) = hdivn(:,:,:)
304      ENDIF
305
306      !!sm: TO BE MOVED IN NEW SURFACE MODULE...
307
308#if defined key_ice_lim
309      ! Louvain La Neuve Sea Ice Model
310      IF( iom_varid( numror, 'nfice', ldstop = .FALSE. ) > 0 ) then
311         CALL iom_get( numror             , 'nfice'  , znfice  )   ! ice computation frequency
312         CALL iom_get( numror, jpdom_autoglo, 'sst_io' , sst_io  )
313         CALL iom_get( numror, jpdom_autoglo, 'sss_io' , sss_io  )
314         CALL iom_get( numror, jpdom_autoglo, 'u_io'   , u_io    )
315         CALL iom_get( numror, jpdom_autoglo, 'v_io'   , v_io    )
316# if defined key_coupled
317         CALL iom_get( numror, jpdom_autoglo, 'alb_ice', alb_ice )
318# endif
319         IF( znfice /= REAL( nfice, wp ) ) THEN      ! if nfice changed between 2 runs
320            zcoef = REAL( nfice-1, wp ) / znfice
321            sst_io(:,:) = zcoef * sst_io(:,:)
322            sss_io(:,:) = zcoef * sss_io(:,:)
323            u_io  (:,:) = zcoef * u_io  (:,:)
324            v_io  (:,:) = zcoef * v_io  (:,:)
325         ENDIF
326      ELSE
327         IF(lwp) WRITE(numout,*)
328         IF(lwp) WRITE(numout,*) 'rst_read :  LLN sea Ice Model => Ice initialization'
329         IF(lwp) WRITE(numout,*)
330         zcoef = REAL( nfice-1, wp )
331         sst_io(:,:) = zcoef *( tn(:,:,1) + rt0 )          !!bug a explanation is needed here!
332         sss_io(:,:) = zcoef *  sn(:,:,1)
333         zcoef = 0.5 * REAL( nfice-1, wp )
334         DO jj = 2, jpj
335            DO ji = fs_2, jpi   ! vector opt.
336               u_io(ji,jj) = zcoef * ( un(ji-1,jj  ,1) + un(ji-1,jj-1,1) )
337               v_io(ji,jj) = zcoef * ( vn(ji  ,jj-1,1) + vn(ji-1,jj-1,1) )
338            END DO
339         END DO
340# if defined key_coupled
341         alb_ice(:,:) = 0.8 * tmask(:,:,1)
342# endif
343      ENDIF
344#endif
345#if defined key_flx_bulk_monthly || defined key_flx_bulk_daily || defined key_flx_core 
346      ! Louvain La Neuve Sea Ice Model
347      IF( iom_varid( numror, 'nfbulk', ldstop = .FALSE. ) > 0 ) THEN
348         CALL iom_get( numror             , 'nfbulk', znfbulk )   ! bulk computation frequency
349         CALL iom_get( numror, jpdom_autoglo, 'gsst'  , gsst    )
350         IF( znfbulk /= REAL(nfbulk, wp) ) THEN      ! if you change nfbulk between 2 runs
351            zcoef = REAL( nfbulk-1, wp ) / znfbulk
352            gsst(:,:) = zcoef * gsst(:,:)
353         ENDIF
354      ELSE
355         IF(lwp) WRITE(numout,*)
356         IF(lwp) WRITE(numout,*) 'rst_read :  LLN sea Ice Model => Ice initialization'
357         IF(lwp) WRITE(numout,*)
358         gsst(:,:) = REAL( nfbulk - 1, wp )*( tn(:,:,1) + rt0 )
359      ENDIF
360#endif
361
362      !!sm: end of TO BE MOVED IN NEW SURFACE MODULE...
363
364      IF( iom_varid( numror, 'rhd', ldstop = .FALSE. ) > 0 ) THEN
365         CALL iom_get( numror, jpdom_autoglo, 'rhd' , rhd  )
366         CALL iom_get( numror, jpdom_autoglo, 'rhop', rhop )
367      ELSE
368         CALL eos( tb, sb, rhd, rhop )        ! before potential and in situ densities
369      ENDIF
370      IF( ln_zps .AND. .NOT. lk_cfg_1d ) THEN
371         IF( iom_varid( numror, 'gtu', ldstop = .FALSE. ) > 0 ) THEN
372            CALL iom_get( numror, jpdom_autoglo, 'gtu' , gtu )
373            CALL iom_get( numror, jpdom_autoglo, 'gsu' , gsu )
374            CALL iom_get( numror, jpdom_autoglo, 'gru' , gru )
375            CALL iom_get( numror, jpdom_autoglo, 'gtv' , gtv )
376            CALL iom_get( numror, jpdom_autoglo, 'gsv' , gsv )
377            CALL iom_get( numror, jpdom_autoglo, 'grv' , grv )
378         ELSE
379            CALL zps_hde( nit000, tb , sb , rhd,   &  ! Partial steps: before Horizontal DErivative
380               &                  gtu, gsu, gru,   &  ! of t, s, rd at the bottom ocean level
381               &                  gtv, gsv, grv )
382         ENDIF
383      ENDIF
384      !
385   END SUBROUTINE rst_read
386
387
388   !!=====================================================================
389END MODULE restart
Note: See TracBrowser for help on using the repository browser.