New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 508 for trunk/NEMO/OPA_SRC – NEMO

Changeset 508 for trunk/NEMO/OPA_SRC


Ignore:
Timestamp:
2006-10-03T17:58:55+02:00 (18 years ago)
Author:
opalod
Message:

nemo_v1_update_071:RB: add iom for restart and reorganization of restart

Location:
trunk/NEMO/OPA_SRC
Files:
15 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMO/OPA_SRC/DIA/diaptr.F90

    r460 r508  
    55   !!                 (please no more than 2 lines) 
    66   !!===================================================================== 
     7   !! History :  9.0  !  03-09  (C. Talandir, G. Madec)  Original code 
     8   !!            9.0  !  06-01  (A. Biastoch)  Allow sub-basins computation 
     9   !!---------------------------------------------------------------------- 
     10 
    711   !!---------------------------------------------------------------------- 
    812   !!   dia_ptr      : Poleward Transport Diagnostics module 
     
    1418   !!                : flux array; Generic interface: ptr_vj_3d, ptr_vj_2d 
    1519   !!---------------------------------------------------------------------- 
    16    !! History : 
    17    !!   9.0  !  03-09  (C. Talandir, G. Madec)  Original code 
    18    !!   9.0  !  06-01  (A. Biastoch)  Allow sub-basins computation 
    19    !!---------------------------------------------------------------------- 
    20    !! * Modules used 
    2120   USE oce           ! ocean dynamics and active tracers 
    2221   USE dom_oce       ! ocean space and time domain 
     
    2625   USE dianam 
    2726   USE phycst 
    28    USE ioipsl          ! NetCDF IPSL library 
     27   USE iom 
     28   USE ioipsl          
    2929   USE daymod 
    3030 
     
    3636   END INTERFACE 
    3737 
    38    !! *  Routine accessibility 
    39    PUBLIC dia_ptr_init   ! call in opa module 
    40    PUBLIC dia_ptr        ! call in step module 
    41    PUBLIC ptr_vj         ! call by tra_ldf & tra_adv routines 
    42    PUBLIC ptr_vjk        ! call by tra_ldf & tra_adv routines 
    43  
    44    !! * Share Module variables 
    45    LOGICAL, PUBLIC ::       & !!! ** init namelist (namptr) ** 
    46       ln_diaptr = .FALSE.,  &  !: Poleward transport flag (T) or not (F) 
    47       ln_subbas = .FALSE.      !: Atlantic/Pacific/Indian basins calculation 
    48    INTEGER, PUBLIC ::       & !!: ** ptr namelist (namptr) ** 
    49       nf_ptr = 15              !: frequency of ptr computation 
    50    REAL(wp), PUBLIC, DIMENSION(jpj) ::   &   !!: poleward transport 
    51       pht_adv, pst_adv,     &  !: heat and salt: advection 
    52       pht_ove, pst_ove,     &  !: heat and salt: overturning 
    53       pht_ldf, pst_ldf,     &  !: heat and salt: lateral diffusion 
    54 #if defined key_diaeiv 
    55       pht_eiv, pst_eiv,     &  !: heat and salt: bolus advection 
    56 #endif 
    57       ht_atl,ht_ind,ht_pac, &  !: heat 
    58       st_atl,st_ind,st_pac     !: salt 
    59    REAL(wp),DIMENSION(jpi,jpj) ::   & 
    60       abasin,pbasin,ibasin     !: return function value 
     38   PUBLIC   dia_ptr_init   ! call in opa module 
     39   PUBLIC   dia_ptr        ! call in step module 
     40   PUBLIC   ptr_vj         ! call by tra_ldf & tra_adv routines 
     41   PUBLIC   ptr_vjk        ! call by tra_ldf & tra_adv routines 
     42 
     43   !!! ** init namelist (namptr) 
     44   LOGICAL , PUBLIC                 ::   ln_diaptr = .FALSE.   !: Poleward transport flag (T) or not (F) 
     45   LOGICAL , PUBLIC                 ::   ln_subbas = .FALSE.   !: Atlantic/Pacific/Indian basins calculation 
     46   INTEGER , PUBLIC                 ::   nf_ptr = 15           !: frequency of ptr computation 
     47 
     48   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_adv, pst_adv      !: heat and salt poleward transport: advection 
     49   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ove, pst_ove      !: heat and salt poleward transport: overturning 
     50   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_ldf, pst_ldf      !: heat and salt poleward transport: lateral diffusion 
     51#if defined key_diaeiv 
     52   REAL(wp), PUBLIC, DIMENSION(jpj) ::   pht_eiv, pst_eiv      !: heat and salt poleward transport: bolus advection 
     53#endif 
     54   REAL(wp), PUBLIC, DIMENSION(jpj) ::   ht_atl,ht_ind,ht_pac  !: heat 
     55   REAL(wp), PUBLIC, DIMENSION(jpj) ::   st_atl,st_ind,st_pac  !: salt 
     56 
    6157      
    6258 
    63    !! Module variables 
    64    REAL(wp), DIMENSION(jpj,jpk) ::   &   
    65       tn_jk  , sn_jk  ,  &  !: "zonal" mean temperature and salinity 
    66       v_msf_atl       ,  &  !: "meridional" Stream-Function 
    67       v_msf_glo       ,  &  !: "meridional" Stream-Function 
    68       v_msf_ipc       ,  &  !: "meridional" Stream-Function 
    69 #if defined key_diaeiv 
    70       v_msf_eiv       ,  &  !: bolus "meridional" Stream-Function 
    71 #endif 
    72       surf_jk_r             !: inverse of the ocean "zonal" section surface 
     59   REAL(wp), DIMENSION(jpj,jpk) ::   tn_jk  , sn_jk  ,  &  !: "zonal" mean temperature and salinity 
     60      &                              v_msf_atl       ,  &  !: "meridional" Stream-Function 
     61      &                              v_msf_glo       ,  &  !: "meridional" Stream-Function 
     62      &                              v_msf_ipc       ,  &  !: "meridional" Stream-Function 
     63      &                              surf_jk_r             !: inverse of the ocean "zonal" section surface 
     64#if defined key_diaeiv 
     65   REAL(wp), DIMENSION(jpj,jpk) ::   v_msf_eiv                  !: bolus "meridional" Stream-Function 
     66#endif 
     67   REAL(wp), DIMENSION(jpi,jpj) ::   abasin, pbasin, ibasin     !: return function value 
    7368 
    7469   !! * Substitutions 
     
    7873   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    7974   !! $Header$  
    80    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     75   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    8176   !!---------------------------------------------------------------------- 
    8277 
     
    9489      !! 
    9590      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    96       !! 
    97       !!---------------------------------------------------------------------- 
    98       !! * arguments 
    99       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   & 
    100          pva                         ! mask flux array at V-point 
    101  
    102       !! * local declarations 
    103       INTEGER  ::   ji, jj, jk        ! dummy loop arguments 
    104       INTEGER  ::   ijpj             ! ??? 
    105       REAL(wp),DIMENSION(jpj) ::   & 
    106          p_fval                       ! function value 
     91      !!---------------------------------------------------------------------- 
     92      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point 
     93      !! 
     94      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments 
     95      INTEGER                  ::   ijpj         ! ??? 
     96      REAL(wp), DIMENSION(jpj) ::   p_fval       ! function value 
    10797      !!-------------------------------------------------------------------- 
    108  
     98      ! 
    10999      ijpj = jpj 
    110100      p_fval(:) = 0.e0 
     
    116106         END DO 
    117107      END DO 
    118  
     108      ! 
    119109      IF( lk_mpp )   CALL mpp_sum( p_fval, ijpj )     !!bug  I presume 
    120  
     110      ! 
    121111   END FUNCTION ptr_vj_3d 
    122  
    123112 
    124113 
     
    134123      !! 
    135124      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    136       !! 
    137       !!---------------------------------------------------------------------- 
    138       !! * arguments 
    139       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   & 
    140          pva                         ! mask flux array at V-point 
    141  
    142       !! * local declarations 
    143       INTEGER  ::   ji,jj             ! dummy loop arguments 
    144       INTEGER  ::   ijpj             ! ??? 
    145       REAL(wp),DIMENSION(jpj) ::   & 
    146          p_fval                       ! function value 
     125      !!---------------------------------------------------------------------- 
     126      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pva   ! mask flux array at V-point 
     127      !! 
     128      INTEGER                  ::   ji,jj    ! dummy loop arguments 
     129      INTEGER                  ::   ijpj     ! ??? 
     130      REAL(wp), DIMENSION(jpj) ::   p_fval   ! function value 
    147131      !!-------------------------------------------------------------------- 
    148        
     132      !  
    149133      ijpj = jpj 
    150134      p_fval(:) = 0.e0 
     
    154138         END DO 
    155139      END DO 
    156  
     140      ! 
    157141      IF( lk_mpp )   CALL mpp_sum( p_fval, ijpj )     !!bug  I presume 
    158   
    159     END FUNCTION ptr_vj_2d 
    160  
     142      !  
     143   END FUNCTION ptr_vj_2d 
    161144 
    162145 
     
    171154      !! 
    172155      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    173       !! 
    174       !!---------------------------------------------------------------------- 
    175       !! * arguments 
    176       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   & 
    177          pva                         ! mask flux array at V-point 
    178  
    179       !! * local declarations 
    180       INTEGER  ::   ji, jj, jk        ! dummy loop arguments 
    181       INTEGER, DIMENSION (1) :: ish 
    182       INTEGER, DIMENSION (2) :: ish2 
    183       REAL(wp),DIMENSION(jpj*jpk) ::   & 
    184          zwork                        ! temporary vector for mpp_sum 
    185       REAL(wp),DIMENSION(jpj,jpk) ::   & 
    186          p_fval                       ! return function value 
     156      !!---------------------------------------------------------------------- 
     157      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point 
     158      !! 
     159      INTEGER                      ::   ji, jj, jk   ! dummy loop arguments 
     160      INTEGER , DIMENSION (1)      ::   ish 
     161      INTEGER , DIMENSION (2)      ::   ish2 
     162      REAL(wp), DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum 
     163      REAL(wp), DIMENSION(jpj,jpk) ::   p_fval       ! return function value 
    187164      !!-------------------------------------------------------------------- 
    188   
     165      !  
    189166      p_fval(:,:) = 0.e0 
    190  
     167      ! 
    191168      DO jk = 1, jpkm1 
    192169         DO jj = 2, jpjm1 
     
    197174         END DO 
    198175      END DO 
    199  
     176      ! 
    200177      IF(lk_mpp) THEN 
    201178         ish(1) = jpj*jpk  ;  ish2(1) = jpj  ;  ish2(2) = jpk 
     
    204181         p_fval(:,:)= RESHAPE( zwork, ish2 ) 
    205182      END IF 
    206  
     183      ! 
    207184   END FUNCTION ptr_vjk 
     185 
    208186 
    209187   FUNCTION ptr_vtjk( pva )   RESULT ( p_fval ) 
     
    218196      !! 
    219197      !! ** Action  : - p_fval: i-k-mean poleward flux of pva 
    220       !! 
    221       !!---------------------------------------------------------------------- 
    222       !! * arguments 
    223       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   & 
    224          pva                         ! mask flux array at V-point 
    225   
    226       !! * local declarations 
    227       INTEGER  ::   ji, jj, jk        ! dummy loop arguments 
    228       INTEGER, DIMENSION (1) :: ish 
    229       INTEGER, DIMENSION (2) :: ish2 
    230       REAL(wp),DIMENSION(jpj*jpk) ::   & 
    231          zwork                        ! temporary vector for mpp_sum 
    232       REAL(wp),DIMENSION(jpj,jpk) ::   & 
    233          p_fval                       ! return function value 
     198      !!---------------------------------------------------------------------- 
     199      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point 
     200      !! 
     201      INTEGER                     ::   ji, jj, jk   ! dummy loop arguments 
     202      INTEGER, DIMENSION (1)      ::   ish 
     203      INTEGER, DIMENSION (2)      ::   ish2 
     204      REAL(wp),DIMENSION(jpj*jpk) ::   zwork        ! temporary vector for mpp_sum 
     205      REAL(wp),DIMENSION(jpj,jpk) ::   p_fval       ! return function value 
    234206      !!--------------------------------------------------------------------  
    235  
     207      ! 
    236208      p_fval(:,:) = 0.e0 
    237209      DO jk = 1, jpkm1 
     
    251223         p_fval(:,:)= RESHAPE(zwork,ish2) 
    252224      END IF 
    253  
     225      ! 
    254226   END FUNCTION ptr_vtjk 
    255227 
     
    259231      !!                  ***  ROUTINE dia_ptr  *** 
    260232      !!---------------------------------------------------------------------- 
    261       !! * Moudules used 
    262       USE ioipsl 
    263  
    264       !! * Argument 
    265233      INTEGER, INTENT(in) ::   kt   ! ocean time step index 
    266  
    267       !! * Local variables 
    268       INTEGER ::   jk,jj,ji               ! dummy loop 
    269       REAL(wp) ::    & 
    270          zsverdrup,  &              ! conversion from m3/s to Sverdrup 
    271          zpwatt,     &              ! conversion from W    to PW 
    272          zggram                     ! conversion from g    to Pg 
     234      !! 
     235      INTEGER  ::   jk, jj, ji               ! dummy loop 
     236      REAL(wp) ::   zsverdrup,  &              ! conversion from m3/s to Sverdrup 
     237         &          zpwatt,     &              ! conversion from W    to PW 
     238         &          zggram                     ! conversion from g    to Pg 
    273239 
    274240      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  & 
     
    277243         vs_atl, vs_pac, vs_ind,           & 
    278244         zv_eiv 
    279       CHARACTER (len=32) ::   & 
    280          clnam = 'subbasins.nc'                 
    281       INTEGER ::  itime,inum,ipi,ipj,ipk       ! temporary integer 
    282       INTEGER, DIMENSION (1) ::   istep 
    283       REAL(wp) ::    zdate0,zsecond,zdt        ! temporary scalars 
    284       REAL(wp), DIMENSION(jpidta,jpjdta) ::   & 
    285          zlamt, zphit, zdta             ! temporary workspace (NetCDF read) 
    286       REAL(wp), DIMENSION(jpk) ::   & 
    287          zdept                          ! temporary workspace (NetCDF read) 
     245      INTEGER ::  inum       ! temporary logical unit 
    288246      !!---------------------------------------------------------------------- 
    289247 
     
    293251         zpwatt    = 1.e-15 
    294252         zggram    = 1.e-6 
    295          ipi       = jpidta 
    296          ipj       = jpjdta 
    297          ipk       = 1 
    298          itime     = 1 
    299          zsecond   = 0.e0 
    300          zdate0    = 0.e0 
    301253    
    302254# if defined key_diaeiv 
     
    315267         IF( ln_subbas ) THEN              ! Basins computation 
    316268 
    317             IF( kt == nit000 ) THEN                ! load basin mask 
    318                itime = 1 
    319                ipi   = jpidta 
    320                ipj   = jpjdta 
    321                ipk   = 1 
    322                zdt   = 0.e0 
    323                istep = 0 
    324                clnam = 'subbasins.nc' 
    325  
    326                CALL flinopen(clnam,1,jpidta,1,jpjdta,.FALSE.,ipi,ipj, & 
    327                   &          ipk,zlamt,zphit,zdept,itime,istep,zdate0,zdt,inum) 
    328  
    329                ! get basins: 
    330                abasin (:,:) = 0.e0 
    331                pbasin (:,:) = 0.e0 
    332                ibasin (:,:) = 0.e0 
    333  
    334                ! Atlantic basin 
    335                CALL flinget(inum,'atlmsk',jpidta,jpjdta,1,itime,1,   & 
    336                   &         0,1,jpidta,1,jpjdta,zdta(:,:)) 
    337                DO jj = 1, nlcj                                 ! interior values 
    338                   DO ji = 1, nlci 
    339                      abasin (ji,jj) = zdta( mig(ji), mjg(jj) ) 
    340                   END DO 
    341                END DO 
    342  
    343                ! Pacific basin 
    344                CALL flinget(inum,'pacmsk',jpidta,jpjdta,1,itime,1,   & 
    345                   &         0,1,jpidta,1,jpjdta,zdta(:,:)) 
    346                DO jj = 1, nlcj                                 ! interior values 
    347                   DO ji = 1, nlci 
    348                      pbasin (ji,jj) = zdta( mig(ji), mjg(jj) ) 
    349                   END DO 
    350                END DO 
    351  
    352                ! Indian basin 
    353                CALL flinget(inum,'indmsk',jpidta,jpjdta,1,itime,1,   & 
    354                   &         0,1,jpidta,1,jpjdta,zdta(:,:)) 
    355                DO jj = 1, nlcj                                 ! interior values 
    356                   DO ji = 1, nlci 
    357                      ibasin (ji,jj) = zdta( mig(ji), mjg(jj) ) 
    358                   END DO 
    359                END DO 
    360  
    361                CALL flinclo(inum) 
    362  
     269            IF( kt == nit000 ) THEN                ! load sub-basin mask 
     270               CALL iom_open( 'subbasins', inum ) 
     271               CALL iom_get( inum, jpdom_data, 'atlmsk', abasin )      ! Atlantic basin 
     272               CALL iom_get( inum, jpdom_data, 'pacmsk', pbasin )      ! Pacific basin 
     273               CALL iom_get( inum, jpdom_data, 'indmsk', ibasin )      ! Indian basin 
     274               CALL iom_close( inum ) 
    363275            ENDIF 
    364276 
     
    396308#endif 
    397309         IF( ln_subbas ) THEN 
    398             v_msf_atl(:,:) = ptr_vjk( v_atl(:,:,:) )  
    399             v_msf_ipc(:,:) = ptr_vjk( v_ipc(:,:,:) )  
    400             ht_atl(:) = SUM(ptr_vjk( vt_atl(:,:,:)),2 ) 
    401             ht_pac(:) = SUM(ptr_vjk( vt_pac(:,:,:)),2 ) 
    402             ht_ind(:) = SUM(ptr_vjk( vt_ind(:,:,:)),2 ) 
    403             st_atl(:) = SUM(ptr_vjk( vs_atl(:,:,:)),2 ) 
    404             st_pac(:) = SUM(ptr_vjk( vs_pac(:,:,:)),2 ) 
    405             st_ind(:) = SUM(ptr_vjk( vs_ind(:,:,:)),2 ) 
     310            v_msf_atl(:,:) = ptr_vjk( v_atl (:,:,:) )  
     311            v_msf_ipc(:,:) = ptr_vjk( v_ipc (:,:,:) )  
     312            ht_atl(:) = SUM( ptr_vjk( vt_atl(:,:,:)), 2 ) 
     313            ht_pac(:) = SUM( ptr_vjk( vt_pac(:,:,:)), 2 ) 
     314            ht_ind(:) = SUM( ptr_vjk( vt_ind(:,:,:)), 2 ) 
     315            st_atl(:) = SUM( ptr_vjk( vs_atl(:,:,:)), 2 ) 
     316            st_pac(:) = SUM( ptr_vjk( vs_pac(:,:,:)), 2 ) 
     317            st_ind(:) = SUM( ptr_vjk( vs_ind(:,:,:)), 2 ) 
    406318         ENDIF 
    407319 
     
    466378      ! Close the file 
    467379      IF( kt == nitend ) CALL histclo( numptr ) 
    468  
     380      ! 
    469381   END SUBROUTINE dia_ptr 
    470382 
     
    475387      !!                    
    476388      !! ** Purpose :   Initialization, namelist read 
    477       !! 
    478389      !!---------------------------------------------------------------------- 
    479390      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   z_1         ! temporary workspace 
     
    485396      REWIND ( numnam ) 
    486397      READ   ( numnam, namptr ) 
    487  
    488398 
    489399      ! Control print 
     
    513423      !! 
    514424      !! ** Method  :   NetCDF file 
    515       !! 
    516       !!---------------------------------------------------------------------- 
    517       !! * Arguments 
     425      !!---------------------------------------------------------------------- 
    518426      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    519  
    520       !! * Save variables    
     427      !! 
    521428      INTEGER, SAVE ::   nhoridz, ndepidzt, ndepidzw, ndex(1) 
    522429 
    523       !! * Local variables 
    524       CHARACTER (len=40) ::   & 
    525          clhstnam, clop             ! temporary names 
    526       INTEGER ::   iline, it, ji    ! 
    527       REAL(wp) ::   & 
    528          zsto, zout, zdt, zmax, &   ! temporary scalars 
    529          zjulian 
     430      CHARACTER (len=40)       ::   clhstnam, clop                   ! temporary names 
     431      INTEGER                  ::   iline, it, ji                    ! 
     432      REAL(wp)                 ::   zsto, zout, zdt, zmax, zjulian   ! temporary scalars 
    530433      REAL(wp), DIMENSION(jpj) ::   zphi, zfoo 
    531434      !!---------------------------------------------------------------------- 
     
    720623  
    721624      ENDIF 
    722  
     625      ! 
    723626   END SUBROUTINE dia_ptr_wri 
    724627 
  • trunk/NEMO/OPA_SRC/DYN/dynspg_flt.F90

    r474 r508  
    44   !! Ocean dynamics:  surface pressure gradient trend 
    55   !!====================================================================== 
    6 #if ( defined key_dynspg_flt && ! defined key_mpp_omp )   ||   defined key_esopa 
     6   !! History    8.0  !  98-05  (G. Roullet)  free surface 
     7   !!                 !  98-10  (G. Madec, M. Imbard)  release 8.2 
     8   !!            8.5  !  02-08  (G. Madec)  F90: Free form and module 
     9   !!            " "  !  02-11  (C. Talandier, A-M Treguier) Open boundaries 
     10   !!            9.0  !  04-08  (C. Talandier) New trends organization 
     11   !!            " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
     12   !!            " "  !  06-07  (S. Masson)  distributed restart using iom 
     13   !!---------------------------------------------------------------------- 
     14#if ( defined key_dynspg_flt  && ! defined key_mpp_omp )  ||   defined key_esopa   
    715   !!---------------------------------------------------------------------- 
    816   !!   'key_dynspg_flt'                              filtered free surface 
    917   !!   NOT 'key_mpp_omp'                          k-j-i loop (vector opt.) 
     18   !!---------------------------------------------------------------------- 
    1019   !!---------------------------------------------------------------------- 
    1120   !!   dyn_spg_flt  : update the momentum trend with the surface pressure 
    1221   !!                  gradient in the filtered free surface case with 
    1322   !!                  vector optimization 
    14    !!---------------------------------------------------------------------- 
    15    !! * Modules used 
     23   !!   flt_rst      : read/write the time-splitting restart fields in the ocean restart file 
     24   !!---------------------------------------------------------------------- 
    1625   USE oce             ! ocean dynamics and tracers  
    1726   USE dom_oce         ! ocean space and time domain  
     
    2130   USE flxrnf          ! ocean runoffs 
    2231   USE sol_oce         ! ocean elliptic solver 
     32   USE solver          ! solver initialization 
    2333   USE solpcg          ! preconditionned conjugate gradient solver 
    2434   USE solsor          ! Successive Over-relaxation solver 
     
    3242   USE cla_dynspg      ! cross land advection 
    3343   USE prtctl          ! Print control 
    34    USE in_out_manager  ! I/O manager 
    3544   USE solmat          ! matrix construction for elliptic solvers 
    3645   USE agrif_opa_interp 
     46   USE in_out_manager  ! I/O manager 
     47   USE iom 
     48   USE restart         ! only for lrst_oce 
    3749 
    3850   IMPLICIT NONE 
    3951   PRIVATE 
    4052 
    41    !! * Accessibility 
    4253   PUBLIC dyn_spg_flt  ! routine called by step.F90 
    4354 
     
    4859   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    4960   !! $Header$  
    50    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     61   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    5162   !!---------------------------------------------------------------------- 
    5263 
     
    96107      !! ** Action : - Update (ua,va) with the surf. pressure gradient trend 
    97108      !! 
    98       !! References : 
    99       !!      Roullet and Madec 1999, JGR. 
    100       !! 
    101       !! History : 
    102       !!        !  98-05 (G. Roullet)  Original code 
    103       !!        !  98-10 (G. Madec, M. Imbard)  release 8.2 
    104       !!   8.5  !  02-08 (G. Madec)  F90: Free form and module 
    105       !!        !  02-11 (C. Talandier, A-M Treguier) Open boundaries 
    106       !!   9.0  !  04-08 (C. Talandier) New trends organization 
    107       !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization 
     109      !! References : Roullet and Madec 1999, JGR. 
    108110      !!--------------------------------------------------------------------- 
    109       !! * Arguments 
    110111      INTEGER, INTENT( in )  ::   kt         ! ocean time-step index 
    111       INTEGER, INTENT( out ) ::   kindic     ! solver convergence flag 
    112                                              ! if the solver doesn t converge 
    113                                              ! the flag is < 0 
    114       !! * Local declarations 
    115       INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    116       REAL(wp) ::                         &  
    117          z2dt, z2dtg, zraur, znugdt,      &  ! temporary scalars 
    118          znurau, zssha, zgcb, zbtd,       &  !   "          " 
    119          ztdgu, ztdgv                        !   "          " 
     112      INTEGER, INTENT( out ) ::   kindic     ! solver convergence flag (<0 if not converge) 
     113      !!                                    
     114      INTEGER  ::   ji, jj, jk                          ! dummy loop indices 
     115      REAL(wp) ::   z2dt, z2dtg, zraur, znugdt,      &  ! temporary scalars 
     116         &          znurau, zssha, zgcb, zbtd,       &  !   "          " 
     117         &          ztdgu, ztdgv                        !   "          " 
    120118      !!---------------------------------------------------------------------- 
    121  
     119      ! 
    122120      IF( kt == nit000 ) THEN 
    123121         IF(lwp) WRITE(numout,*) 
     
    128126         spgu(:,:) = 0.e0                     ! surface pressure gradient (i-direction) 
    129127         spgv(:,:) = 0.e0                     ! surface pressure gradient (j-direction) 
     128         CALL solver_init( nit000 )           ! Elliptic solver initialisation 
     129 
     130         ! read filtered free surface arrays in restart file 
     131         CALL flt_rst( nit000, 'READ' )       ! read or initialize the following fields: 
     132         !                                    ! gcx, gcxb, sshb, sshn 
    130133      ENDIF 
    131134 
     
    168171 
    169172#if defined key_obc 
    170       ! Update velocities on each open boundary with the radiation algorithm 
    171       CALL obc_dyn( kt ) 
    172       ! Correction of the barotropic componant velocity to control the volume of the system 
    173       CALL obc_vol( kt ) 
     173      CALL obc_dyn( kt )      ! Update velocities on each open boundary with the radiation algorithm 
     174      CALL obc_vol( kt )      ! Correction of the barotropic componant velocity to control the volume of the system 
    174175#endif 
    175176#if defined key_agrif 
    176       ! Update velocities on each coarse/fine interfaces 
    177  
    178       CALL Agrif_dyn( kt ) 
     177      CALL Agrif_dyn( kt )    ! Update velocities on each coarse/fine interfaces  
    179178#endif 
    180179#if defined key_orca_r2 
     
    243242      IF( .NOT. AGRIF_ROOT() ) THEN 
    244243         ! add contribution of gradient of after barotropic transport divergence  
    245          IF( (nbondi == -1) .OR. (nbondi == 2) ) gcb(3,:) = gcb(3,:) & 
    246             &            -znugdt * z2dt*laplacu(2,:)*gcdprc(3,:)*hu(2,:)*e2u(2,:) 
    247          IF( (nbondi ==  1) .OR. (nbondi == 2) )  gcb(nlci-2,:) = gcb(nlci-2,:) & 
    248             &           +znugdt * z2dt*laplacu(nlci-2,:)*gcdprc(nlci-2,:)*hu(nlci-2,:)*e2u(nlci-2,:) 
    249          IF( (nbondj == -1) .OR. (nbondj == 2) ) gcb(:,3) = gcb(:,3) & 
    250             &           -znugdt * z2dt*laplacv(:,2)*gcdprc(:,3)*hv(:,2)*e1v(:,2) 
    251          IF( (nbondj ==  1) .OR. (nbondj == 2) )  gcb(:,nlcj-2) = gcb(:,nlcj-2) & 
    252             &           +znugdt * z2dt*laplacv(:,nlcj-2)*gcdprc(:,nlcj-2)*hv(:,nlcj-2)*e1v(:,nlcj-2) 
     244         IF( nbondi == -1 .OR. nbondi == 2 )   gcb(3     ,:) =  & 
     245            &    gcb(3     ,:) - znugdt * z2dt * laplacu(2     ,:) * gcdprc(3     ,:) * hu(2     ,:) * e2u(2     ,:) 
     246         IF( nbondi ==  1 .OR. nbondi == 2 )   gcb(nlci-2,:) =  & 
     247            &    gcb(nlci-2,:) + znugdt * z2dt * laplacu(nlci-2,:) * gcdprc(nlci-2,:) * hu(nlci-2,:) * e2u(nlci-2,:) 
     248         IF( nbondj == -1 .OR. nbondj == 2 )   gcb(:     ,3) =  & 
     249            &    gcb(:,3     ) - znugdt * z2dt * laplacv(:,2     ) * gcdprc(:,3     ) * hv(:,2     ) * e1v(:,2     ) 
     250         IF( nbondj ==  1 .OR. nbondj == 2 )   gcb(:,nlcj-2) =  & 
     251            &    gcb(:,nlcj-2) + znugdt * z2dt * laplacv(:,nlcj-2) * gcdprc(:,nlcj-2) * hv(:,nlcj-2) * e1v(:,nlcj-2) 
    253252      ENDIF 
    254253#endif 
     
    263262      epsr = eps * eps * rnorme 
    264263      ncut = 0 
    265       ! if rnorme is 0, the solution is 0, the solver isn't called 
     264      ! if rnorme is 0, the solution is 0, the solver is not called 
    266265      IF( rnorme == 0.e0 ) THEN 
    267266         gcx(:,:) = 0.e0 
     
    313312      IF( .NOT. Agrif_Root() ) THEN 
    314313         ! caution : grad D (fine) = grad D (coarse) at coarse/fine interface 
    315          IF( (nbondi == -1) .OR. (nbondi == 2) ) spgu(2,:) = znugdt * z2dt * laplacu(2,:) * umask(2,:,1) 
    316          IF( (nbondi ==  1) .OR. (nbondi == 2) ) spgu(nlci-2,:) = znugdt * z2dt * laplacu(nlci-2,:) * umask(nlci-2,:,1) 
    317          IF( (nbondj == -1) .OR. (nbondj == 2) ) spgv(:,2) = znugdt * z2dt * laplacv(:,2) * vmask(:,2,1) 
    318          IF( (nbondj ==  1) .OR. (nbondj == 2) ) spgv(:,nlcj-2) = znugdt * z2dt * laplacv(:,nlcj-2) * vmask(:,nlcj-2,1) 
     314         IF( nbondi == -1 .OR. nbondi == 2 ) spgu(2     ,:) = znugdt * z2dt * laplacu(2     ,:) * umask(2     ,:,1) 
     315         IF( nbondi ==  1 .OR. nbondi == 2 ) spgu(nlci-2,:) = znugdt * z2dt * laplacu(nlci-2,:) * umask(nlci-2,:,1) 
     316         IF( nbondj == -1 .OR. nbondj == 2 ) spgv(:,2     ) = znugdt * z2dt * laplacv(:,2     ) * vmask(:     ,2,1) 
     317         IF( nbondj ==  1 .OR. nbondj == 2 ) spgv(:,nlcj-2) = znugdt * z2dt * laplacv(:,nlcj-2) * vmask(:,nlcj-2,1) 
    319318      ENDIF 
    320319#endif       
     
    323322      !     ( c a u t i o n : (ua,va) here are the after velocity not the 
    324323      !                       trend, the leap-frog time stepping will not 
    325       !                       be done in dynnxt.F routine) 
     324      !                       be done in dynnxt.F90 routine) 
    326325      DO jk = 1, jpkm1 
    327326         DO jj = 2, jpjm1 
     
    332331         END DO 
    333332      END DO 
    334  
    335333 
    336334      ! Sea surface elevation time stepping 
     
    358356      ENDIF 
    359357 
    360       !                       ! print sum trends (used for debugging) 
    361       IF(ln_ctl)   CALL prt_ctl( tab2d_1=sshn, clinfo1=' spg  - ssh: ', mask1=tmask ) 
    362  
     358      ! write filtered free surface arrays in restart file 
     359      ! -------------------------------------------------- 
     360      IF( lrst_oce ) CALL flt_rst( kt, 'WRITE' ) 
     361 
     362      ! print sum trends (used for debugging) 
     363      IF(ln_ctl)     CALL prt_ctl( tab2d_1=sshn, clinfo1=' spg  - ssh: ', mask1=tmask ) 
     364      ! 
    363365   END SUBROUTINE dyn_spg_flt 
     366 
     367 
     368   SUBROUTINE flt_rst( kt, cdrw ) 
     369     !!--------------------------------------------------------------------- 
     370     !!                   ***  ROUTINE ts_rst  *** 
     371     !! 
     372     !! ** Purpose : Read or write filtered free surface arrays in restart file 
     373     !!---------------------------------------------------------------------- 
     374     INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
     375     CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
     376     !!---------------------------------------------------------------------- 
     377 
     378     IF( TRIM(cdrw) == 'READ' ) THEN 
     379        IF( iom_varid( numror, 'gcx' ) > 0 ) THEN 
     380! Caution : extra-hallow 
     381! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     382           CALL iom_get( numror, jpdom_local, 'gcx' , gcx (1:jpi,1:jpj) ) 
     383           CALL iom_get( numror, jpdom_local, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     384           CALL iom_get( numror, jpdom_local, 'sshb', sshb(:,:)         ) 
     385           CALL iom_get( numror, jpdom_local, 'sshn', sshn(:,:)         ) 
     386           IF( neuler == 0 ) THEN 
     387              sshb(:,:) = sshn(:,:) 
     388              gcxb(:,:) = gcx (:,:) 
     389           ENDIF 
     390        ELSE 
     391           gcx (:,:) = 0.e0 
     392           gcxb(:,:) = 0.e0 
     393           sshb(:,:) = 0.e0 
     394           sshn(:,:) = 0.e0 
     395        ENDIF 
     396     ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
     397! Caution : extra-hallow 
     398! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     399        CALL iom_rstput( kt, nitrst, numrow, 'gcx' , gcx( 1:jpi,1:jpj) ) 
     400        CALL iom_rstput( kt, nitrst, numrow, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     401        CALL iom_rstput( kt, nitrst, numrow, 'sshb', sshb(:,:)         ) 
     402        CALL iom_rstput( kt, nitrst, numrow, 'sshn', sshn(:,:)         ) 
     403     ENDIF 
     404     ! 
     405   END SUBROUTINE flt_rst 
    364406 
    365407#else 
  • trunk/NEMO/OPA_SRC/DYN/dynspg_flt_jki.F90

    r503 r508  
    2525   USE solfet          ! FETI solver 
    2626   USE solsor_e        ! Successive Over-relaxation solver with MPP optimization 
     27   USE solver 
    2728   USE obc_oce         ! Lateral open boundary condition 
    2829   USE obcdyn          ! ocean open boundary condition (obc_dyn routines) 
     
    3536   USE solmat          ! matrix construction for elliptic solvers 
    3637   USE agrif_opa_interp 
     38   USE restart         ! only for lrst_oce 
     39   USE iom 
    3740 
    3841   IMPLICIT NONE 
     
    112115         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~  (free surface constant volume, autotasking case)' 
    113116 
    114          ! set to zero free surface specific arrays  
    115          spgu(:,:) = 0.e0      ! surface pressure gradient (i-direction) 
    116          spgv(:,:) = 0.e0      ! surface pressure gradient (j-direction) 
     117         ! set to zero free surface specific arrays 
     118         spgu(:,:) = 0.e0                     ! surface pressure gradient (i-direction) 
     119          
     120         spgv(:,:) = 0.e0                     ! surface pressure gradient (j-direction) 
     121          
     122         CALL solver_init( nit000 )           ! Elliptic solver initialisation 
     123 
     124         ! read filtered free surface arrays in restart file 
     125         CALL flt_rst( nit000, 'READ' )       ! read or initialize the following fields: 
     126         !                                    ! gcx, gcxb, sshb, sshn 
    117127      ENDIF 
    118128 
     
    354364      CALL lbc_lnk( sshn, 'T', 1. ) 
    355365 
     366      ! write filtered free surface arrays in restart file 
     367      ! -------------------------------------------------- 
     368      IF( lrst_oce ) CALL flt_rst( kt, 'WRITE' ) 
     369 
    356370      IF(ln_ctl) THEN         ! print sum trends (used for debugging) 
    357371         CALL prt_ctl( tab3d_1=ua  , clinfo1=' spg  - Ua : ', mask1=umask,   & 
     
    362376   END SUBROUTINE dyn_spg_flt_jki 
    363377 
     378   SUBROUTINE flt_rst( kt, cdrw ) 
     379     !!--------------------------------------------------------------------- 
     380     !!                   ***  ROUTINE ts_rst  *** 
     381     !! 
     382     !! ** Purpose : Read or write filtered free surface arrays in restart file 
     383     !!---------------------------------------------------------------------- 
     384     INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
     385     CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
     386     !!---------------------------------------------------------------------- 
     387 
     388     IF( TRIM(cdrw) == 'READ' ) THEN 
     389        IF( iom_varid( numror, 'gcx' ) > 0 ) THEN 
     390! Caution : extra-hallow 
     391! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     392           CALL iom_get( numror, jpdom_local, 'gcx' , gcx (1:jpi,1:jpj) ) 
     393           CALL iom_get( numror, jpdom_local, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     394           CALL iom_get( numror, jpdom_local, 'sshb', sshb(:,:)         ) 
     395           CALL iom_get( numror, jpdom_local, 'sshn', sshn(:,:)         ) 
     396           IF( neuler == 0 ) THEN 
     397              sshb(:,:) = sshn(:,:) 
     398              gcxb(:,:) = gcx (:,:) 
     399           ENDIF 
     400        ELSE 
     401           gcx (:,:) = 0.e0 
     402           gcxb(:,:) = 0.e0 
     403           sshb(:,:) = 0.e0 
     404           sshn(:,:) = 0.e0 
     405        ENDIF 
     406     ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
     407! Caution : extra-hallow 
     408! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     409        CALL iom_rstput( kt, nitrst, numrow, 'gcx' , gcx( 1:jpi,1:jpj) ) 
     410        CALL iom_rstput( kt, nitrst, numrow, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     411        CALL iom_rstput( kt, nitrst, numrow, 'sshb', sshb(:,:)         ) 
     412        CALL iom_rstput( kt, nitrst, numrow, 'sshn', sshn(:,:)         ) 
     413     ENDIF 
     414     ! 
     415   END SUBROUTINE flt_rst 
     416 
    364417#else 
    365418   !!---------------------------------------------------------------------- 
  • trunk/NEMO/OPA_SRC/DYN/dynspg_rl.F90

    r474 r508  
    44   !! Ocean dynamics:  surface pressure gradient trend 
    55   !!====================================================================== 
     6   !! History :  7.0  !  96-05  (G. Madec, M. Imbard, M. Guyon)  rewitting in 1 
     7   !!                           routine, without pointers, and with the same matrix 
     8   !!                           for sor and pcg, mpp exchange, and symmetric conditions 
     9   !!            " "  !  96-07  (A. Weaver)  Euler forward step 
     10   !!            " "  !  96-11  (A. Weaver)  correction to preconditioning: 
     11   !!            8.0  !  98-02  (M. Guyon)  FETI method 
     12   !!            " "  !  97-09  (J.-M. Molines)  Open boundaries 
     13   !!            8.5  !  02-08  (G. Madec)  F90: Free form and module 
     14   !!                 !  02-11  (C. Talandier, A-M Treguier) Open boundaries 
     15   !!            9.0  !  04-08  (C. Talandier)  New trends organization 
     16   !!            " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
     17   !!            " "  !  06-07  (S. Masson)  distributed restart using iom 
     18   !!--------------------------------------------------------------------- 
    619#if   defined key_dynspg_rl   ||   defined key_esopa 
    720   !!---------------------------------------------------------------------- 
     
    1023   !!   dyn_spg_rl   : update the momentum trend with the surface pressure 
    1124   !!                  for the rigid-lid case. 
    12    !!---------------------------------------------------------------------- 
    13    !! * Modules used 
     25   !!   rl_rst       : read/write the rigid-lid restart fields in the ocean restart file 
     26   !!---------------------------------------------------------------------- 
    1427   USE oce             ! ocean dynamics and tracers 
    1528   USE dom_oce         ! ocean space and time domain 
     
    1932   USE zdf_oce         ! ocean vertical physics 
    2033   USE sol_oce         ! ocean elliptic solver 
     34   USE solver          ! solver initialization 
    2135   USE solpcg          ! preconditionned conjugate gradient solver 
    2236   USE solsor          ! Successive Over-relaxation solver 
     
    2842   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2943   USE in_out_manager  ! I/O manager 
     44   USE iom 
     45   USE restart         ! only for lrst_oce 
    3046 
    3147   IMPLICIT NONE 
    3248   PRIVATE 
    3349 
    34    !! * Accessibility 
    3550   PUBLIC dyn_spg_rl   ! called by step.F90 
    3651 
     
    4257   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    4358   !! $Header$  
    44    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     59   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    4560   !!---------------------------------------------------------------------- 
    4661 
     
    7691      !! ** Action : - Update (ua,va) with the surf. pressure gradient trend 
    7792      !! 
    78       !! References : 
    79       !!      Madec et al. 1988, ocean modelling, issue 78, 1-6. 
    80       !! 
    81       !! History : 
    82       !!        !  96-05  (G. Madec, M. Imbard, M. Guyon)  rewitting in 1 
    83       !!                  routine, without pointers, and with the same matrix 
    84       !!                  for sor and pcg, mpp exchange, and symmetric conditions 
    85       !!        !  96-07  (A. Weaver)  Euler forward step 
    86       !!        !  96-11  (A. Weaver)  correction to preconditioning: 
    87       !!        !  98-02  (M. Guyon)  FETI method 
    88       !!        !  98-05  (G. Roullet)  free surface 
    89       !!        !  97-09  (J.-M. Molines)  Open boundaries 
    90       !!   8.5  !  02-08  (G. Madec)  F90: Free form and module 
    91       !!        !  02-11  (C. Talandier, A-M Treguier) Open boundaries 
    92       !!   9.0  !  04-08  (C. Talandier)  New trends organization 
    93       !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization 
     93      !! References : Madec et al. 1988, ocean modelling, issue 78, 1-6. 
    9494      !!--------------------------------------------------------------------- 
    95       !! * Arguments 
    9695      INTEGER, INTENT( in  ) ::   kt       ! ocean time-step index 
    97       INTEGER, INTENT( out ) ::   kindic   ! solver flag, take a negative value 
    98       !                                    ! when the solver doesnot converge 
    99       !! * Local declarations 
     96      INTEGER, INTENT( out ) ::   kindic   ! solver flag (<0 when the solver does not converge) 
    10097      INTEGER ::   ji, jj, jk    ! dummy loop indices 
    10198      REAL(wp) ::  zbsfa, zgcx, z2dt 
     
    114111 
    115112         ! set to zero rigid-lid specific arrays 
    116          spgu(:,:) = 0.e0      ! surface pressure gradient (i-direction)  
    117          spgv(:,:) = 0.e0      ! surface pressure gradient (j-direction) 
    118       ENDIF 
    119  
    120       ! 0. Initializations: 
    121       ! ------------------- 
    122 # if defined key_obc 
    123       ! space index on boundary arrays 
    124       ib = 1 
    125       ibm = 2 
    126       ibm2 = 3 
    127       ! time index on boundary arrays 
    128       it = 1 
    129       itm = 2 
    130       itm2 = 3 
    131 # endif 
    132  
     113         spgu(:,:) = 0.e0                   ! surface pressure gradient (i-direction)  
     114         spgv(:,:) = 0.e0                   ! surface pressure gradient (j-direction) 
     115 
     116         CALL solver_init( nit000 )         ! Elliptic solver initialisation 
     117 
     118         ! read rigid lid arrays in restart file 
     119         CALL rl_rst( nit000, 'READ' )      ! read or initialize the following fields: 
     120         !                                  ! gcx, gcxb, bsfb, bsfn, bsfd 
     121      ENDIF 
     122 
     123      !  Vertically averaged momentum trend 
     124      ! ------------------------------------ 
    133125      !                                                ! =============== 
    134126      DO jj = 2, jpjm1                                 !  Vertical slab 
    135127         !                                             ! =============== 
    136  
    137          ! 1. Vertically averaged momentum trend 
    138          ! ------------------------------------- 
    139          ! initialization to zero 
    140          spgu(:,jj) = 0. 
     128          
     129         spgu(:,jj) = 0.                          ! initialization to zero 
    141130         spgv(:,jj) = 0. 
    142          ! vertical sum 
    143          DO jk = 1, jpk 
     131         DO jk = 1, jpk                           ! vertical sum 
    144132            DO ji = 2, jpim1 
    145133               spgu(ji,jj) = spgu(ji,jj) + ua(ji,jj,jk) * fse3u(ji,jj,jk) * umask(ji,jj,jk) 
     
    147135            END DO  
    148136         END DO  
    149          ! divide by the depth 
    150          spgu(:,jj) = spgu(:,jj) * hur(:,jj) 
     137         spgu(:,jj) = spgu(:,jj) * hur(:,jj)      ! divide by the depth  
    151138         spgv(:,jj) = spgv(:,jj) * hvr(:,jj) 
    152139 
     
    155142      !                                                ! =============== 
    156143 
    157       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    158  
    159144      ! Boundary conditions on (spgu,spgv) 
    160145      CALL lbc_lnk( spgu, 'U', -1. ) 
    161146      CALL lbc_lnk( spgv, 'V', -1. ) 
    162147 
    163       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    164  
    165       ! 2. Barotropic streamfunction trend (bsfd) 
     148      !  Barotropic streamfunction trend (bsfd) 
    166149      ! ---------------------------------- 
    167  
    168       ! Curl of the vertically averaged velocity 
    169       DO jj = 2, jpjm1 
     150      DO jj = 2, jpjm1                            ! Curl of the vertically averaged velocity  
    170151         DO ji = 2, jpim1 
    171152            gcb(ji,jj) = -gcdprc(ji,jj)   & 
    172                        *(  ( e2v(ji+1,jj  )*spgv(ji+1,jj  ) - e2v(ji,jj)*spgv(ji,jj) )   & 
    173                           -( e1u(ji  ,jj+1)*spgu(ji  ,jj+1) - e1u(ji,jj)*spgu(ji,jj) )   )  
     153               &       *(  ( e2v(ji+1,jj  )*spgv(ji+1,jj  ) - e2v(ji,jj)*spgv(ji,jj) )   & 
     154               &          -( e1u(ji  ,jj+1)*spgu(ji  ,jj+1) - e1u(ji,jj)*spgu(ji,jj) )   )  
    174155         END DO 
    175156      END DO 
    176157 
    177158# if defined key_obc 
    178       ! Open boundary contribution 
    179       DO jj = 2, jpjm1 
     159      DO jj = 2, jpjm1                            ! Open boundary contribution  
    180160         DO ji = 2, jpim1 
    181161           gcb(ji,jj) = gcb(ji,jj) - gcdprc(ji,jj) * gcbob(ji,jj) 
     
    198178      ! applied the lateral boundary conditions 
    199179      IF( nsolv == 4)   CALL lbc_lnk_e( gcb, c_solver_pt, 1. )    
    200  
    201       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    202180 
    203181      ! Relative precision (computation on one processor) 
     
    234212      ENDIF 
    235213 
    236  
    237214      ! bsf trend update 
    238215      ! ---------------- 
    239  
    240216      bsfd(1:nlci,1:nlcj) = gcx(1:nlci,1:nlcj) 
    241  
    242217       
    243218      ! update bsf trend with islands trend 
    244219      ! ----------------------------------- 
    245  
    246220      IF( lk_isl )   CALL isl_dyn_spg         ! update bsfd 
    247  
    248221 
    249222# if defined key_obc 
     
    337310      ! 4. Barotropic stream function and array swap 
    338311      ! -------------------------------------------- 
    339  
    340312      ! Leap-frog time scheme, time filter and array swap 
    341313      IF( neuler == 0 .AND. kt == nit000 ) THEN 
     
    362334 
    363335# if defined key_obc 
     336      ib   = 1      ! space index on boundary arrays 
     337      ibm  = 2 
     338      ibm2 = 3 
     339      it   = 1      ! time index on boundary arrays 
     340      itm  = 2 
     341      itm2 = 3 
     342 
    364343      ! Swap of boundary arrays 
    365344      IF( lp_obc_east ) THEN 
     
    499478      ENDIF 
    500479# endif 
    501       ! 
    502       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    503480       
    504481      !  add the surface pressure trend to the general trend 
    505482      ! ----------------------------------------------------- 
    506        
    507483      DO jj = 2, jpjm1 
    508  
    509484         ! update the surface pressure gradient with the barotropic trend 
    510485         DO ji = 2, jpim1 
     
    519494            END DO 
    520495         END DO 
    521  
    522       END DO 
     496      END DO 
     497 
     498      ! write rigid lid arrays in restart file 
     499      ! -------------------------------------- 
     500      IF( lrst_oce ) CALL rl_rst( kt, 'WRITE' ) 
    523501 
    524502   END SUBROUTINE dyn_spg_rl 
     503 
     504 
     505   SUBROUTINE rl_rst( kt, cdrw ) 
     506     !!--------------------------------------------------------------------- 
     507     !!                   ***  ROUTINE rl_rst  *** 
     508     !! 
     509     !! ** Purpose : Read or write rigid-lid arrays in ocean restart file 
     510     !!---------------------------------------------------------------------- 
     511     INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
     512     CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
     513     !!---------------------------------------------------------------------- 
     514     ! 
     515     IF( TRIM(cdrw) == 'READ' ) THEN 
     516        IF( iom_varid( numror, 'gcx' ) > 0 ) THEN 
     517     ! Caution : extra-hallow 
     518     ! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     519           CALL iom_get( numror, jpdom_local, 'gcx' , gcx (1:jpi,1:jpj) ) 
     520           CALL iom_get( numror, jpdom_local, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     521           CALL iom_get( numror, jpdom_local, 'bsfb', bsfb(:,:)         ) 
     522           CALL iom_get( numror, jpdom_local, 'bsfn', bsfn(:,:)         ) 
     523           CALL iom_get( numror, jpdom_local, 'bsfd', bsfd(:,:)         ) 
     524           IF( neuler == 0 ) THEN 
     525              gcxb(:,:) = gcx (:,:) 
     526              bsfb(:,:) = bsfn(:,:) 
     527           ENDIF 
     528        ELSE 
     529           gcx (:,:) = 0.e0 
     530           gcxb(:,:) = 0.e0 
     531           bsfb(:,:) = 0.e0 
     532           bsfn(:,:) = 0.e0 
     533           bsfd(:,:) = 0.e0 
     534        ENDIF 
     535     ELSEIF(  TRIM(cdrw) == 'WRITE' ) THEN 
     536        ! Caution : extra-hallow, gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
     537        CALL iom_rstput( kt, nitrst, numrow, 'gcx' , gcx (1:jpi,1:jpj) ) 
     538        CALL iom_rstput( kt, nitrst, numrow, 'gcxb', gcxb(1:jpi,1:jpj) ) 
     539        CALL iom_rstput( kt, nitrst, numrow, 'bsfb', bsfb(:,:)         ) 
     540        CALL iom_rstput( kt, nitrst, numrow, 'bsfn', bsfn(:,:)         ) 
     541        CALL iom_rstput( kt, nitrst, numrow, 'bsfd', bsfd(:,:)         ) 
     542     ENDIF 
     543     ! 
     544   END SUBROUTINE rl_rst 
    525545 
    526546#else 
  • trunk/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r455 r508  
    44   !! Ocean dynamics:  surface pressure gradient trend 
    55   !!====================================================================== 
     6   !! History :   9.0  !  04-12  (L. Bessieres, G. Madec)  Original code 
     7   !!             " "  !  05-11  (V. Garnier, G. Madec)  optimization 
     8   !!             9.0  !  06-08  (S. Masson)  distributed restart using iom 
     9   !!--------------------------------------------------------------------- 
    610#if ( defined key_dynspg_ts && ! defined key_mpp_omp ) ||   defined key_esopa 
    711   !!---------------------------------------------------------------------- 
     
    913   !!   NOT 'key_mpp_omp'                          k-j-i loop (vector opt.) 
    1014   !!---------------------------------------------------------------------- 
     15   !!---------------------------------------------------------------------- 
    1116   !!   dyn_spg_ts  : compute surface pressure gradient trend using a time- 
    1217   !!                 splitting scheme and add to the general trend  
     18   !!   ts_rst      : read/write the time-splitting restart fields in the ocean restart file 
    1319   !!---------------------------------------------------------------------- 
    1420   !! * Modules used 
     
    2733   USE dynspg_oce      ! surface pressure gradient variables 
    2834   USE in_out_manager  ! I/O manager 
     35   USE iom 
     36   USE restart         ! only for lrst_oce 
    2937 
    3038   IMPLICIT NONE 
    3139   PRIVATE 
    3240 
    33    !! * Accessibility 
    3441   PUBLIC dyn_spg_ts  ! routine called by step.F90 
     42 
     43   REAL(wp), DIMENSION(jpi,jpj) ::  ftnw, ftne,   &  ! triad of coriolis parameter 
     44      &                             ftsw, ftse       ! (only used with een vorticity scheme) 
     45 
    3546 
    3647   !! * Substitutions 
     
    7485      !! ** Action : - Update (ua,va) with the surf. pressure gradient trend 
    7586      !! 
    76       !! References : 
    77       !!   Griffies et al., (2003): A technical guide to MOM4. NOAA/GFDL 
    78       !! 
    79       !! History : 
    80       !!   9.0  !  04-12  (L. Bessieres, G. Madec)  Original code 
    81       !!        !  05-11  (V. Garnier, G. Madec)  optimization 
     87      !! References : Griffies et al., (2003): A technical guide to MOM4. NOAA/GFDL 
    8288      !!--------------------------------------------------------------------- 
    83       !! * Arguments 
    8489      INTEGER, INTENT( in )  ::   kt           ! ocean time-step index 
    8590 
     
    97102         zsshb_e, zub_e, zvb_e,             &  !     "        " 
    98103         zun_e, zvn_e                          !     "        " 
    99       REAL(wp), DIMENSION(jpi,jpj),SAVE ::  & 
    100          ztnw, ztne, ztsw, ztse 
    101104      !!---------------------------------------------------------------------- 
    102105 
     
    109112 
    110113      IF( kt == nit000 ) THEN 
    111  
     114         ! 
    112115         IF(lwp) WRITE(numout,*) 
    113116         IF(lwp) WRITE(numout,*) 'dyn_spg_ts : surface pressure gradient trend' 
     
    115118         IF(lwp) WRITE(numout,*) ' Number of sub cycle in 1 time-step (2 rdt) : icycle = ', FLOOR( 2*rdt/rdtbt ) 
    116119 
    117          IF( .NOT. ln_rstart ) THEN 
    118             ! initialize barotropic specific arrays 
    119             sshb_b(:,:) = sshb(:,:) 
    120             sshn_b(:,:) = sshn(:,:) 
    121             un_b(:,:)   = 0.e0 
    122             vn_b(:,:)   = 0.e0 
    123             ! vertical sum 
    124             IF( lk_vopt_loop ) THEN          ! vector opt., forced unroll 
    125                DO jk = 1, jpkm1 
    126                   DO ji = 1, jpij 
    127                      un_b(ji,1) = un_b(ji,1) + fse3u(ji,1,jk) * un(ji,1,jk) 
    128                      vn_b(ji,1) = vn_b(ji,1) + fse3v(ji,1,jk) * vn(ji,1,jk) 
    129                   END DO 
    130                END DO 
    131             ELSE                             ! No  vector opt. 
    132                DO jk = 1, jpkm1 
    133                   un_b(:,:) = un_b(:,:) + fse3u(:,:,jk) * un(:,:,jk) 
    134                   vn_b(:,:) = vn_b(:,:) + fse3v(:,:,jk) * vn(:,:,jk) 
    135                END DO 
    136             ENDIF 
    137          ENDIF 
     120         CALL ts_rst( nit000, 'READ' )   ! read or initialize the following fields: 
     121         !                               ! sshb, sshn, sshb_b, sshn_b, un_b, vn_b 
     122 
    138123         ssha_e(:,:) = sshn(:,:) 
    139124         ua_e(:,:)   = un_b(:,:) 
     
    141126 
    142127         IF( ln_dynvor_een ) THEN 
    143             ztne(1,:) = 0.e0   ;   ztnw(1,:) = 0.e0   ;   ztse(1,:) = 0.e0   ;   ztsw(1,:) = 0.e0 
     128            ftne(1,:) = 0.e0   ;   ftnw(1,:) = 0.e0   ;   ftse(1,:) = 0.e0   ;   ftsw(1,:) = 0.e0 
    144129            DO jj = 2, jpj 
    145130               DO ji = fs_2, jpi   ! vector opt. 
    146                   ztne(ji,jj) = ( ff(ji-1,jj  ) + ff(ji  ,jj  ) + ff(ji  ,jj-1) ) / 3. 
    147                   ztnw(ji,jj) = ( ff(ji-1,jj-1) + ff(ji-1,jj  ) + ff(ji  ,jj  ) ) / 3. 
    148                   ztse(ji,jj) = ( ff(ji  ,jj  ) + ff(ji  ,jj-1) + ff(ji-1,jj-1) ) / 3. 
    149                   ztsw(ji,jj) = ( ff(ji  ,jj-1) + ff(ji-1,jj-1) + ff(ji-1,jj  ) ) / 3. 
     131                  ftne(ji,jj) = ( ff(ji-1,jj  ) + ff(ji  ,jj  ) + ff(ji  ,jj-1) ) / 3. 
     132                  ftnw(ji,jj) = ( ff(ji-1,jj-1) + ff(ji-1,jj  ) + ff(ji  ,jj  ) ) / 3. 
     133                  ftse(ji,jj) = ( ff(ji  ,jj  ) + ff(ji  ,jj-1) + ff(ji-1,jj-1) ) / 3. 
     134                  ftsw(ji,jj) = ( ff(ji  ,jj-1) + ff(ji-1,jj-1) + ff(ji-1,jj  ) ) / 3. 
    150135               END DO 
    151136            END DO 
    152137         ENDIF 
    153  
     138         ! 
    154139      ENDIF 
    155      
     140 
    156141      ! Local constant initialization 
    157142      ! -------------------------------- 
     
    216201            END DO 
    217202         END DO 
    218  
     203         ! 
    219204      ELSEIF ( ln_dynvor_ens ) THEN                    ! enstrophy conserving scheme 
    220205         DO jj = 2, jpjm1 
     
    228213            END DO 
    229214         END DO 
    230  
     215         ! 
    231216      ELSEIF ( ln_dynvor_een ) THEN                    ! enstrophy and energy conserving scheme 
    232217         zfac25 = 0.25 
     
    241226            END DO 
    242227         END DO 
    243  
     228         ! 
    244229      ENDIF 
    245230 
     
    300285      DO jit = 1, icycle                                   !  sub-time-step loop  ! 
    301286         !                                                 ! ==================== ! 
    302  
    303287         z2dt_e = 2. * rdtbt 
    304288         IF ( jit == 1 )   z2dt_e = rdtbt 
     
    360344               END DO 
    361345            END DO 
    362  
     346            ! 
    363347         ELSEIF ( ln_dynvor_ens ) THEN                    ! enstrophy conserving scheme 
    364348            DO jj = 2, jpjm1 
     
    379363               END DO 
    380364            END DO 
    381  
     365            ! 
    382366         ELSEIF ( ln_dynvor_een ) THEN                    ! energy and enstrophy conserving scheme 
    383367            zfac25 = 0.25 
     
    397381               END DO 
    398382            END DO 
     383            !  
    399384         ENDIF 
    400385 
     
    504489      END DO 
    505490 
    506       IF(ln_ctl) THEN         ! print sum trends (used for debugging) 
    507          CALL prt_ctl(tab2d_1=sshn, clinfo1=' ssh      : ', mask1=tmask) 
     491      ! write filtered free surface arrays in restart file 
     492      ! -------------------------------------------------- 
     493      IF( lrst_oce )   CALL ts_rst( kt, 'WRITE' ) 
     494 
     495      ! print sum trends (used for debugging) 
     496      IF( ln_ctl )     CALL prt_ctl( tab2d_1=sshn, clinfo1=' ssh      : ', mask1=tmask ) 
     497      ! 
     498   END SUBROUTINE dyn_spg_ts 
     499 
     500 
     501   SUBROUTINE ts_rst( kt, cdrw ) 
     502      !!--------------------------------------------------------------------- 
     503      !!                   ***  ROUTINE ts_rst  *** 
     504      !! 
     505      !! ** Purpose : Read or write time-splitting arrays in restart file 
     506      !!---------------------------------------------------------------------- 
     507      INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
     508      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
     509      ! 
     510      INTEGER ::  ji, jk        ! dummy loop indices 
     511      !!---------------------------------------------------------------------- 
     512      ! 
     513      IF( TRIM(cdrw) == 'READ' ) THEN 
     514         IF( iom_varid( numror, 'sshn' ) > 0 ) THEN 
     515            CALL iom_get( numror, jpdom_local, 'sshb'  , sshb(:,:)   ) 
     516            CALL iom_get( numror, jpdom_local, 'sshn'  , sshn(:,:)   ) 
     517            IF( neuler == 0 ) sshb(:,:) = sshn(:,:) 
     518         ELSE 
     519            sshb(:,:) = 0.e0 
     520            sshn(:,:) = 0.e0 
     521         ENDIF 
     522         IF( iom_varid( numror, 'sshn_b' ) > 0 ) THEN 
     523            CALL iom_get( numror, jpdom_local, 'sshb_b', sshb_b(:,:) )   ! free surface issued 
     524            CALL iom_get( numror, jpdom_local, 'sshn_b', sshn_b(:,:) )   ! from time-splitting loop 
     525            CALL iom_get( numror, jpdom_local, 'un_b'  , un_b  (:,:) )   ! horizontal transports issued 
     526            CALL iom_get( numror, jpdom_local, 'vn_b'  , vn_b  (:,:) )   ! from barotropic loop 
     527            IF( neuler == 0 ) sshb_b(:,:) = sshn_b(:,:) 
     528         ELSE 
     529            sshb_b(:,:) = sshb(:,:) 
     530            sshn_b(:,:) = sshn(:,:) 
     531            un_b  (:,:) = 0.e0 
     532            vn_b  (:,:) = 0.e0 
     533            ! vertical sum 
     534            IF( lk_vopt_loop ) THEN          ! vector opt., forced unroll 
     535               DO jk = 1, jpkm1 
     536                  DO ji = 1, jpij 
     537                     un_b(ji,1) = un_b(ji,1) + fse3u(ji,1,jk) * un(ji,1,jk) 
     538                     vn_b(ji,1) = vn_b(ji,1) + fse3v(ji,1,jk) * vn(ji,1,jk) 
     539                  END DO 
     540               END DO 
     541            ELSE                             ! No  vector opt. 
     542               DO jk = 1, jpkm1 
     543                  un_b(:,:) = un_b(:,:) + fse3u(:,:,jk) * un(:,:,jk) 
     544                  vn_b(:,:) = vn_b(:,:) + fse3v(:,:,jk) * vn(:,:,jk) 
     545               END DO 
     546            ENDIF 
     547         ENDIF 
     548      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
     549         CALL iom_rstput( kt, nitrst, numrow, 'sshb'  , sshb  (:,:) ) 
     550         CALL iom_rstput( kt, nitrst, numrow, 'sshn'  , sshn  (:,:) ) 
     551         CALL iom_rstput( kt, nitrst, numrow, 'sshb_b', sshb_b(:,:) )   ! free surface issued 
     552         CALL iom_rstput( kt, nitrst, numrow, 'sshn_b', sshn_b(:,:) )   ! from barotropic loop 
     553         CALL iom_rstput( kt, nitrst, numrow, 'un_b'  , un_b  (:,:) )   ! horizontal transports issued 
     554         CALL iom_rstput( kt, nitrst, numrow, 'vn_b'  , vn_b  (:,:) )   ! from barotropic loop 
    508555      ENDIF 
    509        
    510    END SUBROUTINE dyn_spg_ts 
     556      ! 
     557   END SUBROUTINE ts_rst 
     558 
    511559#else 
    512560   !!---------------------------------------------------------------------- 
  • trunk/NEMO/OPA_SRC/SOL/solmat.F90

    r413 r508  
    44   !! solver       : construction of the matrix  
    55   !!====================================================================== 
     6   !! History :   1.0  !  88-04  (G. Madec)  Original code 
     7   !!                  !  93-03  (M. Guyon)  symetrical conditions 
     8   !!                  !  93-06  (M. Guyon)  suppress pointers 
     9   !!                  !  96-05  (G. Madec)  merge sor and pcg formulations 
     10   !!                  !  96-11  (A. Weaver)  correction to preconditioning 
     11   !!                  !  98-02  (M. Guyon)  FETI method 
     12   !!             8.5  !  02-08  (G. Madec)  F90: Free form 
     13   !!                  !  02-11  (C. Talandier, A-M. Treguier) Free surface & Open boundaries 
     14   !!             9.0  !  05-09  (R. Benshila)  add sol_exd for extra outer halo 
     15   !!             9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
     16   !!             9.0  !  06-07  (S. Masson)  distributed restart using iom 
     17   !!---------------------------------------------------------------------- 
    618 
    719   !!---------------------------------------------------------------------- 
     
    2941   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    3042   !! $Header$  
    31    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     43   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    3244   !!---------------------------------------------------------------------- 
    3345 
     
    5567      !!              - gcdmat : preconditioning matrix (diagonal elements) 
    5668      !!              - gcdprc : inverse of the preconditioning matrix 
    57       !! 
    58       !! History : 
    59       !!   1.0  !  88-04  (G. Madec)  Original code 
    60       !!        !  91-11  (G. Madec) 
    61       !!        !  93-03  (M. Guyon)  symetrical conditions 
    62       !!        !  93-06  (M. Guyon)  suppress pointers 
    63       !!        !  96-05  (G. Madec)  merge sor and pcg formulations 
    64       !!        !  96-11  (A. Weaver)  correction to preconditioning 
    65       !!        !  98-02  (M. Guyon)  FETI method 
    66       !!   8.5  !  02-08  (G. Madec)  F90: Free form 
    67       !!        !  02-11  (C. Talandier, A-M. Treguier) Free surface & Open boundaries 
    68       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    6969      !!---------------------------------------------------------------------- 
    7070      !! * Arguments 
  • trunk/NEMO/OPA_SRC/ZDF/zdftke.F90

    r474 r508  
    55   !!                 turbulent closure parameterization 
    66   !!===================================================================== 
     7   !! History :   6.0  !  91-03  (b. blanke)  Original code 
     8   !!             7.0  !  91-11  (G. Madec)   bug fix 
     9   !!             7.1  !  92-10  (G. Madec)   new mixing length and eav 
     10   !!             7.2  !  93-03  (M. Guyon)   symetrical conditions 
     11   !!             7.3  !  94-08  (G. Madec, M. Imbard)   npdl flag 
     12   !!             7.5  !  96-01  (G. Madec)   s-coordinates 
     13   !!             8.0  !  97-07  (G. Madec)   lbc 
     14   !!             8.1  !  99-01  (E. Stretta) new option for the mixing length 
     15   !!             8.5  !  02-06  (G. Madec) add zdf_tke_init routine 
     16   !!             8.5  !  02-08  (G. Madec)  ri_c and Free form, F90 
     17   !!             9.0  !  04-10  (C. Ethe )  1D configuration 
     18   !!             9.0  !  02-08  (G. Madec)  autotasking optimization 
     19   !!             9.0  !  06-07  (S. Masson)  distributed restart using iom 
     20   !!---------------------------------------------------------------------- 
    721#if defined key_zdftke   ||   defined key_esopa 
    822   !!---------------------------------------------------------------------- 
    9    !!   'key_zdftke'                                             TKE scheme 
     23   !!   'key_zdftke'                                   TKE vertical physics 
     24   !!---------------------------------------------------------------------- 
    1025   !!---------------------------------------------------------------------- 
    1126   !!   zdf_tke      : update momentum and tracer Kz from a tke scheme 
    1227   !!   zdf_tke_init : initialization, namelist read, and parameters control 
     28   !!   tke_rst      : read/write tke restart in ocean restart file 
    1329   !!---------------------------------------------------------------------- 
    14    !! * Modules used 
    1530   USE oce             ! ocean dynamics and active tracers  
    1631   USE dom_oce         ! ocean space and time domain 
    1732   USE zdf_oce         ! ocean vertical physics 
    18    USE in_out_manager  ! I/O manager 
    19    USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2033   USE phycst          ! physical constants 
    2134   USE taumod          ! surface stress 
     35   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2236   USE prtctl          ! Print control 
     37   USE in_out_manager  ! I/O manager 
     38   USE iom 
     39   USE restart         ! only for lrst_oce 
    2340 
    2441   IMPLICIT NONE 
    2542   PRIVATE 
    2643 
    27    !! * Routine accessibility 
    28    PUBLIC zdf_tke        ! routine called in step module 
    29    PUBLIC zdf_tke_init   ! routine called in zdftke_jki module 
    30  
    31    !! * Share Module variables 
    32    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftke = .TRUE.    !: TKE vertical mixing flag 
    33    LOGICAL, PUBLIC ::         & !!: ** tke namelist (namtke) ** 
    34      ln_rstke = .FALSE.          !: =T restart with tke from a run without tke with  
    35      !                           !  a none zero initial value for en 
    36    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &   !: 
    37       en                         !: now turbulent kinetic energy 
    38  
    39    INTEGER, PUBLIC ::         & !!! ** tke namelist (namtke) ** 
    40       nitke = 50 ,            &  ! number of restart iterative loops 
    41       nmxl  =  2 ,            &  ! = 0/1/2/3 flag for the type of mixing length used 
    42       npdl  =  1 ,            &  ! = 0/1/2 flag on prandtl number on vert. eddy coeff. 
    43       nave  =  1 ,            &  ! = 0/1 flag for horizontal average on avt, avmu, avmv 
    44       navb  =  0                 ! = 0/1 flag for constant or profile background avt 
    45    REAL(wp), PUBLIC ::        & !!! ** tke namlist (namtke) ** 
    46       ediff = 0.1_wp       ,  &  ! coeff. for vertical eddy coef.; avt=ediff*mxl*sqrt(e) 
    47       ediss = 0.7_wp       ,  &  ! coef. of the Kolmogoroff dissipation  
    48       ebb   = 3.75_wp      ,  &  ! coef. of the surface input of tke 
    49       efave = 1._wp        ,  &  ! coef. for the tke vert. diff. coeff.; avtke=efave*avm 
    50       emin  = 0.7071e-6_wp ,  &  ! minimum value of tke (m2/s2) 
    51       emin0 = 1.e-4_wp     ,  &  ! surface minimum value of tke (m2/s2) 
    52       ri_c  = 2._wp / 9._wp      ! critic Richardson number 
    53    REAL(wp), PUBLIC ::        & 
    54       eboost                     ! multiplicative coeff of the shear product. 
    55  
    56    !! caution vectopt_memory change the solution (last digit of the solver stat) 
     44   PUBLIC   zdf_tke        ! routine called in step module 
     45   PUBLIC   zdf_tke_init   ! routine also called in zdftke_jki module 
     46   PUBLIC   tke_rst        ! routine also called in zdftke_jki module 
     47 
     48   LOGICAL , PUBLIC, PARAMETER              ::   lk_zdftke = .TRUE.  !: TKE vertical mixing flag 
     49   REAL(wp), PUBLIC                         ::   eboost              !: multiplicative coeff of the shear product. 
     50   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   en                  !: now turbulent kinetic energy 
    5751# if defined key_vectopt_memory 
    58    REAL(wp), DIMENSION(jpi,jpj,jpk), PUBLIC ::   & 
    59       etmean,    &  ! coefficient used for horizontal smoothing 
    60       eumean,    &  ! at t-, u- and v-points 
    61       evmean        ! 
     52   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   etmean              !: coefficient used for horizontal smoothing 
     53   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   eumean, evmean      !: at t-, u- and v-points 
    6254# endif 
    6355 
     56   !! * Namelist (namtke) 
     57   LOGICAL , PUBLIC ::   ln_rstke = .FALSE.         !: =T restart with tke from a run without tke with  
     58     !                                              !  a none zero initial value for en 
     59   INTEGER , PUBLIC ::   nitke = 50 ,            &  !: number of restart iterative loops 
     60      &                  nmxl  =  2 ,            &  !: = 0/1/2/3 flag for the type of mixing length used 
     61      &                  npdl  =  1 ,            &  !: = 0/1/2 flag on prandtl number on vert. eddy coeff. 
     62      &                  nave  =  1 ,            &  !: = 0/1 flag for horizontal average on avt, avmu, avmv 
     63      &                  navb  =  0                 !: = 0/1 flag for constant or profile background avt 
     64   REAL(wp), PUBLIC ::   ediff = 0.1_wp       ,  &  !: coeff. for vertical eddy coef.; avt=ediff*mxl*sqrt(e) 
     65      &                  ediss = 0.7_wp       ,  &  !: coef. of the Kolmogoroff dissipation  
     66      &                  ebb   = 3.75_wp      ,  &  !: coef. of the surface input of tke 
     67      &                  efave = 1._wp        ,  &  !: coef. for the tke vert. diff. coeff.; avtke=efave*avm 
     68      &                  emin  = 0.7071e-6_wp ,  &  !: minimum value of tke (m2/s2) 
     69      &                  emin0 = 1.e-4_wp     ,  &  !: surface minimum value of tke (m2/s2) 
     70      &                  ri_c  = 2._wp / 9._wp      !: critic Richardson number 
     71   NAMELIST/namtke/ ln_rstke, ediff, ediss, ebb, efave, emin, emin0,   & 
     72      &             ri_c, nitke, nmxl, npdl, nave, navb 
     73 
    6474# if defined key_cfg_1d 
    65    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &    
    66       e_dis,    &   ! dissipation turbulent lengh scale 
    67       e_mix,    &   ! mixing turbulent lengh scale 
    68       e_pdl,    &   ! prandl number 
    69       e_ric         ! local Richardson number 
     75   !                                                                   ! 1D cfg only 
     76   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   e_dis, e_mix,      &  ! dissipation and mixing turbulent lengh scales 
     77      &                                          e_pdl, e_ric          ! prandl and local Richardson numbers 
    7078#endif 
    7179 
     
    7482#  include "vectopt_loop_substitute.h90" 
    7583   !!---------------------------------------------------------------------- 
    76    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
     84   !!   OPA 9.0 , LOCEAN-IPSL (2006)  
     85   !! $Header$  
     86   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    7787   !!---------------------------------------------------------------------- 
    7888 
    7989CONTAINS 
    8090 
    81    SUBROUTINE zdf_tke ( kt ) 
     91   SUBROUTINE zdf_tke( kt ) 
    8292      !!---------------------------------------------------------------------- 
    8393      !!                   ***  ROUTINE zdf_tke  *** 
     
    136146      !!                update avt, avmu, avmv (before vertical eddy coef.) 
    137147      !! 
    138       !! References : 
    139       !!      Gaspar et al., jgr, 95, 1990, 
    140       !!      Blanke and Delecluse, jpo, 1991 
    141       !! History : 
    142       !!   6.0  !  91-03  (b. blanke)  Original code 
    143       !!   7.0  !  91-11  (G. Madec)   bug fix 
    144       !!   7.1  !  92-10  (G. Madec)   new mixing length and eav 
    145       !!   7.2  !  93-03  (M. Guyon)   symetrical conditions 
    146       !!   7.3  !  94-08  (G. Madec, M. Imbard)   npdl flag 
    147       !!   7.5  !  96-01  (G. Madec)   s-coordinates 
    148       !!   8.0  !  97-07  (G. Madec)   lbc 
    149       !!   8.1  !  99-01  (E. Stretta) new option for the mixing length 
    150       !!   8.5  !  02-08  (G. Madec)  ri_c and Free form, F90 
    151       !!   9.0  !  04-10  (C. Ethe )  1D configuration 
     148      !! References : Gaspar et al., jgr, 95, 1990, 
     149      !!              Blanke and Delecluse, jpo, 1991 
    152150      !!---------------------------------------------------------------------- 
    153       !! * Modules used 
    154151      USE oce     , zwd   => ua,  &  ! use ua as workspace 
    155152         &          zmxlm => ta,  &  ! use ta as workspace 
    156153         &          zmxld => sa      ! use sa as workspace 
    157  
    158       !! * arguments 
    159       INTEGER, INTENT( in  ) ::   kt ! ocean time step 
    160  
    161       !! * local declarations 
    162       INTEGER ::   ji, jj, jk        ! dummy loop arguments 
    163       REAL(wp) ::   & 
    164          zmlmin, zbbrau,          &  ! temporary scalars 
    165          zfact1, zfact2, zfact3,  &  ! 
    166          zrn2, zesurf,            &  ! 
    167          ztx2, zty2, zav,         &  ! 
    168          zcoef, zcof, zsh2,       &  ! 
    169          zdku, zdkv, zpdl, zri,   &  ! 
    170          zsqen, zesh2,            &  ! 
    171          zemxl, zemlm, zemlp 
     154      ! 
     155      INTEGER, INTENT(in) ::   kt ! ocean time step 
     156      ! 
     157      INTEGER  ::   ji, jj, jk                  ! dummy loop arguments 
     158      REAL(wp) ::   zmlmin, zbbrau,          &  ! temporary scalars 
     159         &          zfact1, zfact2, zfact3,  &  ! 
     160         &          zrn2, zesurf,            &  ! 
     161         &          ztx2, zty2, zav,         &  ! 
     162         &          zcoef, zcof, zsh2,       &  ! 
     163         &          zdku, zdkv, zpdl, zri,   &  ! 
     164         &          zsqen, zesh2,            &  ! 
     165         &          zemxl, zemlm, zemlp 
    172166      !!-------------------------------------------------------------------- 
    173167 
    174       ! Initialization (first time-step only) 
    175       ! -------------- 
    176       IF( kt == nit000  )   CALL zdf_tke_init 
    177  
    178       ! Local constant initialization 
     168      IF( kt == nit000  )   CALL zdf_tke_init      ! Initialization (first time-step only) 
     169 
     170      !                                            ! Local constant initialization 
    179171      zmlmin = 1.e-8 
    180172      zbbrau =  .5 * ebb / rau0 
     
    183175      zfact3 = 0.5 * rdt * ediss 
    184176 
    185  
    186177      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    187178      ! I.  Mixing length 
    188179      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    189  
    190180 
    191181      ! Buoyancy length scale: l=sqrt(2*e/n**2) 
     
    204194         END DO 
    205195      END DO 
    206  
    207196 
    208197      ! Physical limits for the mixing length 
     
    291280# endif 
    292281 
    293  
    294282      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    295283      ! II  Tubulent kinetic energy time stepping 
    296284      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    297  
    298285 
    299286      ! 1. Vertical eddy viscosity on tke (put in zmxlm) and first estimate of avt 
     
    475462      CALL lbc_lnk( en , 'W', 1. )   ;   CALL lbc_lnk( avt, 'W', 1. ) 
    476463 
    477  
    478464      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    479465      ! III.  Before vertical eddy vicosity and diffusivity coefficients 
     
    601587      ! ------------------------------===== 
    602588      CALL lbc_lnk( avt, 'W', 1. ) 
     589 
     590      ! write en in restart file 
     591      ! ------------------------ 
     592      IF( lrst_oce )   CALL tke_rst( kt, 'WRITE' ) 
    603593 
    604594      IF(ln_ctl) THEN 
     
    624614      !! ** Action  :   Increase by 1 the nstop flag is setting problem encounter 
    625615      !! 
    626       !! history : 
    627       !!  8.5  ! 02-06 (G. Madec) original code 
    628616      !!---------------------------------------------------------------------- 
    629       !! * Module used 
    630617      USE dynzdf_exp 
    631618      USE trazdf_exp 
    632  
    633       !! * local declarations 
    634       !! caution vectopt_memory change the solution (last digit of the solver stat) 
     619      ! 
    635620# if defined key_vectopt_memory 
    636       INTEGER ::   ji, jj, jk, jit   ! dummy loop indices 
     621      ! caution vectopt_memory change the solution (last digit of the solver stat) 
     622      INTEGER ::   ji, jj, jk   ! dummy loop indices 
    637623# else 
    638       INTEGER ::           jk, jit   ! dummy loop indices 
     624      INTEGER ::           jk   ! dummy loop indices 
    639625# endif 
    640  
    641       NAMELIST/namtke/ ln_rstke, ediff, ediss, ebb, efave, emin, emin0,   & 
    642          ri_c, nitke, nmxl, npdl, nave, navb 
    643626      !!---------------------------------------------------------------------- 
    644627 
     
    681664      ! Check nmxl and npdl values 
    682665      IF( nmxl < 0 .OR. nmxl > 3 ) CALL ctl_stop( '          bad flag: nmxl is < 0 or > 3 ' ) 
    683       IF ( npdl < 0 .OR. npdl > 1 ) CALL ctl_stop( '          bad flag: npdl is < 0 or > 1 ' ) 
     666      IF( npdl < 0 .OR. npdl > 1 ) CALL ctl_stop( '          bad flag: npdl is < 0 or > 1 ' ) 
    684667 
    685668      ! Horizontal average : initialization of weighting arrays  
     
    691674         IF(lwp) WRITE(numout,*) '          no horizontal average on avt, avmu, avmv' 
    692675         IF(lwp) WRITE(numout,*) '          only in very high horizontal resolution !' 
    693 !! caution vectopt_memory change the solution (last digit of the solver stat) 
    694676# if defined key_vectopt_memory 
     677         ! caution vectopt_memory change the solution (last digit of the solver stat) 
    695678         ! weighting mean arrays etmean, eumean and evmean 
    696679         !           ( 1  1 )                                          ( 1 ) 
     
    720703      CASE ( 1 )                ! horizontal average  
    721704         IF(lwp) WRITE(numout,*) '          horizontal average on avt, avmu, avmv' 
    722 !! caution vectopt_memory change the solution (last digit of the solver stat) 
    723705# if defined key_vectopt_memory 
     706         ! caution vectopt_memory change the solution (last digit of the solver stat) 
    724707         ! weighting mean arrays etmean, eumean and evmean 
    725708         !           ( 1  1 )              ( 1/2  1/2 )             ( 1/2  1  1/2 ) 
     
    790773 
    791774 
    792       ! Initialization of turbulent kinetic energy ( en ) 
     775      ! read or initialize turbulent kinetic energy ( en ) 
    793776      ! ------------------------------------------------- 
    794       IF( ln_rstart ) THEN 
    795          ! no en field in the restart file, en set by iterative loop 
    796          IF( ln_rstke ) THEN 
    797             en (:,:,:) = emin * tmask(:,:,:) 
    798             DO jit = 2, nitke+1 
    799                CALL zdf_tke( jit ) 
    800             END DO 
    801          ENDIF 
    802          ! otherwise en is already read in dtrlec called by inidtr 
    803       ELSE 
    804          ! no restart: en set to emin 
    805          en(:,:,:) = emin * tmask(:,:,:) 
    806       ENDIF 
    807  
     777      CALL tke_rst( nit000, 'READ' ) 
     778      ! 
    808779   END SUBROUTINE zdf_tke_init 
     780 
     781 
     782   SUBROUTINE tke_rst( kt, cdrw ) 
     783     !!--------------------------------------------------------------------- 
     784     !!                   ***  ROUTINE ts_rst  *** 
     785     !!                      
     786     !! ** Purpose : Read or write filtered free surface arrays in restart file 
     787     !! 
     788     !! ** Method  :  
     789     !! 
     790     !!---------------------------------------------------------------------- 
     791     INTEGER         , INTENT(in) ::   kt         ! ocean time-step 
     792     CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
     793     ! 
     794     INTEGER ::   jit   ! dummy loop indices 
     795     !!---------------------------------------------------------------------- 
     796     ! 
     797     IF( TRIM(cdrw) == 'READ' ) THEN 
     798        IF( ln_rstart ) THEN 
     799           IF( iom_varid( numror, 'en' ) > 0 .AND. .NOT.(ln_rstke) ) THEN  
     800              CALL iom_get( numror, jpdom_local, 'en', en ) 
     801           ELSE 
     802              IF(lwp .AND. iom_varid(numror,'en') > 0 ) WRITE(numout,*) ' ===>>>> : previous run without tke scheme' 
     803              IF(lwp .AND. ln_rstke ) WRITE(numout,*) ' ===>>>> : We do not use en from the restart file' 
     804              IF(lwp) WRITE(numout,*) ' ===>>>> : en set by iterative loop' 
     805              IF(lwp) WRITE(numout,*) ' =======             =========' 
     806              en (:,:,:) = emin * tmask(:,:,:) 
     807              DO jit = 2, nitke+1 
     808                 CALL zdf_tke( jit ) 
     809              END DO 
     810           ENDIF 
     811        ELSE 
     812           en(:,:,:) = emin * tmask(:,:,:)      ! no restart: en set to emin 
     813        ENDIF 
     814     ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
     815        CALL iom_rstput( kt, nitrst, numrow, 'en', en ) 
     816     ENDIF 
     817     ! 
     818   END SUBROUTINE tke_rst 
    809819 
    810820#else 
     
    812822   !!   Dummy module :                                        NO TKE scheme 
    813823   !!---------------------------------------------------------------------- 
    814    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftke = .FALSE.   !: TKE flag 
     824   PUBLIC, LOGICAL, PARAMETER ::   lk_zdftke = .FALSE.   !: TKE flag 
    815825CONTAINS 
    816826   SUBROUTINE zdf_tke( kt )          ! Empty routine 
  • trunk/NEMO/OPA_SRC/ZDF/zdftke_jki.F90

    r463 r508  
    55   !!                 turbulent closure parameterization 
    66   !!===================================================================== 
     7   !! History : 
     8   !!   9.0  !  02-08  (G. Madec)  autotasking optimization 
     9   !!---------------------------------------------------------------------- 
    710#if defined key_zdftke   ||   defined key_esopa 
    811   !!---------------------------------------------------------------------- 
     
    2326   USE taumod          ! surface stress 
    2427   USE prtctl          ! Print control 
     28   USE restart         ! only for lrst_oce 
    2529 
    2630   IMPLICIT NONE 
     
    99103      !!      Gaspar et al., jgr, 95, 1990, 
    100104      !!      Blanke and Delecluse, jpo, 1991 
    101       !! History : 
    102       !!   9.0  !  02-08  (G. Madec)  autotasking optimization 
    103105      !!---------------------------------------------------------------------- 
    104106      !! * Modules used 
     
    518520      CALL lbc_lnk( avt, 'W', 1. ) 
    519521 
     522      ! write en in restart file 
     523      ! ------------------------ 
     524      IF( lrst_oce )   CALL tke_rst( kt, 'WRITE' ) 
     525 
    520526      IF(ln_ctl) THEN 
    521527         CALL prt_ctl(tab3d_1=en  , clinfo1=' tke  - e: ', tab3d_2=avt , clinfo2=' t: ', ovlap=1, kdim=jpk) 
  • trunk/NEMO/OPA_SRC/in_out_manager.F90

    r472 r508  
    11MODULE in_out_manager    
     2   !!====================================================================== 
     3   !!                       ***  MODULE  in_out_manager  *** 
     4   !! Ocean physics:  vertical mixing coefficient compute from the tke  
     5   !!                 turbulent closure parameterization 
     6   !!===================================================================== 
     7   !! History :   8.5  !  02-06  (G. Madec)  original code 
     8   !!             9.0  !  06-07  (S. Masson)  iom, add ctl_stop, ctl_warn 
     9   !!---------------------------------------------------------------------- 
    210 
    3    USE lib_print         ! formated print library 
     11   !!---------------------------------------------------------------------- 
     12   !!   ctl_stop   : update momentum and tracer Kz from a tke scheme 
     13   !!   ctl_warn   : initialization, namelist read, and parameters control 
     14   !!---------------------------------------------------------------------- 
    415   USE par_kind 
    516   USE par_oce 
     17   USE lib_print         ! formated print library 
    618 
    719   PUBLIC 
    820 
    921   !!---------------------------------------------------------------------- 
    10    !! namelist parameters 
    11    !! ------------------------------------- 
    12    ! namrun:  parameters of the run 
    13    CHARACTER (len=16) ::    &   !: 
    14       cexper = "exp0"           !: experiment name used for output filename 
    15     
    16    LOGICAL ::   &              !!: * namelist namrun * 
    17       ln_rstart = .FALSE. ,  &  !: start from (F) rest or (T) a restart file 
    18       ln_ctl    = .FALSE.       !: run control for debugging 
    19     
    20    INTEGER ::                & !!: * namelist namrun * 
    21       no     = 0        ,    &  !: job number 
    22       nrstdt = 0        ,    &  !: control of the time step (0, 1 or 2) 
    23       nit000 = 1        ,    &  !: index of the first time step 
    24       nitend = 10       ,    &  !: index of the last time step 
    25       ndate0 = 961115   ,    &  !: initial calendar date aammjj 
    26       nleapy = 0        ,    &  !: Leap year calendar flag (0/1 or 30) 
    27       ninist = 0        ,    &  !: initial state output flag (0/1) 
    28       nbench = 0                !: benchmark parameter (0/1) 
     22   !!                   namrun namelist parameters 
     23   !!---------------------------------------------------------------------- 
     24   CHARACTER (len=16) ::   cexper    = "exp0"        !: experiment name used for output filename 
     25   LOGICAL            ::   ln_rstart = .FALSE. ,  &  !: start from (F) rest or (T) a restart file 
     26      &                    ln_ctl    = .FALSE.       !: run control for debugging 
     27   INTEGER            ::   no     = 0        ,    &  !: job number 
     28      &                    nrstdt = 0        ,    &  !: control of the time step (0, 1 or 2) 
     29      &                    nit000 = 1        ,    &  !: index of the first time step 
     30      &                    nitend = 10       ,    &  !: index of the last time step 
     31      &                    ndate0 = 961115   ,    &  !: initial calendar date aammjj 
     32      &                    nleapy = 0        ,    &  !: Leap year calendar flag (0/1 or 30) 
     33      &                    ninist = 0        ,    &  !: initial state output flag (0/1) 
     34      &                    nbench = 0                !: benchmark parameter (0/1) 
    2935    
    3036   !!---------------------------------------------------------------------- 
    31    !! output monitoring 
    32    !! ----------------------------------- 
    33  
    34    INTEGER ::                &  !: 
    35       nstock =   10 ,        &  !: restart file frequency 
    36       nprint =    0 ,        &  !: level of print (0 no print) 
    37       nwrite =   10 ,        &  !: restart file frequency 
    38       nictls =    0 ,        &  !: Start i indice for the SUM control 
    39       nictle =    0 ,        &  !: End   i indice for the SUM control 
    40       njctls =    0 ,        &  !: Start j indice for the SUM control 
    41       njctle =    0 ,        &  !: End   j indice for the SUM control 
    42       isplt  =    1 ,        &  !: number of processors following i 
    43       jsplt  =    1 ,        &  !: number of processors following j 
    44       ijsplt =    1             !: nb of local domain = nb of processors 
     37   !!                    output monitoring 
     38   !!---------------------------------------------------------------------- 
     39   INTEGER ::   nstock =   10 ,        &  !: restart file frequency 
     40      &         nprint =    0 ,        &  !: level of print (0 no print) 
     41      &         nwrite =   10 ,        &  !: restart file frequency 
     42      &         nictls =    0 ,        &  !: Start i indice for the SUM control 
     43      &         nictle =    0 ,        &  !: End   i indice for the SUM control 
     44      &         njctls =    0 ,        &  !: Start j indice for the SUM control 
     45      &         njctle =    0 ,        &  !: End   j indice for the SUM control 
     46      &         isplt  =    1 ,        &  !: number of processors following i 
     47      &         jsplt  =    1 ,        &  !: number of processors following j 
     48      &         ijsplt =    1             !: nb of local domain = nb of processors 
    4549 
    4650   !!---------------------------------------------------------------------- 
    47    !! logical units 
    48    !! ------------------------------ 
    49    INTEGER ::                &  !: 
    50       numstp     =  1 ,      &  !: logical unit for time step 
    51       numout     =  2 ,      &  !: logical unit for output print 
    52       numnam     =  3 ,      &  !: logical unit for namelist 
    53       numnam_ice =  4 ,      &  !: logical unit for ice namelist 
    54       numevo_ice = 17 ,      &  !: logical unit for ice variables (temp. evolution) 
    55       numice_dmp = 18 ,      &  !: logical unit for ice variables (damping) 
    56       numsol     = 25 ,      &  !: logical unit for solver statistics 
    57       numwri     = 40 ,      &  !: logical unit for output write 
    58       numisp     = 41 ,      &  !: logical unit for island statistics 
    59       numgap     = 45 ,      &  !: logical unit for differences diagnostic 
    60       numbol     = 67 ,      &  !: logical unit for "bol" diagnostics 
    61       numptr     = 68 ,      &  !: logical unit for Poleward TRansports 
    62       numflo     = 69           !: logical unit for drifting floats 
    63       !                         !: * coupled units 
     51   !!                        logical units 
     52   !!---------------------------------------------------------------------- 
     53   INTEGER ::   numstp     =  1 ,      &  !: logical unit for time step 
     54      &         numout     =  2 ,      &  !: logical unit for output print 
     55      &         numnam     =  3 ,      &  !: logical unit for namelist 
     56      &         numnam_ice =  4 ,      &  !: logical unit for ice namelist 
     57      &         numevo_ice = 17 ,      &  !: logical unit for ice variables (temp. evolution) 
     58      &         numsol     = 25 ,      &  !: logical unit for solver statistics 
     59      &         numwri     = 40 ,      &  !: logical unit for output write 
     60      &         numisp     = 41 ,      &  !: logical unit for island statistics 
     61      &         numgap     = 45 ,      &  !: logical unit for differences diagnostic 
     62      &         numbol     = 67 ,      &  !: logical unit for "bol" diagnostics 
     63      &         numptr     = 68 ,      &  !: logical unit for Poleward TRansports 
     64      &         numflo     = 69           !: logical unit for drifting floats 
    6465 
    6566   !!---------------------------------------------------------------------- 
     
    6768   !!---------------------------------------------------------------------- 
    6869    
    69    INTEGER ::                &  !: 
    70       nstop = 0 ,            &  !: e r r o r  flag (=number of reason for a 
    71       !                         !                   prematurely stop the run) 
    72       nwarn = 0                 !: w a r n i n g  flag (=number of warning 
    73       !                         !                       found during the run) 
    74  
    75     
    76    CHARACTER(LEN=100) :: ctmp1, ctmp2, ctmp3    ! temporary character 
    77    CHARACTER (len=64) ::        &                                                    !: 
    78       cform_err="(/,' ===>>> : E R R O R',     /,'         ===========',/)"    ,   & !: 
    79       cform_war="(/,' ===>>> : W A R N I N G', /,'         ===============',/)"      !: 
    80    LOGICAL ::   &               !: 
    81       lwp                ,   &  !: boolean : true on the 1st processor only 
    82       lsp_area = .TRUE.         !: to make a control print over a specific area 
     70   INTEGER            ::   nstop = 0 ,           &  !: error flag (=number of reason for a premature stop run) 
     71      &                    nwarn = 0                !: warning flag (=number of warning found during the run) 
     72   CHARACTER(LEN=100) ::   ctmp1, ctmp2, ctmp3      !: temporary character 
     73   CHARACTER (len=64) ::   cform_err="(/,' ===>>> : E R R O R',     /,'         ===========',/)"    ,   &  !: 
     74      &                    cform_war="(/,' ===>>> : W A R N I N G', /,'         ===============',/)"       !: 
     75   LOGICAL            ::   lwp               ,   &  !: boolean : true on the 1st processor only 
     76      &                    lsp_area = .TRUE.        !: to make a control print over a specific area 
    8377   !!---------------------------------------------------------------------- 
    8478   !!  OPA 9.0 , LOCEAN-IPSL (2005)  
    8579   !! $Header$  
    86    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     80   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    8781   !!---------------------------------------------------------------------- 
    8882 
    89  
    9083CONTAINS 
    91  
    9284 
    9385   SUBROUTINE ctl_stop( cd1, cd2, cd3, cd4, cd5,   & 
     
    9688      !!                  ***  ROUTINE  stop_opa  *** 
    9789      !! 
    98       !! ** Purpose : ??? 
    99       !! 
     90      !! ** Purpose : ??? blah blah.... 
    10091      !!----------------------------------------------------------------------- 
    101       CHARACTER(len=*),INTENT(in),OPTIONAL ::  cd1, cd2, cd3, cd4, cd5, cd6, cd7, cd8, cd9, cd10 
     92      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5,   & 
     93         &                                       cd6, cd7, cd8, cd9, cd10 
    10294      !!----------------------------------------------------------------------- 
    103        
     95      ! 
    10496      nstop = nstop + 1  
    10597      IF(lwp) THEN 
     
    117109      ENDIF 
    118110      CALL FLUSH(numout) 
    119  
     111      ! 
    120112   END SUBROUTINE ctl_stop 
    121113 
     
    124116      &                 cd6, cd7, cd8, cd9, cd10 ) 
    125117      !!----------------------------------------------------------------------- 
    126       !!                  ***  ROUTINE  stop_opa  *** 
     118      !!                  ***  ROUTINE  stop_warn  *** 
    127119      !! 
    128       !! ** Purpose : ??? 
    129       !! 
     120      !! ** Purpose : ???  blah blah.... 
    130121      !!----------------------------------------------------------------------- 
    131       CHARACTER(len=*),INTENT(in),OPTIONAL ::  cd1, cd2, cd3, cd4, cd5, cd6, cd7, cd8, cd9, cd10 
     122      CHARACTER(len=*), INTENT(in), OPTIONAL ::  cd1, cd2, cd3, cd4, cd5,   & 
     123         &                                       cd6, cd7, cd8, cd9, cd10 
    132124      !!----------------------------------------------------------------------- 
    133        
     125      !  
    134126      nwarn = nwarn + 1  
    135127      IF(lwp) THEN 
     
    147139      ENDIF 
    148140      CALL FLUSH(numout) 
    149  
     141      ! 
    150142   END SUBROUTINE ctl_warn 
    151143 
     144   !!===================================================================== 
    152145END MODULE in_out_manager 
  • trunk/NEMO/OPA_SRC/iom.F90

    r485 r508  
    22   !!===================================================================== 
    33   !!                    ***  MODULE  iom *** 
    4    !! 
    54   !! Input/Output manager :  Library to read input files 
    6    !! 
    7    !! Ongoing work : This code is here to help discussions about I/O 
    8    !!                library in the NEMO system 
    95   !!==================================================================== 
     6   !! History :  9.0  ! 05 12  (J. Belier) Original code 
     7   !!            9.0  ! 06 02  (S. Masson) Adaptation to NEMO 
     8   !!-------------------------------------------------------------------- 
     9   !!gm  caution add !DIR nec: improved performance to be checked as well as no result changes 
     10 
    1011   !!-------------------------------------------------------------------- 
    1112   !!   iom_open       : open a file read only 
    1213   !!   iom_close      : close a file or all files opened by iom 
    13    !!   iom_get        : read a field : interface to several routines 
     14   !!   iom_get        : read a field (interfaced to several routines) 
     15   !!   iom_gettime    : read the time axis cdvar in the file               !!gm : never call ?????? 
    1416   !!   iom_varid      : get the id of a variable in a file 
    15    !!   iom_get_gblatt : ??? 
     17   !!   iom_rstput     : write a field in a restart file (interfaced to several routines) 
    1618   !!-------------------------------------------------------------------- 
    17    !! History :  9.0  ! 05 12  (J. Belier) Original code 
    18    !!            9.0  ! 06 02  (S. Masson) Adaptation to NEMO 
    19    !!-------------------------------------------------------------------- 
    20    !! * Modules used 
    2119   USE in_out_manager  ! I/O manager 
    2220   USE dom_oce         ! ocean space and time domain 
    23    USE lbclnk          ! ??? 
    24    USE ioipsl          ! ??? 
     21   USE lbclnk          ! lateal boundary condition / mpp exchanges 
     22   USE ioipsl          ! IOIPSL library 
    2523 
    2624   IMPLICIT NONE 
    2725   PRIVATE 
    2826 
    29    PUBLIC iom_open, iom_close, iom_get, iom_varid, iom_get_gblatt 
    30  
    31    !! * Interfaces 
     27   PUBLIC iom_open, iom_close, iom_get, iom_varid, iom_rstput, iom_gettime 
     28 
    3229   INTERFACE iom_get 
    33       MODULE PROCEDURE iom_get_r_1d, iom_get_r_2d, iom_get_r_3d 
     30      MODULE PROCEDURE iom_get_r_0d, iom_get_r_1d, iom_get_r_2d, iom_get_r_3d 
    3431   END INTERFACE 
    35  
    36    !! * Share module variables 
    37    INTEGER, PARAMETER, PUBLIC ::        &  !: 
    38       jpdom_data                  = 1,  &  !: ( 1  :jpidta, 1  :jpjdta) 
    39       jpdom_global                = 2,  &  !: ( 1  :jpiglo, 1  :jpjglo) 
    40       jpdom_local                 = 3,  &  !: One of the 3 following cases 
    41       jpdom_local_full            = 4,  &  !: ( 1  :jpi   , 1  :jpi   ) 
    42       jpdom_local_noextra         = 5,  &  !: ( 1  :nlci  , 1  :nlcj  ) 
    43       jpdom_local_noovlap         = 6,  &  !: (nldi:nlei  ,nldj:nlej  ) 
    44       jpdom_unknown               = 7      !: No dimension checking 
    45  
    46    !! * Module variables 
    47    INTEGER, PARAMETER ::    & 
    48       jpmax_vars    = 50,   &  ! maximum number of variables in one file 
    49       jpmax_dims    =  5,   &  ! maximum number of dimensions for one variable 
    50       jpmax_digits  =  5       ! maximum number of digits in the file name to reference the cpu number 
    51   
     32   INTERFACE iom_rstput 
     33      MODULE PROCEDURE iom_rstput_0d, iom_rstput_1d, iom_rstput_2d, iom_rstput_3d 
     34   END INTERFACE 
     35 
     36   INTEGER, PARAMETER, PUBLIC ::   jpdom_data          = 1   !: ( 1  :jpidta, 1  :jpjdta) 
     37   INTEGER, PARAMETER, PUBLIC ::   jpdom_global        = 2   !: ( 1  :jpiglo, 1  :jpjglo) 
     38   INTEGER, PARAMETER, PUBLIC ::   jpdom_local         = 3   !: One of the 3 following cases 
     39   INTEGER, PARAMETER, PUBLIC ::   jpdom_local_full    = 4   !: ( 1  :jpi   , 1  :jpi   ) 
     40   INTEGER, PARAMETER, PUBLIC ::   jpdom_local_noextra = 5   !: ( 1  :nlci  , 1  :nlcj  ) 
     41   INTEGER, PARAMETER, PUBLIC ::   jpdom_local_noovlap = 6   !: (nldi:nlei  ,nldj:nlej  ) 
     42   INTEGER, PARAMETER, PUBLIC ::   jpdom_unknown       = 7   !: No dimension checking 
     43 
     44   INTEGER, PARAMETER ::   jpmax_vars   = 60,   &  ! maximum number of variables in one file 
     45      &                    jpmax_dims   =  4,   &  ! maximum number of dimensions for one variable 
     46      &                    jpmax_digits =  5       ! maximum number of digits for the cpu number in the file name 
    5247!$AGRIF_DO_NOT_TREAT 
    53    INTEGER :: iom_init = 0 
    54  
    55    TYPE :: flio_file 
     48   INTEGER ::   iom_init = 0 
     49   TYPE    ::   flio_file 
    5650      CHARACTER(LEN=240)                        ::   name     ! name of the file 
    57       INTEGER                                   ::   iopen    ! 1/0 is the file is open/not open 
     51      INTEGER                                   ::   iopen    ! 1(0) if the file is open(not open) 
    5852      INTEGER                                   ::   nvars    ! number of identified varibles in the file 
    5953      INTEGER                                   ::   iduld    ! id of the unlimited dimension 
    6054      CHARACTER(LEN=16), DIMENSION(jpmax_vars)  ::   cn_var   ! names of the variables 
    6155      INTEGER, DIMENSION(jpmax_vars)            ::   ndims    ! number of dimensions of the variables 
    62       LOGICAL, DIMENSION(jpmax_vars)            ::   luld     ! variable including unlimited dimension 
    63       INTEGER, DIMENSION(jpmax_dims,jpmax_vars) ::   dimsz    ! size of the dimensions of the variables 
     56      LOGICAL, DIMENSION(jpmax_vars)            ::   luld     ! variable using the unlimited dimension 
     57      INTEGER, DIMENSION(jpmax_dims,jpmax_vars) ::   dimsz    ! size of variables dimensions  
    6458      REAL(kind=wp), DIMENSION(jpmax_vars)      ::   scf      ! scale_factor of the variables 
    6559      REAL(kind=wp), DIMENSION(jpmax_vars)      ::   ofs      ! add_offset of the variables 
    6660   END TYPE flio_file 
    67    TYPE(flio_file), DIMENSION(flio_max_files)   :: iom_file ! array containing the info for all opened files 
     61   TYPE(flio_file), DIMENSION(flio_max_files)   ::   iom_file ! array containing the info for all opened files 
    6862!$AGRIF_END_DO_NOT_TREAT 
    6963 
     
    7670CONTAINS 
    7771 
    78    SUBROUTINE iom_open( cdname, knumfl, ldimg ) 
     72   SUBROUTINE iom_open( cdname, knumfl, ldwrt, kdom, ldimg ) 
    7973      !!--------------------------------------------------------------------- 
    8074      !!                   ***  SUBROUTINE  iom_open  *** 
    8175      !! 
    8276      !! ** Purpose :  open an input file read only (return 0 if not found) 
    83       !! 
    84       !! ** Method : 
    85       !! 
    8677      !!--------------------------------------------------------------------- 
    87       CHARACTER(len=*), INTENT(in )  ::   cdname   ! File name 
    88       INTEGER, INTENT(out)           ::   knumfl   ! Identifier of the opened file 
    89       LOGICAL, INTENT(in ), OPTIONAL ::   ldimg    ! flg to specify that we use dimg format 
    90  
    91       CHARACTER(LEN=100) :: clname   ! the name of the file based on cdname [[+clcpu]+clcpu] 
    92       CHARACTER(LEN=10) :: clsuffix  ! ".nc" or ".dimg" 
    93       CHARACTER(LEN=10) :: clcpu     ! the cpu number (max jpmax_digits digits) 
    94       LOGICAL :: llok                ! check the existence  
    95       INTEGER :: icnt                ! counter for digits in clcpu (max = jpmax_digits) 
     78      CHARACTER(len=*), INTENT(in   )           ::   cdname   ! File name 
     79      INTEGER         , INTENT(  out)           ::   knumfl   ! Identifier of the opened file 
     80      LOGICAL         , INTENT(in   ), OPTIONAL ::   ldwrt    ! read or write the file? 
     81      INTEGER         , INTENT(in   ), OPTIONAL ::   kdom     ! Type of domain to be written 
     82      LOGICAL         , INTENT(in   ), OPTIONAL ::   ldimg    ! use dimg format? 
     83 
     84      CHARACTER(LEN=100)    ::   clname    ! the name of the file based on cdname [[+clcpu]+clcpu] 
     85      CHARACTER(LEN=100)    ::   cltmpn    ! tempory name to store clname (in writting mode) 
     86      CHARACTER(LEN=10)     ::   clsuffix  ! ".nc" or ".dimg" 
     87      CHARACTER(LEN=10)     ::   clcpu     ! the cpu number (max jpmax_digits digits) 
     88      LOGICAL               ::   llok      ! check the existence  
     89      LOGICAL               ::   llwrt     !  
     90      INTEGER               ::   icnt      ! counter for digits in clcpu (max = jpmax_digits) 
     91      INTEGER               ::   iln, ils  ! lengths of character 
     92      INTEGER               ::   idom      ! type of domain 
     93      INTEGER               ::   ifliodom  ! model domain identifier (see flio_dom_set) 
     94      INTEGER, DIMENSION(2) ::   iszl      ! local number of points for x,y dimensions 
     95      INTEGER, DIMENSION(2) ::   ifst      ! position of first local point for x,y dimensions 
     96      INTEGER, DIMENSION(2) ::   ilst      ! position of last local point for x,y dimensions 
     97      INTEGER, DIMENSION(2) ::   ihst      ! start halo size for x,y dimensions 
     98      INTEGER, DIMENSION(2) ::   ihnd      ! end halo size for x,y dimensions 
    9699      !--------------------------------------------------------------------- 
    97  
    98       ! find the file 
     100      ! if iom_open is called for the first time: initialize iom_file(:)%iopen to 0 
     101      ! (could be done when defining iom_file in f95 but not in f90) 
     102      IF( iom_init == 0 ) THEN 
     103         iom_file(:)%iopen = 0 
     104         iom_init = 1 
     105      ENDIF 
     106      ! do we read or write the file? 
     107      ! ============= 
     108      IF( PRESENT(ldwrt) ) THEN   ;   llwrt = ldwrt 
     109      ELSE                        ;   llwrt = .FALSE. 
     110      ENDIF 
     111      ! create the file name by added, if needed, TRIM(Agrif_CFixed()) and TRIM(clsuffix) 
    99112      ! ============= 
    100113      clname   = trim(cdname) 
    101114#if defined key_agrif 
    102115      if ( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname) 
    103 #endif                 
     116#endif     
     117      ! which suffix should we use? 
    104118      clsuffix = ".nc" 
    105       IF( PRESENT(ldimg) ) THEN 
    106          IF ( ldimg ) clsuffix = ".dimg" 
    107       ENDIF 
    108       ! 
     119      IF( PRESENT(ldimg) ) THEN   ;   IF( ldimg )   clsuffix = ".dimg"   ;   ENDIF 
     120      ! Add the suffix if needed 
     121      iln = LEN_TRIM(clname) 
     122      ils = LEN_TRIM(clsuffix) 
     123      IF( iln <= ils) clname = clname(1:iln)//TRIM(clsuffix) 
     124      IF( clname(iln-ils+1:iln) /= TRIM(clsuffix) )   clname = clname(1:iln)//TRIM(clsuffix) 
     125      cltmpn = clname   ! store this name 
     126      ! try to find if the file to be opened already exist 
    109127      INQUIRE( FILE = clname, EXIST = llok ) 
    110       IF( .NOT.llok ) THEN                         ! try to complete the name with the suffix only 
    111          clname = TRIM(cdname)//TRIM(clsuffix) 
    112 #if defined key_agrif 
    113          if ( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname) 
    114 #endif                 
    115          INQUIRE( FILE = clname, EXIST = llok ) 
    116          IF( .NOT.llok ) THEN                      ! try to complete the name with both cpu number and suffix 
    117             WRITE(clcpu,*) narea-1 
    118             clcpu  = trim(adjustl(clcpu)) 
    119             clname = trim(cdname)//"_" 
    120 #if defined key_agrif 
    121             if ( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname) 
    122 #endif                 
    123             icnt = 0 
    124             INQUIRE( FILE = trim(clname)//trim(clcpu)//trim(clsuffix), EXIST = llok )  
    125             DO WHILE( .NOT.llok .AND. icnt < jpmax_digits )   ! we try fifferent formats for the cpu number by adding 0 
    126                clname = trim(clname)//"0" 
    127                INQUIRE( FILE = trim(clname)//trim(clcpu)//trim(clsuffix), EXIST = llok ) 
    128                icnt = icnt + 1 
    129             END DO 
    130             IF( .NOT.llok ) THEN                   ! no way to find the files... 
    131                CALL ctl_stop( 'iom_open: file '//trim(clname)//'... not found' ) 
     128      IF( .NOT.llok ) THEN 
     129         ! we try to add the cpu number to the name 
     130         WRITE(clcpu,*) narea-1 
     131         clcpu  = TRIM(ADJUSTL(clcpu)) 
     132         iln = INDEX(clname,TRIM(clsuffix)) 
     133         clname = clname(1:iln-1)//'_'//TRIM(clcpu)//TRIM(clsuffix) 
     134         icnt = 0 
     135         INQUIRE( FILE = clname, EXIST = llok )  
     136         ! we try different formats for the cpu number by adding 0 
     137         DO WHILE( .NOT.llok .AND. icnt < jpmax_digits ) 
     138            clcpu  = "0"//trim(clcpu) 
     139            clname = clname(1:iln-1)//'_'//TRIM(clcpu)//TRIM(clsuffix) 
     140            INQUIRE( FILE = clname, EXIST = llok ) 
     141            icnt = icnt + 1 
     142         END DO 
     143      ENDIF 
     144      ! 
     145      IF( llok ) THEN      ! Open the file 
     146         !                 ! ============= 
     147         IF( llwrt ) THEN  
     148            IF(lwp) WRITE(numout,*) '          iom_open ~~~  open existing file: '//TRIM(clname)//' in WRITE mode' 
     149            CALL flioopfd( TRIM(clname), knumfl, "WRITE" ) 
     150         ELSE 
     151            IF(lwp) WRITE(numout,*) '          iom_open ~~~  open existing file: '//TRIM(clname)//' in READ mode' 
     152            CALL flioopfd( TRIM(clname), knumfl ) 
     153         ENDIF 
     154      ELSE                 ! no way to find the file... 
     155         !                 ! ======================= 
     156         IF( llwrt ) THEN  
     157            ! file opened in write mode 
     158            ! the file does not exist, we must create it... 
     159            ! ============= 
     160            llok = .TRUE. 
     161            ! on which domain must the file be written?? 
     162            ! check the domain definition 
     163            idom = jpdom_local_noovlap   ! default definition 
     164            IF( PRESENT(kdom) )   idom = kdom 
     165            ! create the domain informations 
     166            ! ============= 
     167            SELECT CASE (idom) 
     168            CASE (jpdom_local_full) 
     169               iszl = (/ jpi             , jpj              /) 
     170               ifst = (/ nimpp           , njmpp            /) 
     171               ilst = (/ nimpp + jpi - 1 , njmpp + jpj - 1  /) 
     172               ihst = (/ nldi - 1        , nldj - 1         /) 
     173               ihnd = (/ jpi - nlei      , jpj - nlej       /) 
     174            CASE (jpdom_local_noextra) 
     175               iszl = (/ nlci            , nlcj             /) 
     176               ifst = (/ nimpp           , njmpp            /) 
     177               ilst = (/ nimpp + nlci - 1, njmpp + nlcj - 1 /) 
     178               ihst = (/ nldi - 1        , nldj - 1         /) 
     179               ihnd = (/ nlci - nlei     , nlcj - nlej      /) 
     180            CASE (jpdom_local_noovlap) 
     181               iszl = (/ nlei - nldi + 1 , nlej - nldj + 1  /) 
     182               ifst = (/ nimpp + nldi - 1, njmpp + nldj - 1 /) 
     183               ilst = (/ nimpp + nlei - 1, njmpp + nlej - 1 /) 
     184               ihst = (/ 0               , 0                /) 
     185               ihnd = (/ 0               , 0                /) 
     186            CASE DEFAULT 
     187               llok = .FALSE. 
     188               CALL ctl_stop( 'iom_open: wrong value of kdom, only jpdom_local* cases are accepted' ) 
     189            END SELECT 
     190            IF( llok ) THEN 
     191               CALL flio_dom_set( jpnij, narea-1, (/1, 2/), (/jpiglo, jpjglo/)   & 
     192                    &                  , iszl, ifst, ilst, ihst, ihnd, 'BOX', ifliodom )         
     193               ! create the file 
     194               ! ============= 
     195               ! Note that fliocrfd may change the value of clname (add the cpu number...) 
     196               clname = cltmpn   ! get back the file name without the cpu number in it 
     197               IF(lwp) WRITE(numout,*) '          iom_open ~~~  create new file: '//trim(clname)//' in WRITE mode' 
     198               CALL fliocrfd( clname, (/'x'    , 'y'    , 'z', 't'/)   & 
     199                    &               , (/iszl(1), iszl(2), jpk, -1 /)   & 
     200                    &               , knumfl, ifliodom ) 
    132201            ENDIF 
    133             clname = trim(clname)//trim(clcpu)//trim(clsuffix) 
     202         ELSE 
     203            ! the file is open for read-only, it must exist... 
     204            iln = INDEX( cltmpn,TRIM(clsuffix) ) 
     205            CALL ctl_stop( 'iom_open: file '//cltmpn(1:iln-1)//'* not found' ) 
    134206         ENDIF 
    135207      ENDIF 
    136  
    137       ! Open the file 
     208      ! start to fill the information of opened files 
    138209      ! ============= 
    139210      IF( llok ) THEN                          
    140          IF (lwp) WRITE(numout,*) 'iom_open ~~~  open file: '//trim(clname) 
    141          CALL flioopfd( trim(clname), knumfl ) 
    142          IF( iom_init == 0 ) THEN 
    143             iom_file(:)%iopen = 0 
    144             iom_init = 1 
    145          ENDIF 
    146211         iom_file(knumfl)%iopen      = 1 
    147212         iom_file(knumfl)%name       = TRIM(clname) 
     
    152217         ! does the file contain time axis (that must be unlimitted) ? 
    153218         CALL flioinqf( knumfl, id_uld = iom_file(knumfl)%iduld ) 
     219         IF(lwp) WRITE(numout,*) '                   ---> OK' 
    154220      ELSE 
    155          knumfl = 0 
    156       ENDIF 
    157    
     221         knumfl = 0      ! return error flag 
     222      ENDIF 
     223      ! 
    158224   END SUBROUTINE iom_open 
    159225 
     
    164230      !! 
    165231      !! ** Purpose : close an input file, or all files opened by iom 
    166       !! 
    167       !! ** Method : 
    168       !! 
    169232      !!-------------------------------------------------------------------- 
    170       INTEGER,INTENT(in), OPTIONAL ::   knumfl   ! Identifier of the file to be closed 
    171       !                                          ! If this argument is not present, 
    172       !                                          ! all the files opened by iom are closed. 
    173  
    174       INTEGER           ::   jf            ! dummy loop indices 
    175       INTEGER           ::   i_s, i_e      ! temporary integer 
     233      INTEGER, INTENT(in), OPTIONAL ::   knumfl   ! Identifier of the file to be closed 
     234      !                                           ! No argument : all the files opened by iom are closed 
     235 
     236      INTEGER ::   jf         ! dummy loop indices 
     237      INTEGER ::   i_s, i_e   ! temporary integer 
    176238      !--------------------------------------------------------------------- 
    177  
     239      ! 
    178240      IF( PRESENT(knumfl) ) THEN 
    179241         i_s = knumfl 
     
    183245         i_e = flio_max_files 
    184246      ENDIF 
    185       IF ( i_s > 0 ) THEN 
     247       
     248      IF( i_s > 0 ) THEN 
    186249         DO jf = i_s, i_e 
    187250            IF( iom_file(jf)%iopen > 0 ) THEN 
    188251               CALL flioclo( jf ) 
     252               IF(lwp) WRITE(numout,*) '          iom_close, close file: '//TRIM(iom_file(knumfl)%name)//' ok' 
    189253               iom_file(jf)%iopen      = 0 
    190254               iom_file(jf)%name       = 'NONE' 
     
    200264         END DO 
    201265      ENDIF 
    202            
     266      !     
    203267   END SUBROUTINE iom_close 
    204268  
    205269 
    206270   !!---------------------------------------------------------------------- 
    207    !!                   INTERFACE iom_u_getv 
     271   !!                   INTERFACE iom_get_123d 
    208272   !!---------------------------------------------------------------------- 
    209    SUBROUTINE iom_get_r_1d( knumfl, kdom , cdvar , pvar  ,   & 
    210       &                             ktime, kstart, kcount ) 
    211       INTEGER               , INTENT(in )           ::   knumfl    ! Identifier of the file 
    212       INTEGER               , INTENT(in )           ::   kdom      ! Type of domain to be read 
    213       CHARACTER(len=*)      , INTENT(in )           ::   cdvar     ! Name of the variable 
    214       REAL(wp), DIMENSION(:), INTENT(out)           ::   pvar      ! read field 
    215       INTEGER               , INTENT(in ) ,OPTIONAL ::   ktime     ! record number 
    216       INTEGER , DIMENSION(:), INTENT(in ), OPTIONAL ::   kstart    ! start position of the reading in each axis  
    217       INTEGER , DIMENSION(:), INTENT(in ), OPTIONAL ::   kcount    ! number of points to be read in each axis 
    218 ! 
    219       CHARACTER(LEN=100) :: clinfo                    ! info character 
    220 ! 
    221       clinfo = 'iom_get_r_1d, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    222       IF( PRESENT(kstart) ) THEN 
    223          IF ( SIZE(kstart) /= 1 ) CALL ctl_stop( trim(clinfo), 'kstart must be a 1 element vector' ) 
    224       ENDIF 
    225       IF( PRESENT(kcount) ) THEN 
    226          IF ( SIZE(kcount) /= 1 ) CALL ctl_stop( trim(clinfo), 'kcount must be a 1 element vector' ) 
    227       ENDIF 
    228       IF ( knumfl > 0 ) CALL iom_u_getv( knumfl, kdom       , cdvar        , pv_r1d=pvar,   & 
    229            &                                     ktime=ktime, kstart=kstart, kcount=kcount ) 
     273   SUBROUTINE iom_get_r_0d( knumfl, cdvar, pvar ) 
     274      INTEGER         , INTENT(in   )                 ::   knumfl    ! Identifier of the file 
     275      CHARACTER(len=*), INTENT(in   )                 ::   cdvar     ! Name of the variable 
     276      REAL(wp)        , INTENT(  out)                 ::   pvar      ! read field 
     277      ! 
     278      IF( knumfl > 0 .AND. iom_varid( knumfl, cdvar ) > 0 )   CALL fliogetv( knumfl, cdvar, pvar ) 
     279   END SUBROUTINE iom_get_r_0d 
     280 
     281   SUBROUTINE iom_get_r_1d( knumfl, kdom, cdvar, pvar, ktime, kstart, kcount ) 
     282      INTEGER         , INTENT(in   )                         ::   knumfl    ! Identifier of the file 
     283      INTEGER         , INTENT(in   )                         ::   kdom      ! Type of domain to be read 
     284      CHARACTER(len=*), INTENT(in   )                         ::   cdvar     ! Name of the variable 
     285      REAL(wp)        , INTENT(  out), DIMENSION(:)           ::   pvar      ! read field 
     286      INTEGER         , INTENT(in   )              , OPTIONAL ::   ktime     ! record number 
     287      INTEGER         , INTENT(in   ), DIMENSION(1), OPTIONAL ::   kstart    ! start axis position of the reading  
     288      INTEGER         , INTENT(in   ), DIMENSION(1), OPTIONAL ::   kcount    ! number of points in each axis 
     289      ! 
     290      IF( knumfl > 0 ) CALL iom_get_123d( knumfl, kdom       , cdvar        , pv_r1d=pvar,   & 
     291         &                                        ktime=ktime, kstart=kstart, kcount=kcount ) 
    230292   END SUBROUTINE iom_get_r_1d 
    231    SUBROUTINE iom_get_r_2d( knumfl, kdom , cdvar , pvar  ,   & 
    232       &                             ktime, kstart, kcount ) 
    233       INTEGER,INTENT(in) :: knumfl 
    234       INTEGER,INTENT(in) :: kdom 
    235       CHARACTER(len=*),INTENT(in) :: cdvar 
    236       REAL(wp),INTENT(out),DIMENSION(:,:) :: pvar 
    237       INTEGER,INTENT(in),OPTIONAL :: ktime 
    238       INTEGER,DIMENSION(:),INTENT(in),OPTIONAL :: kstart 
    239       INTEGER,DIMENSION(:),INTENT(in),OPTIONAL :: kcount 
    240 ! 
    241       CHARACTER(LEN=100) :: clinfo                    ! info character 
    242 ! 
    243       clinfo = 'iom_get_r_2d, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    244       IF( PRESENT(kstart) ) THEN 
    245          IF ( size(kstart) /= 2 ) CALL ctl_stop(trim(clinfo), 'kstart must be a 2 element vector') 
    246       ENDIF 
    247       IF( PRESENT(kcount) ) THEN 
    248          IF ( size(kcount) /= 2 ) CALL ctl_stop(trim(clinfo), 'kcount must be a 2 element vector') 
    249       ENDIF 
    250       IF ( knumfl > 0 ) CALL iom_u_getv( knumfl, kdom       , cdvar        , pv_r2d=pvar,   & 
    251          &                                       ktime=ktime, kstart=kstart, kcount=kcount ) 
     293 
     294   SUBROUTINE iom_get_r_2d( knumfl, kdom, cdvar, pvar, ktime, kstart, kcount ) 
     295      INTEGER         , INTENT(in   )                           ::   knumfl    ! Identifier of the file 
     296      INTEGER         , INTENT(in   )                           ::   kdom      ! Type of domain to be read 
     297      CHARACTER(len=*), INTENT(in   )                           ::   cdvar     ! Name of the variable 
     298      REAL(wp)        , INTENT(  out), DIMENSION(:,:)           ::   pvar      ! read field 
     299      INTEGER         , INTENT(in   )                , OPTIONAL ::   ktime     ! record number 
     300      INTEGER         , INTENT(in   ), DIMENSION(2)  , OPTIONAL ::   kstart    ! start axis position of the reading  
     301      INTEGER         , INTENT(in   ), DIMENSION(2)  , OPTIONAL ::   kcount    ! number of points in each axis 
     302      ! 
     303      IF( knumfl > 0 ) CALL iom_get_123d( knumfl, kdom       , cdvar        , pv_r2d=pvar,   & 
     304         &                                        ktime=ktime, kstart=kstart, kcount=kcount ) 
    252305   END SUBROUTINE iom_get_r_2d 
    253    SUBROUTINE iom_get_r_3d( knumfl, kdom , cdvar , pvar  ,   & 
    254       &                             ktime, kstart, kcount ) 
    255       INTEGER,INTENT(in) :: knumfl 
    256       INTEGER,INTENT(in) :: kdom 
    257       CHARACTER(len=*),INTENT(in) :: cdvar 
    258       REAL(wp),INTENT(out),DIMENSION(:,:,:) :: pvar 
    259       INTEGER,INTENT(in),OPTIONAL :: ktime 
    260       INTEGER,DIMENSION(:),INTENT(in),OPTIONAL :: kstart 
    261       INTEGER,DIMENSION(:),INTENT(in),OPTIONAL :: kcount 
    262 ! 
    263       CHARACTER(LEN=100) :: clinfo                    ! info character 
    264 ! 
    265       clinfo = 'iom_get_r_3d, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    266       IF ( PRESENT(kstart) ) THEN 
    267          IF ( size(kstart) /= 3 ) CALL ctl_stop(trim(clinfo), 'kstart must be a 3 element vector') 
    268       ENDIF 
    269       IF ( PRESENT(kcount) ) THEN 
    270          IF ( size(kcount) /= 3 ) CALL ctl_stop(trim(clinfo), 'kcount must be a 3 element vector') 
    271       ENDIF 
    272       IF ( knumfl > 0 ) CALL iom_u_getv( knumfl, kdom       , cdvar        , pv_r3d=pvar,   & 
    273         &                                   ktime=ktime, kstart=kstart, kcount=kcount ) 
     306 
     307   SUBROUTINE iom_get_r_3d( knumfl, kdom, cdvar, pvar, ktime, kstart, kcount ) 
     308      INTEGER         , INTENT(in   )                             ::   knumfl    ! Identifier of the file 
     309      INTEGER         , INTENT(in   )                             ::   kdom      ! Type of domain to be read 
     310      CHARACTER(len=*), INTENT(in   )                             ::   cdvar     ! Name of the variable 
     311      REAL(wp)        , INTENT(  out), DIMENSION(:,:,:)           ::   pvar      ! read field 
     312      INTEGER         , INTENT(in   )                  , OPTIONAL ::   ktime     ! record number 
     313      INTEGER         , INTENT(in   ), DIMENSION(3)    , OPTIONAL ::   kstart    ! start axis position of the reading  
     314      INTEGER         , INTENT(in   ), DIMENSION(3)    , OPTIONAL ::   kcount    ! number of points in each axis 
     315      ! 
     316      IF( knumfl > 0 ) CALL iom_get_123d( knumfl, kdom       , cdvar        , pv_r3d=pvar,   & 
     317         &                                        ktime=ktime, kstart=kstart, kcount=kcount ) 
    274318   END SUBROUTINE iom_get_r_3d 
    275319   !!---------------------------------------------------------------------- 
    276320 
    277  
    278    SUBROUTINE iom_u_getv( knumfl, kdom  , cdvar , & 
    279         &                   pv_r1d, pv_r2d, pv_r3d, & 
    280         &                   ktime , kstart, kcount ) 
     321   SUBROUTINE iom_get_123d( knumfl, kdom  , cdvar ,   & 
     322        &                   pv_r1d, pv_r2d, pv_r3d,   & 
     323        &                   ktime , kstart, kcount  ) 
    281324     !!----------------------------------------------------------------------- 
    282      !!                  ***  ROUTINE  iom_u_getv  *** 
     325     !!                  ***  ROUTINE  iom_get_123d  *** 
    283326     !! 
    284327     !! ** Purpose : read a 1D/2D/3D variable 
    285328     !! 
    286      !! ** Method : read ONE time step at each CALL 
    287      !! 
     329     !! ** Method : read ONE record at each CALL 
    288330     !!----------------------------------------------------------------------- 
    289      INTEGER,                      INTENT(in )           ::   knumfl     ! Identifier of the file 
    290      INTEGER,                      INTENT(in )           ::   kdom       ! Type of domain to be read 
    291      CHARACTER(len=*),             INTENT(in )           ::   cdvar      ! Name of the variable 
    292      REAL(wp), DIMENSION(:)      , INTENT(out), OPTIONAL ::   pv_r1d     ! read field (1D case) 
    293      REAL(wp), DIMENSION(:,:)    , INTENT(out), OPTIONAL ::   pv_r2d     ! read field (2D case) 
    294      REAL(wp), DIMENSION(:,:,:)  , INTENT(out), OPTIONAL ::   pv_r3d     ! read field (3D case) 
    295      INTEGER                     , INTENT(in ), OPTIONAL ::   ktime      ! record number 
    296      INTEGER , DIMENSION(:)      , INTENT(in ), OPTIONAL ::   kstart     ! start position of the reading in each axis  
    297      INTEGER , DIMENSION(:)      , INTENT(in ), OPTIONAL ::   kcount     ! number of points to be read in each axis 
    298  
    299      INTEGER                        :: jl            ! loop on number of dimension  
    300      INTEGER                        :: idom,    &    ! type of domain 
    301           &                            idvar,   &    ! id of the variable 
    302           &                            inbdim,  &    ! number of dimensions of the variable 
    303           &                            idmspc,  &    ! number of spatial dimensions  
    304           &                            itime,   &    ! record number 
    305           &                            istop         ! temporary value of nstop 
    306      INTEGER, DIMENSION(jpmax_dims) :: istart,  &    ! starting point to read for each axis 
    307           &                            icnt,    &    ! number of value to read along each axis  
    308           &                            idimsz        ! size of the dimensions of the variable 
    309      REAL(wp)                       :: zscf, zofs    ! sacle_factor and add_offset 
    310      INTEGER                        ::  itmp         ! temporary integer 
    311      CHARACTER(LEN=100)             :: clinfo        ! info character 
     331     INTEGER                    , INTENT(in  )           ::   knumfl     ! Identifier of the file 
     332     INTEGER                    , INTENT(in  )           ::   kdom       ! Type of domain to be read 
     333     CHARACTER(len=*)           , INTENT(in  )           ::   cdvar      ! Name of the variable 
     334     REAL(wp), DIMENSION(:)     , INTENT(  out), OPTIONAL ::   pv_r1d     ! read field (1D case) 
     335     REAL(wp), DIMENSION(:,:)   , INTENT(  out), OPTIONAL ::   pv_r2d     ! read field (2D case) 
     336     REAL(wp), DIMENSION(:,:,:) , INTENT(  out), OPTIONAL ::   pv_r3d     ! read field (3D case) 
     337     INTEGER                    , INTENT(in  ), OPTIONAL ::   ktime      ! record number 
     338     INTEGER , DIMENSION(:)     , INTENT(in  ), OPTIONAL ::   kstart     ! start position of the reading in each axis  
     339     INTEGER , DIMENSION(:)     , INTENT(in  ), OPTIONAL ::   kcount     ! number of points to be read in each axis 
     340     ! 
     341     INTEGER                        ::   jl          ! loop on number of dimension  
     342     INTEGER                        ::   idom,    &  ! type of domain 
     343          &                              idvar,   &  ! id of the variable 
     344          &                              inbdim,  &  ! number of dimensions of the variable 
     345          &                              idmspc,  &  ! number of spatial dimensions  
     346          &                              itime,   &  ! record number 
     347          &                              istop       ! temporary value of nstop 
     348     INTEGER, DIMENSION(jpmax_dims) ::   istart,  &  ! starting point to read for each axis 
     349          &                              icnt,    &  ! number of value to read along each axis  
     350          &                              idimsz      ! size of the dimensions of the variable 
     351     REAL(wp)                       ::   zscf, zofs  ! sacle_factor and add_offset 
     352     INTEGER                        ::   itmp        ! temporary integer 
     353     CHARACTER(LEN=100)             ::   clinfo      ! info character 
    312354     !--------------------------------------------------------------------- 
    313      clinfo = 'iom_u_getv, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
     355     ! 
     356     clinfo = '          iom_get_123d, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    314357     ! local definition of the domain ? 
    315358     idom = kdom 
    316359     ! check kcount and kstart optionals parameters... 
    317      IF( PRESENT(kcount) .AND. (.NOT. PRESENT(kstart)) ) & 
     360     IF( PRESENT(kcount) .AND. (.NOT. PRESENT(kstart)) )   & 
    318361          CALL ctl_stop( trim(clinfo), 'kcount present needs kstart present' ) 
    319      IF( PRESENT(kstart) .AND. (.NOT. PRESENT(kcount)) ) & 
     362     IF( PRESENT(kstart) .AND. (.NOT. PRESENT(kcount)) )   & 
    320363          CALL ctl_stop( trim(clinfo), 'kstart present needs kcount present' ) 
    321      IF( PRESENT(kstart) .AND. idom /= jpdom_unknown ) & 
     364     IF( PRESENT(kstart) .AND. idom /= jpdom_unknown   )  & 
    322365          CALL ctl_stop( trim(clinfo), 'kstart present needs kdom = jpdom_unknown' ) 
    323366 
    324367     ! Search for the variable in the data base (eventually actualize data) 
    325      !- 
    326368     istop = nstop 
    327369     idvar = iom_varid( knumfl, cdvar ) 
    328370     ! 
    329      IF ( idvar > 0 ) THEN 
     371     IF( idvar > 0 ) THEN 
    330372        ! to write iom_file(knumfl)%dimsz in a shorter way ! 
    331373        idimsz(:) = iom_file(knumfl)%dimsz(:, idvar)  
    332         inbdim = iom_file(knumfl)%ndims(idvar)! number of dimensions in the file 
    333         idmspc = inbdim ! number of spatial dimensions in the file 
    334         IF( iom_file(knumfl)%luld(idvar) ) idmspc = inbdim - 1 
    335         IF( idmspc > 3 ) CALL ctl_stop(trim(clinfo), 'the file has more than 3 spatial dimensions',   & 
    336              &                    'this case is not coded...')  
    337         ! Identify the domain in case of jpdom_local definition 
    338         !- 
    339         IF( idom == jpdom_local ) THEN 
     374        inbdim = iom_file(knumfl)%ndims(idvar)            ! number of dimensions in the file 
     375        idmspc = inbdim                                   ! number of spatial dimensions in the file 
     376        IF( iom_file(knumfl)%luld(idvar) )   idmspc = inbdim - 1 
     377        IF( idmspc > 3 )   CALL ctl_stop(trim(clinfo),   & 
     378           &                    'the file has more than 3 spatial dimensions this case is not coded...' )  
     379        IF( idom == jpdom_local ) THEN        ! Identify the domain in case of jpdom_local definition 
    340380           IF( idimsz(1) == jpi .AND. idimsz(2) == jpj ) THEN  
    341381              idom = jpdom_local_full 
     
    348388           ENDIF 
    349389        ENDIF 
    350  
     390        ! 
    351391        ! definition of istart and icnt 
    352         !- 
     392        ! 
    353393        ! initializations 
    354394        istart(:) = 1 
    355         icnt(:) = 1 
     395        icnt  (:) = 1 
    356396        itime = 1 
    357397        IF( PRESENT(ktime) ) itime = ktime 
     
    383423           CASE (2) 
    384424              ! data is 2d array (+ maybe a temporal dimension) 
    385               IF ( PRESENT(kstart) ) THEN 
     425              IF( PRESENT(kstart) ) THEN 
    386426                 istart(1:3) = (/ kstart(1:2), itime /) 
    387427                 icnt(1:2) = kcount(1:2) 
     
    404444              ENDIF 
    405445           CASE DEFAULT 
    406               IF ( itime == 1 .AND. idimsz(3) == 1 .AND. idmspc == 3 ) THEN 
    407                  CALL ctl_warn( trim(clinfo), '2D array but 3 spatial dimensions for the data...', & 
    408                       &         'As the size of the z dimension is 1 and as we try to read the first reacord, ',     & 
    409                       &         'we accept this case even if there is a possible mix-up between z and time dimension...')            
    410                  IF ( PRESENT(kstart) ) THEN 
     446              IF( itime == 1 .AND. idimsz(3) == 1 .AND. idmspc == 3 ) THEN 
     447                 CALL ctl_warn( trim(clinfo), '2D array but 3 spatial dimensions for the data...',                & 
     448                      &         'As the size of the z dimension is 1 and as we try to read the first record, ',   & 
     449                      &         'we accept this case even if there is a possible mix-up between z and time dimension' )            
     450                 IF( PRESENT(kstart) ) THEN 
    411451                    istart(1:2) = kstart(1:2) 
    412452                    icnt(1:2) = kcount(1:2) 
     
    428468                 ENDIF 
    429469              ELSE 
    430                  CALL ctl_stop( trim(clinfo), 'to keep iom lisibility, when reading a 2D array,', & 
    431                       &         'we do not accept data with more than 2 spatial dimension',     & 
    432                       &         'Use ncwa -a to suppress the unnecessary dimensions')            
     470                 CALL ctl_stop( trim(clinfo), 'to keep iom lisibility, when reading a 2D array,',           & 
     471                      &                       'we do not accept data with more than 2 spatial dimension',   & 
     472                      &                       'Use ncwa -a to suppress the unnecessary dimensions' ) 
    433473              ENDIF 
    434474           END SELECT 
    435475        ELSEIF( PRESENT(pv_r3d) ) THEN 
    436476           SELECT CASE (idmspc) 
    437            CASE (1) 
    438               CALL ctl_stop(trim(clinfo), 'the file has only 1 spatial dimension',   & 
    439                    &        'it is impossible to read a 3d array from this file...') 
    440            CASE (2) 
    441               CALL ctl_stop(trim(clinfo), 'the file has only 2 spatial dimension',   & 
    442                    &        'it is impossible to read a 3d array from this file...') 
    443            CASE (3) 
     477           CASE( 1 ) 
     478              CALL ctl_stop( trim(clinfo), 'the file has only 1 spatial dimension',            & 
     479                   &                       'it is impossible to read a 3d array from this file...' ) 
     480           CASE( 2 ) 
     481              CALL ctl_stop( trim(clinfo), 'the file has only 2 spatial dimension',            & 
     482                   &                       'it is impossible to read a 3d array from this file...' ) 
     483           CASE( 3 ) 
    444484              ! data is 3d array (+ maybe a temporal dimension) 
    445485              IF( PRESENT(kstart) ) THEN 
     
    469509              ENDIF 
    470510           CASE DEFAULT 
    471               CALL ctl_stop( trim(clinfo), 'to keep iom lisibility, when reading a 3D array,', & 
    472                    &         'we do not accept data with more than 3 spatial dimension',     & 
    473                    &         'Use ncwa -a to suppress the unnecessary dimensions')            
     511              CALL ctl_stop( trim(clinfo), 'to keep iom lisibility, when reading a 3D array,',   & 
     512                   &         'we do not accept data with more than 3 spatial dimension',         & 
     513                   &         'Use ncwa -a to suppress the unnecessary dimensions' )            
    474514           END SELECT 
    475515        ENDIF 
     
    491531           itmp = size(pv_r1d) 
    492532           WRITE(ctmp1,*) 'size(pv_r1d): ', itmp, ' /= icnt(1): ', icnt(1) 
    493            IF( itmp /= icnt(1) ) CALL ctl_stop( trim(clinfo), ctmp1 ) 
     533           IF( itmp /= icnt(1) )   CALL ctl_stop( trim(clinfo), ctmp1 ) 
    494534        ELSEIF( PRESENT(pv_r2d) ) THEN 
    495535           DO jl = 1, 2 
     
    501541                 WRITE(ctmp1,*) 'size(pv_r2d(nldi:nlei,nldj:nlej), ',jl,'): ',itmp,' /= icnt(',jl,'): ',icnt(jl) 
    502542              ENDIF 
    503               IF( itmp /= icnt(jl) ) CALL ctl_stop( trim(clinfo), ctmp1 ) 
     543              IF( itmp /= icnt(jl) )   CALL ctl_stop( trim(clinfo), ctmp1 ) 
    504544           END DO 
    505545        ELSEIF( PRESENT(pv_r3d) ) THEN 
     
    512552                 WRITE(ctmp1,*) 'size(pv_r3d(nldi:nlei,nldj:nlej), ',jl,'): ',itmp,' /= icnt(',jl,'): ',icnt(jl) 
    513553              ENDIF 
    514               IF( itmp /= icnt(jl) ) CALL ctl_stop( trim(clinfo), ctmp1 ) 
     554              IF( itmp /= icnt(jl) )   CALL ctl_stop( trim(clinfo), ctmp1 ) 
    515555           END DO 
    516556        ENDIF 
     
    520560     !-      
    521561     IF( istop == nstop) THEN ! no additional errors until this point... 
    522         ! 
    523         istop = nstop 
    524562        ! 
    525563        zscf = iom_file(knumfl)%scf(idvar)      ! scale factor 
     
    529567           CALL fliogetv( knumfl, cdvar, pv_r1d(:), start=istart(1:inbdim), count=icnt(1:inbdim) ) 
    530568           !--- Apply scale_factor and offset 
    531            IF( zscf /= 1. ) pv_r1d(:) = pv_r1d(:) * zscf  
    532            IF( zofs /= 0. ) pv_r1d(:) = pv_r1d(:) + zofs 
     569           IF( zscf /= 1. )   pv_r1d(:) = pv_r1d(:) * zscf  
     570           IF( zofs /= 0. )   pv_r1d(:) = pv_r1d(:) + zofs 
    533571        ELSEIF( PRESENT(pv_r2d) ) THEN 
    534572           IF( idom /= jpdom_unknown ) THEN 
    535573              CALL fliogetv( knumfl, cdvar, pv_r2d(nldi:nlei,nldj:nlej), start=istart(1:inbdim), count=icnt(1:inbdim) ) 
    536574              !--- Apply scale_factor and offset 
    537               IF (zscf /= 1.) pv_r2d(nldi:nlei, nldj:nlej) = pv_r2d(nldi:nlei,nldj:nlej) * zscf 
    538               IF (zofs /= 0.) pv_r2d(nldi:nlei, nldj:nlej) = pv_r2d(nldi:nlei,nldj:nlej) + zofs 
     575!CDIR NOUNROLL 
     576              IF( zscf /= 1.)   pv_r2d(nldi:nlei, nldj:nlej) = pv_r2d(nldi:nlei,nldj:nlej) * zscf 
     577!CDIR NOUNROLL 
     578              IF( zofs /= 0.)   pv_r2d(nldi:nlei, nldj:nlej) = pv_r2d(nldi:nlei,nldj:nlej) + zofs 
    539579              !--- Fill the overlap areas and extra hallows (mpp) 
    540580              CALL lbc_lnk (pv_r2d,'Z',-999.,'no0') 
     
    542582              CALL fliogetv( knumfl, cdvar, pv_r2d(:,:), start=istart(1:inbdim), count=icnt(1:inbdim) ) 
    543583              !--- Apply scale_factor and offset 
    544               IF (zscf /= 1.) pv_r2d(:,:) = pv_r2d(:,:) * zscf 
    545               IF (zofs /= 0.) pv_r2d(:,:) = pv_r2d(:,:) + zofs 
     584!CDIR COLLAPSE 
     585              IF( zscf /= 1.)   pv_r2d(:,:) = pv_r2d(:,:) * zscf 
     586!CDIR COLLAPSE 
     587              IF( zofs /= 0.)   pv_r2d(:,:) = pv_r2d(:,:) + zofs 
    546588           ENDIF 
    547589        ELSEIF( PRESENT(pv_r3d) ) THEN 
    548590           IF( idom /= jpdom_unknown ) THEN 
    549               CALL fliogetv( knumfl, cdvar, pv_r3d(nldi:nlei,nldj:nlej,:), start=istart(1:inbdim), count=icnt(1:inbdim) ) 
     591              CALL fliogetv( knumfl, cdvar, pv_r3d(nldi:nlei,nldj:nlej,:), start=istart(1:inbdim),   & 
     592                 &                                                         count=icnt  (1:inbdim) ) 
    550593              !--- Apply scale_factor and offset 
    551               IF( zscf /= 1. ) pv_r3d(nldi:nlei,nldj:nlej,:) = pv_r3d(nldi:nlei,nldj:nlej,:) * zscf 
    552               IF( zofs /= 0. ) pv_r3d(nldi:nlei,nldj:nlej,:) = pv_r3d(nldi:nlei,nldj:nlej,:) + zofs 
     594!CDIR NOUNROLL 
     595              IF( zscf /= 1. )   pv_r3d(nldi:nlei,nldj:nlej,:) = pv_r3d(nldi:nlei,nldj:nlej,:) * zscf 
     596!CDIR NOUNROLL 
     597              IF( zofs /= 0. )   pv_r3d(nldi:nlei,nldj:nlej,:) = pv_r3d(nldi:nlei,nldj:nlej,:) + zofs 
    553598              !--- Fill the overlap areas and extra hallows (mpp) 
    554599              IF( icnt(3) == jpk )   CALL lbc_lnk( pv_r3d,'Z',-999.,'no0' ) ! this if could be removed with the new lbc_lnk ... 
     
    556601              CALL fliogetv( knumfl, cdvar, pv_r3d(:,:,:), start=istart(1:inbdim), count=icnt(1:inbdim) ) 
    557602              !--- Apply scale_factor and offset 
    558               IF (zscf /= 1.)   pv_r3d(:,:,:) = pv_r3d(:,:,:) * zscf 
    559               IF (zofs /= 0.)   pv_r3d(:,:,:) = pv_r3d(:,:,:) + zofs 
     603!CDIR COLLAPSE 
     604              IF( zscf /= 1.)   pv_r3d(:,:,:) = pv_r3d(:,:,:) * zscf 
     605!CDIR COLLAPSE 
     606              IF( zofs /= 0.)   pv_r3d(:,:,:) = pv_r3d(:,:,:) + zofs 
    560607           ENDIF 
    561608        ENDIF 
    562609        ! 
    563         IF( istop == nstop .AND. lwp )  & 
    564              &  WRITE(numout,*) ' read '//trim(cdvar)//' in '//trim(iom_file(knumfl)%name)//' ok' 
     610        IF( istop == nstop .AND. lwp )   & 
     611           WRITE(numout,*) '          read '//trim(cdvar)//' in '//trim(iom_file(knumfl)%name)//' ok' 
    565612     ENDIF 
    566613     ! 
    567    END SUBROUTINE iom_u_getv 
     614   END SUBROUTINE iom_get_123d 
    568615 
    569616    
    570617   SUBROUTINE iom_gettime( knumfl, cdvar, ptime ) 
    571      !!-------------------------------------------------------------------- 
    572      !!                   ***  SUBROUTINE  iom_close  *** 
    573      !! 
    574      !! ** Purpose : read the time axis cdvar in the file  
    575      !! 
    576      !! ** Method : 
    577      !! 
    578      !!-------------------------------------------------------------------- 
    579      INTEGER               , INTENT(in)  ::   knumfl   ! Identifier of the file to be closed 
    580      CHARACTER(len=*)      , INTENT(in)  ::   cdvar    ! time axis name 
    581      REAL(wp), DIMENSION(:), INTENT(out) ::   ptime    ! the time axis 
    582  
    583      INTEGER           :: idvar    ! id of the variable 
    584      CHARACTER(LEN=100) :: clinfo                    ! info character 
    585      !--------------------------------------------------------------------- 
    586      clinfo = 'iom_gettime, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    587      idvar = iom_varid( knumfl, cdvar ) 
    588      ! 
    589      ptime(:) = 0. ! default definition 
    590      IF ( idvar > 0 ) THEN 
    591         IF ( iom_file(knumfl)%ndims(idvar) == 1 ) THEN 
    592            IF ( iom_file(knumfl)%luld(idvar) ) THEN 
    593               IF ( iom_file(knumfl)%dimsz(1,idvar) == size(ptime) ) THEN 
    594                  CALL fliogetv( knumfl, cdvar, ptime(:), start=(/ 1 /), & 
    595                       &                                   count=(/ iom_file(knumfl)%dimsz(1,idvar) /) ) 
    596               ELSE 
    597                  WRITE(ctmp1,*) 'error with the size of ptime ',size(ptime),iom_file(knumfl)%dimsz(1,idvar) 
    598                  CALL ctl_stop( trim(clinfo), trim(ctmp1) ) 
    599               ENDIF 
    600            ELSE 
    601               CALL ctl_stop( trim(clinfo), 'variable dimension is not unlimited... use iom_get' ) 
    602            ENDIF 
    603         ELSE 
    604            CALL ctl_stop( trim(clinfo), 'the variable has more than 1 dimension' ) 
    605         ENDIF 
    606      ELSE 
    607         CALL ctl_stop( trim(clinfo), 'variable not found in '//iom_file(knumfl)%name ) 
    608      ENDIF 
    609  
     618      !!-------------------------------------------------------------------- 
     619      !!                   ***  SUBROUTINE iom_gettime  *** 
     620      !! 
     621      !! ** Purpose : read the time axis cdvar in the file  
     622      !!-------------------------------------------------------------------- 
     623      INTEGER               , INTENT(in   ) ::   knumfl   ! file Identifier 
     624      CHARACTER(len=*)      , INTENT(in   ) ::   cdvar    ! time axis name 
     625      REAL(wp), DIMENSION(:), INTENT(  out) ::   ptime    ! the time axis 
     626      ! 
     627      INTEGER            ::   idvar    ! id of the variable 
     628      CHARACTER(LEN=100) ::   clinfo   ! info character 
     629      !--------------------------------------------------------------------- 
     630      ! 
     631      clinfo = 'iom_gettime, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
     632      idvar = iom_varid( knumfl, cdvar ) 
     633      ! 
     634      ptime(:) = 0. ! default definition 
     635      IF( idvar > 0 ) THEN 
     636         IF( iom_file(knumfl)%ndims(idvar) == 1 ) THEN 
     637            IF( iom_file(knumfl)%luld(idvar) ) THEN 
     638               IF( iom_file(knumfl)%dimsz(1,idvar) == size(ptime) ) THEN 
     639                  CALL fliogetv( knumfl, cdvar, ptime(:), start=(/ 1 /),   & 
     640                     &                                    count=(/ iom_file(knumfl)%dimsz(1,idvar) /) ) 
     641               ELSE 
     642                  WRITE(ctmp1,*) 'error with the size of ptime ',size(ptime),iom_file(knumfl)%dimsz(1,idvar) 
     643                  CALL ctl_stop( trim(clinfo), trim(ctmp1) ) 
     644               ENDIF 
     645            ELSE 
     646               CALL ctl_stop( trim(clinfo), 'variable dimension is not unlimited... use iom_get' ) 
     647            ENDIF 
     648         ELSE 
     649            CALL ctl_stop( trim(clinfo), 'the variable has more than 1 dimension' ) 
     650         ENDIF 
     651      ELSE 
     652         CALL ctl_stop( trim(clinfo), 'variable not found in '//iom_file(knumfl)%name ) 
     653      ENDIF 
     654      ! 
    610655   END SUBROUTINE iom_gettime 
    611656     
    612657     
    613      
    614658   FUNCTION iom_varid ( knumfl, cdvar, kdimsz )   
    615      !!----------------------------------------------------------------------- 
    616      !!                  ***  FUNCTION  iom_varid  *** 
    617      !! 
    618      !! ** Purpose : get the id of a variable in a file (return 0 if not found) 
    619      !! 
    620      !! ** Method : ??? 
    621      !! 
    622      !!----------------------------------------------------------------------- 
    623      INTEGER              , INTENT(in)            ::   knumfl   ! file Identifier 
    624      CHARACTER(len=*)     , INTENT(in)            ::   cdvar    ! name of the variable 
    625      INTEGER, DIMENSION(:), INTENT(out), OPTIONAL ::   kdimsz   ! size of the dimensions 
    626      ! 
    627      INTEGER                        :: iom_varid, idvar, i_nvd, ji 
    628      INTEGER, DIMENSION(jpmax_dims) :: idimid 
    629      LOGICAL                        :: ll_fnd 
    630      CHARACTER(LEN=100)             :: clinfo                   ! info character 
    631      !!----------------------------------------------------------------------- 
    632      clinfo = 'iom_varid, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
    633      iom_varid = 0                         ! default definition 
    634      IF ( PRESENT(kdimsz) ) kdimsz(:) = 0  ! default definition 
    635      ! 
    636      IF ( knumfl > 0 ) THEN 
    637         IF( iom_file(knumfl)%iopen == 0 ) THEN  
    638            CALL ctl_stop( trim(clinfo), 'the file is not open' ) 
    639         ELSE 
    640            ! 
    641            ll_fnd  = .FALSE. 
    642            idvar = 0 
    643            ! 
    644            DO WHILE ( .NOT.ll_fnd .AND. idvar < iom_file(knumfl)%nvars ) 
    645               idvar = idvar + 1 
    646               ll_fnd  = ( TRIM(cdvar) == TRIM(iom_file(knumfl)%cn_var(idvar)) ) 
    647            END DO 
    648            ! 
    649            IF( .NOT.ll_fnd ) THEN 
    650               idvar = idvar + 1 
    651               IF( idvar <= jpmax_vars ) THEN 
    652                  CALL flioinqv( knumfl, cdvar, ll_fnd, nb_dims = i_nvd ) 
    653                  IF( ll_fnd ) THEN 
    654                     IF( i_nvd <= jpmax_dims ) THEN 
    655                        iom_file(knumfl)%nvars           = idvar 
    656                        iom_file(knumfl)%cn_var(idvar) = trim(cdvar) 
    657                        iom_file(knumfl)%ndims(idvar)  = i_nvd 
    658                        CALL flioinqv( knumfl, cdvar, ll_fnd, & 
    659                             &         len_dims = iom_file(knumfl)%dimsz(1:i_nvd,idvar), & 
    660                             &         id_dims  = idimid(1:i_nvd) ) 
    661                        DO ji = 1, i_nvd 
    662                           IF ( idimid(ji) == iom_file(knumfl)%iduld ) iom_file(knumfl)%luld(idvar) = .TRUE. 
    663                        END DO 
    664                        !---------- 
    665                        !---------- Deal with scale_factor and offset 
    666                        CALL flioinqa( knumfl, cdvar, 'scale_factor', ll_fnd ) 
    667                        IF (ll_fnd) THEN 
    668                           CALL fliogeta( knumfl, cdvar, 'scale_factor', iom_file(knumfl)%scf(idvar) ) 
    669                        ELSE 
    670                           iom_file(knumfl)%scf(idvar) = 1. 
    671                        END IF 
    672                        CALL flioinqa( knumfl, cdvar, 'offset', ll_fnd ) 
    673                        IF( ll_fnd ) THEN 
    674                           CALL fliogeta( knumfl, cdvar, 'offset', iom_file(knumfl)%ofs(idvar) ) 
    675                        ELSE 
    676                           iom_file(knumfl)%ofs(idvar) = 0. 
    677                        END IF 
    678                        ! 
    679                        iom_varid = idvar 
    680                        IF ( PRESENT(kdimsz) ) THEN  
    681                           IF ( i_nvd == size(kdimsz) ) THEN 
    682                              kdimsz(:) = iom_file(knumfl)%dimsz(1:i_nvd,idvar) 
    683                           ELSE 
    684                              WRITE(ctmp1,*) i_nvd, size(kdimsz) 
    685                              CALL ctl_stop( trim(clinfo), 'error in kdimsz size'//trim(ctmp1) ) 
    686                           ENDIF 
    687                        ENDIF 
    688                     ELSE 
    689                        CALL ctl_stop( trim(clinfo), 'Too many dimensions in the file '//iom_file(knumfl)%name,   & 
    690                             &           'increase the parameter jpmax_vars') 
    691                     ENDIF 
    692                  ELSE   
    693                     CALL ctl_warn( trim(clinfo), 'Variable '//trim(cdvar)// & 
    694                          & ' is not found in the file '//trim(iom_file(knumfl)%name) ) 
    695                  ENDIF 
    696               ELSE 
    697                  CALL ctl_stop( trim(clinfo), 'Too many variables in the file '//iom_file(knumfl)%name,   & 
    698                       &           'increase the parameter jpmax_vars') 
    699               ENDIF 
    700            ELSE 
    701               iom_varid = idvar 
    702               IF ( PRESENT(kdimsz) ) THEN  
    703                  i_nvd = iom_file(knumfl)%ndims(idvar) 
    704                  IF ( i_nvd == size(kdimsz) ) THEN 
    705                     kdimsz(:) = iom_file(knumfl)%dimsz(1:i_nvd,idvar) 
    706                  ELSE 
    707                     WRITE(ctmp1,*) i_nvd, size(kdimsz) 
    708                     CALL ctl_stop( trim(clinfo), 'error in kdimsz size'//trim(ctmp1) ) 
    709                  ENDIF 
    710               ENDIF 
    711            ENDIF 
    712         ENDIF 
    713      ENDIF 
    714  
     659      !!----------------------------------------------------------------------- 
     660      !!                  ***  FUNCTION  iom_varid  *** 
     661      !! 
     662      !! ** Purpose : get the id of a variable in a file (return 0 if not found) 
     663      !!----------------------------------------------------------------------- 
     664      INTEGER              , INTENT(in   )           ::   knumfl   ! file Identifier 
     665      CHARACTER(len=*)     , INTENT(in   )           ::   cdvar    ! name of the variable 
     666      INTEGER, DIMENSION(:), INTENT(  out), OPTIONAL ::   kdimsz   ! size of the dimensions 
     667      ! 
     668      INTEGER                        ::   ji                       ! dummy loop index 
     669      INTEGER                        ::   iom_varid, idvar, i_nvd 
     670      INTEGER, DIMENSION(jpmax_dims) ::   idimid 
     671      LOGICAL                        ::   ll_fnd 
     672      CHARACTER(LEN=100)             ::   clinfo                   ! info character 
     673      !!----------------------------------------------------------------------- 
     674      clinfo = 'iom_varid, file: '//trim(iom_file(knumfl)%name)//', var: '//trim(cdvar) 
     675      iom_varid = 0                         ! default definition 
     676      IF( PRESENT(kdimsz) ) kdimsz(:) = 0   ! default definition 
     677      ! 
     678      IF( knumfl > 0 ) THEN 
     679         IF( iom_file(knumfl)%iopen == 0 ) THEN  
     680            CALL ctl_stop( trim(clinfo), 'the file is not open' ) 
     681         ELSE 
     682            ll_fnd  = .FALSE. 
     683            idvar = 0 
     684            ! 
     685            DO WHILE ( .NOT.ll_fnd .AND. idvar < iom_file(knumfl)%nvars ) 
     686               idvar = idvar + 1 
     687               ll_fnd  = ( TRIM(cdvar) == TRIM(iom_file(knumfl)%cn_var(idvar)) ) 
     688            END DO 
     689            ! 
     690            IF( .NOT.ll_fnd ) THEN 
     691               idvar = idvar + 1 
     692               IF( idvar <= jpmax_vars ) THEN 
     693                  CALL flioinqv( knumfl, cdvar, ll_fnd, nb_dims = i_nvd ) 
     694                  IF( ll_fnd ) THEN 
     695                     IF( i_nvd <= jpmax_dims ) THEN 
     696                        iom_file(knumfl)%nvars           = idvar 
     697                        iom_file(knumfl)%cn_var(idvar) = trim(cdvar) 
     698                        iom_file(knumfl)%ndims(idvar)  = i_nvd 
     699                        CALL flioinqv( knumfl, cdvar, ll_fnd,   & 
     700                           &           len_dims = iom_file(knumfl)%dimsz(1:i_nvd,idvar), & 
     701                           &           id_dims  = idimid(1:i_nvd) ) 
     702                        DO ji = 1, i_nvd 
     703                           IF( idimid(ji) == iom_file(knumfl)%iduld ) iom_file(knumfl)%luld(idvar) = .TRUE. 
     704                        END DO 
     705                        !---------- 
     706                        !---------- Deal with scale_factor and offset 
     707                        CALL flioinqa( knumfl, cdvar, 'scale_factor', ll_fnd ) 
     708                        IF( ll_fnd) THEN 
     709                           CALL fliogeta( knumfl, cdvar, 'scale_factor', iom_file(knumfl)%scf(idvar) ) 
     710                        ELSE 
     711                           iom_file(knumfl)%scf(idvar) = 1. 
     712                        END IF 
     713                        CALL flioinqa( knumfl, cdvar, 'offset', ll_fnd ) 
     714                        IF( ll_fnd ) THEN 
     715                           CALL fliogeta( knumfl, cdvar, 'offset', iom_file(knumfl)%ofs(idvar) ) 
     716                        ELSE 
     717                           iom_file(knumfl)%ofs(idvar) = 0. 
     718                        END IF 
     719                        ! 
     720                        iom_varid = idvar 
     721                        IF( PRESENT(kdimsz) ) THEN  
     722                           IF( i_nvd == size(kdimsz) ) THEN 
     723                              kdimsz(:) = iom_file(knumfl)%dimsz(1:i_nvd,idvar) 
     724                           ELSE 
     725                              WRITE(ctmp1,*) i_nvd, size(kdimsz) 
     726                              CALL ctl_stop( trim(clinfo), 'error in kdimsz size'//trim(ctmp1) ) 
     727                           ENDIF 
     728                        ENDIF 
     729                     ELSE 
     730                        CALL ctl_stop( trim(clinfo), 'Too many dimensions in the file '//iom_file(knumfl)%name,   & 
     731                           &                        'increase the parameter jpmax_vars') 
     732                     ENDIF 
     733!!$                  ELSE   
     734!!$                     CALL ctl_warn( trim(clinfo), 'Variable '//trim(cdvar)// & 
     735!!$                        &                         ' is not found in the file '//trim(iom_file(knumfl)%name) ) 
     736                  ENDIF 
     737               ELSE 
     738                  CALL ctl_stop( trim(clinfo), 'Too many variables in the file '//iom_file(knumfl)%name,   & 
     739                     &                         'increase the parameter jpmax_vars') 
     740               ENDIF 
     741            ELSE 
     742               iom_varid = idvar 
     743               IF( PRESENT(kdimsz) ) THEN  
     744                  i_nvd = iom_file(knumfl)%ndims(idvar) 
     745                  IF( i_nvd == size(kdimsz) ) THEN 
     746                     kdimsz(:) = iom_file(knumfl)%dimsz(1:i_nvd,idvar) 
     747                  ELSE 
     748                     WRITE(ctmp1,*) i_nvd, size(kdimsz) 
     749                     CALL ctl_stop( trim(clinfo), 'error in kdimsz size'//trim(ctmp1) ) 
     750                  ENDIF 
     751               ENDIF 
     752            ENDIF 
     753         ENDIF 
     754      ENDIF 
     755      ! 
    715756   END FUNCTION iom_varid 
    716757 
    717    
    718    FUNCTION iom_get_gblatt( knumfl, kinfonum ) 
    719       !!----------------------------------------------------------------------- 
    720       !!                  ***  FUNCTION  iom_get_gblatt  *** 
     758   !!---------------------------------------------------------------------- 
     759   !!                   INTERFACE iom_rstput 
     760   !!---------------------------------------------------------------------- 
     761   SUBROUTINE iom_rstput_0d( kt, kwrite, knumfl, cdvar, pvar ) 
     762      INTEGER         , INTENT(in)                         ::   kt       ! ocean time-step 
     763      INTEGER         , INTENT(in)                         ::   kwrite   ! writing time-step 
     764      INTEGER         , INTENT(in)                         ::   knumfl   ! Identifier of the file  
     765      CHARACTER(len=*), INTENT(in)                         ::   cdvar    ! time axis name 
     766      REAL(wp)        , INTENT(in)                         ::   pvar     ! read field 
     767      IF( knumfl > 0 )   CALL iom_rstput_0123d( kt, kwrite, knumfl, cdvar, pv_r0d = pvar ) 
     768   END SUBROUTINE iom_rstput_0d 
     769 
     770   SUBROUTINE iom_rstput_1d( kt, kwrite, knumfl, cdvar, pvar ) 
     771      INTEGER         , INTENT(in)                         ::   kt       ! ocean time-step 
     772      INTEGER         , INTENT(in)                         ::   kwrite   ! writing time-step 
     773      INTEGER         , INTENT(in)                         ::   knumfl   ! Identifier of the file  
     774      CHARACTER(len=*), INTENT(in)                         ::   cdvar    ! time axis name 
     775      REAL(wp)        , INTENT(in), DIMENSION(        jpk) ::   pvar     ! read field 
     776      IF( knumfl > 0 )    CALL iom_rstput_0123d( kt, kwrite, knumfl, cdvar, pv_r1d = pvar ) 
     777   END SUBROUTINE iom_rstput_1d 
     778 
     779   SUBROUTINE iom_rstput_2d( kt, kwrite, knumfl, cdvar, pvar ) 
     780      INTEGER         , INTENT(in)                         ::   kt       ! ocean time-step 
     781      INTEGER         , INTENT(in)                         ::   kwrite   ! writing time-step 
     782      INTEGER         , INTENT(in)                         ::   knumfl   ! Identifier of the file  
     783      CHARACTER(len=*), INTENT(in)                         ::   cdvar    ! time axis name 
     784      REAL(wp)        , INTENT(in), DIMENSION(jpi,jpj    ) ::   pvar     ! read field 
     785      IF( knumfl > 0 )   CALL iom_rstput_0123d( kt, kwrite, knumfl, cdvar, pv_r2d = pvar ) 
     786   END SUBROUTINE iom_rstput_2d 
     787 
     788   SUBROUTINE iom_rstput_3d( kt, kwrite, knumfl, cdvar, pvar ) 
     789      INTEGER         , INTENT(in)                         ::   kt       ! ocean time-step 
     790      INTEGER         , INTENT(in)                         ::   kwrite   ! writing time-step 
     791      INTEGER         , INTENT(in)                         ::   knumfl   ! Identifier of the file  
     792      CHARACTER(len=*), INTENT(in)                         ::   cdvar    ! time axis name 
     793      REAL(wp)        , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pvar     ! read field 
     794      IF( knumfl > 0 )   CALL iom_rstput_0123d( kt, kwrite, knumfl, cdvar, pv_r3d = pvar ) 
     795   END SUBROUTINE iom_rstput_3d 
     796   !!---------------------------------------------------------------------- 
     797 
     798   SUBROUTINE iom_rstput_0123d( kt, kwrite, knumfl, cdvar ,   & 
     799      &                             pv_r0d, pv_r1d, pv_r2d, pv_r3d ) 
     800      !!-------------------------------------------------------------------- 
     801      !!                   ***  SUBROUTINE  iom_rstput  *** 
    721802      !! 
    722       !! ** Purpose : ??? 
    723       !! 
    724       !! ** Method : ??? 
    725       !! 
    726       !!----------------------------------------------------------------------- 
    727       INTEGER,INTENT(in) :: knumfl 
    728       INTEGER,          intent(in) :: kinfonum 
    729       ! 
    730       CHARACTER(LEN=10) :: clinfo 
    731       REAL(wp) :: iom_get_gblatt 
    732       !!----------------------------------------------------------------------- 
    733  
    734       WRITE(clinfo,*) kinfonum 
    735       clinfo = 'info'//trim(adjustl(clinfo)) 
    736       CALL fliogeta (knumfl, "?", clinfo, iom_get_gblatt) 
    737  
    738    END FUNCTION iom_get_gblatt 
    739  
     803      !! ** Purpose : read the time axis cdvar in the file  
     804      !!-------------------------------------------------------------------- 
     805      INTEGER                   , INTENT(in)           ::   kt       ! ocean time-step 
     806      INTEGER                   , INTENT(in)           ::   kwrite   ! writing time-step 
     807      INTEGER                   , INTENT(in)           ::   knumfl   ! Identifier of the file  
     808      CHARACTER(len=*)          , INTENT(in)           ::   cdvar    ! time axis name 
     809      REAL(wp)                  , INTENT(in), OPTIONAL ::   pv_r0d   ! read field 
     810      REAL(wp), DIMENSION(:)    , INTENT(in), OPTIONAL ::   pv_r1d   ! read field 
     811      REAL(wp), DIMENSION(:,:)  , INTENT(in), OPTIONAL ::   pv_r2d   ! read field 
     812      REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL ::   pv_r3d   ! read field 
     813      ! 
     814      INTEGER               :: idims, idvar         
     815      INTEGER               :: ix1, ix2, iy1, iy2        
     816      INTEGER, DIMENSION(4) :: idimsz, idimid      
     817      CHARACTER(LEN=100)    :: clinfo        ! info character 
     818      !--------------------------------------------------------------------- 
     819      ! 
     820      clinfo = '          iom_rstput_0123d, file: '//TRIM(iom_file(knumfl)%name)//', var: '//TRIM(cdvar) 
     821 
     822      ! define dimension variables if it is not already done 
     823      ! ========================== 
     824      IF( iom_file(knumfl)%nvars == 0 ) THEN 
     825         ! define the dimension variables if it is not already done 
     826         CALL fliodefv( knumfl,'nav_lon', (/1,2/), v_t=flio_r4   , axis='X',   & 
     827            &                  long_name="Longitude", units="degrees_east" ) 
     828         CALL fliodefv( knumfl,'nav_lat', (/1,2/), v_t=flio_r4   , axis='Y',   & 
     829            &                  long_name="Latitude", units="degrees_north" ) 
     830         CALL fliodefv( knumfl,'nav_lev', (/3/)  , v_t=flio_i4   , axis='Z',   & 
     831            &                  long_name="Model levels",units="model_levels") 
     832         CALL fliodefv( knumfl,'time_counter', (/4/), v_t=flio_r4, axis='T',   & 
     833            &                  long_name="Time axis", units='seconds since 0001-01-01 00:00:00' ) 
     834         ! update informations structure related the dimension variable we just added... 
     835         iom_file(knumfl)%nvars       = 4 
     836         iom_file(knumfl)%luld(1:4)   = (/ .FALSE., .FALSE., .FALSE., .TRUE. /) 
     837         iom_file(knumfl)%cn_var(1:3) = (/ 'nav_lon', 'nav_lat', 'nav_lev' /) 
     838         iom_file(knumfl)%cn_var(4)   = 'time_counter' 
     839         iom_file(knumfl)%ndims(1:4)  = (/ 2, 2, 1, 1 /)   
     840         IF(lwp) WRITE(numout,*) TRIM(clinfo)//' define dimension variables done' 
     841      ENDIF 
     842 
     843      ! define the data if it is not already done 
     844      ! =============== 
     845      idvar = iom_varid( knumfl, cdvar ) 
     846      IF( idvar <= 0 ) THEN 
     847         ! variable definition 
     848         IF(     PRESENT(pv_r0d) ) THEN   ;   idims = 0 
     849         ELSEIF( PRESENT(pv_r1d) ) THEN   ;   idims = 2   ;   idimid(1:idims) = (/    3,4/) 
     850         ELSEIF( PRESENT(pv_r2d) ) THEN   ;   idims = 3   ;   idimid(1:idims) = (/1,2  ,4/) 
     851         ELSEIF( PRESENT(pv_r3d) ) THEN   ;   idims = 4   ;   idimid(1:idims) = (/1,2,3,4/) 
     852         ENDIF 
     853         IF( PRESENT(pv_r0d) ) THEN   ;   CALL fliodefv (knumfl, cdvar                 , v_t = flio_r8) 
     854         ELSE                         ;   CALL fliodefv (knumfl, cdvar, idimid(1:idims), v_t = flio_r8) 
     855         ENDIF 
     856         ! update informations structure related the new variable we want to add... 
     857         idvar                          = iom_file(knumfl)%nvars + 1 
     858         iom_file(knumfl)%nvars         = idvar 
     859         iom_file(knumfl)%cn_var(idvar) = TRIM(cdvar) 
     860         iom_file(knumfl)%scf(idvar)    = 1. 
     861         iom_file(knumfl)%ofs(idvar)    = 0. 
     862         iom_file(knumfl)%ndims(idvar)  = idims 
     863         IF( .NOT. PRESENT(pv_r0d) ) THEN 
     864            iom_file(knumfl)%luld(idvar) = .TRUE. 
     865            CALL flioinqf( knumfl, ln_dim = idimsz ) 
     866            iom_file(knumfl)%dimsz(1:idims-1,idvar) = idimsz(idimid(1:idims-1)) 
     867         ENDIF 
     868         IF(lwp) WRITE(numout,*) TRIM(clinfo)//' defined ok' 
     869      ENDIF 
     870 
     871      ! time step kwrite : write the variable 
     872      IF( kt == kwrite ) THEN 
     873         ! on what kind of domain must the data be written? 
     874         IF( PRESENT(pv_r2d) .OR. PRESENT(pv_r3d) ) THEN 
     875            idimsz(1:2) = iom_file(knumfl)%dimsz(1:2,idvar) 
     876            IF(     idimsz(1) == (nlei - nldi + 1) .AND. idimsz(2) == (nlej - nldj + 1) ) THEN 
     877               ix1 = nldi   ;   ix2 = nlei   ;   iy1 = nldj   ;   iy2 = nlej 
     878            ELSEIF( idimsz(1) == nlci              .AND. idimsz(2) == nlcj              ) THEN 
     879               ix1 = 1      ;   ix2 = nlci   ;   iy1 = 1      ;   iy2 = nlcj 
     880            ELSEIF( idimsz(1) == jpi               .AND. idimsz(2) == jpj               ) THEN 
     881               ix1 = 1      ;   ix2 = jpi    ;   iy1 = 1      ;   iy2 = jpj 
     882            ELSE  
     883               CALL ctl_stop( 'iom_rstput_0123d: should have been an impossible case...' ) 
     884            ENDIF 
     885 
     886            ! write dimension variables if it is not already done 
     887            ! ============= 
     888            IF( iom_file(knumfl)%dimsz(1, 1) == 0 ) THEN 
     889               CALL flioputv( knumfl, 'nav_lon'     , glamt(ix1:ix2, iy1:iy2) ) 
     890               CALL flioputv( knumfl, 'nav_lat'     , gphit(ix1:ix2, iy1:iy2) ) 
     891               CALL flioputv( knumfl, 'nav_lev'     , gdept_0 ) 
     892               CALL flioputv( knumfl, 'time_counter', kt )   ! +++ WRONG VALUE: to be improved but not really useful... 
     893               ! update informations structure related the dimension variable 
     894               iom_file(knumfl)%dimsz(1:2, 1) = idimsz(1:2) 
     895               iom_file(knumfl)%dimsz(1:2, 2) = idimsz(1:2) 
     896               iom_file(knumfl)%dimsz(1, 3:4) = (/idimsz(3), 1/) 
     897               IF(lwp) WRITE(numout,*) TRIM(clinfo)//' write dimension variables done' 
     898            ENDIF 
     899         ENDIF 
     900 
     901         ! write the data 
     902         ! ============= 
     903         IF(     PRESENT(pv_r0d) ) THEN   ;   CALL flioputv( knumfl, cdvar, pv_r0d                      ) 
     904         ELSEIF( PRESENT(pv_r1d) ) THEN   ;   CALL flioputv( knumfl, cdvar, pv_r1d(                  :) ) 
     905         ELSEIF( PRESENT(pv_r2d) ) THEN   ;   CALL flioputv( knumfl, cdvar, pv_r2d(ix1:ix2, iy1:iy2   ) ) 
     906         ELSEIF( PRESENT(pv_r3d) ) THEN   ;   CALL flioputv( knumfl, cdvar, pv_r3d(ix1:ix2, iy1:iy2, :) ) 
     907         ENDIF 
     908         ! add 1 to the size of the temporal dimension (not really useful...) 
     909         IF( iom_file(knumfl)%luld(idvar) )   iom_file(knumfl)%dimsz(iom_file(knumfl)%ndims(idvar), idvar)    & 
     910            &                               = iom_file(knumfl)%dimsz(iom_file(knumfl)%ndims(idvar), idvar) + 1 
     911         IF(lwp) WRITE(numout,*) TRIM(clinfo)//' written ok' 
     912      ENDIF 
     913      !      
     914   END SUBROUTINE iom_rstput_0123d 
     915    
    740916   !!====================================================================== 
    741917END MODULE iom 
  • trunk/NEMO/OPA_SRC/istate.F90

    r479 r508  
    44   !! Ocean state   :  initial state setting 
    55   !!===================================================================== 
     6   !! History :   4.0  !  89-12  (P. Andrich)  Original code 
     7   !!             5.0  !  91-11  (G. Madec)  rewritting 
     8   !!             6.0  !  96-01  (G. Madec)  terrain following coordinates 
     9   !!             8.0  !  01-09  (M. Levy, M. Ben Jelloul)  istate_eel 
     10   !!             8.0  !  01-09  (M. Levy, M. Ben Jelloul)  istate_uvg 
     11   !!             9.0  !  03-08  (G. Madec)  F90: Free form, modules 
     12   !!             9.0  !  03-09  (G. Madec, C. Talandier)  add EEL R5 
     13   !!             9.0  !  04-05  (A. Koch-Larrouy)  istate_gyre  
     14   !!             9.0  !  06-07  (S. Masson)  distributed restart using iom 
     15   !!---------------------------------------------------------------------- 
    616 
    717   !!---------------------------------------------------------------------- 
     
    1323   !!   istate_uvg    : initial velocity in geostropic balance 
    1424   !!---------------------------------------------------------------------- 
    15    !! * Modules used 
    1625   USE oce             ! ocean dynamics and active tracers  
    1726   USE dom_oce         ! ocean space and time domain  
     
    1928   USE ldftra_oce      ! ocean active tracers: lateral physics 
    2029   USE zdf_oce         ! ocean vertical physics 
    21    USE in_out_manager  ! I/O manager 
    2230   USE phycst          ! physical constants 
    2331   USE wzvmod          ! verctical velocity               (wzv     routine) 
     
    2634   USE restart         ! ocean restart                   (rst_read routine) 
    2735   USE solisl          ! ??? 
    28  
     36   USE in_out_manager  ! I/O manager 
     37   USE iom 
     38    
    2939   IMPLICIT NONE 
    3040   PRIVATE 
    3141 
    32    !! * Routine accessibility 
    33    PUBLIC istate_init   ! routine called by step.F90 
     42   PUBLIC   istate_init   ! routine called by step.F90 
    3443 
    3544   !! * Substitutions 
     
    3746#  include "vectopt_loop_substitute.h90" 
    3847   !!---------------------------------------------------------------------- 
    39    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
     48   !!   OPA 9.0 , LOCEAN-IPSL (2006)  
    4049   !! $Header$  
    41    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     50   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4251   !!---------------------------------------------------------------------- 
    4352 
     
    4857      !!                   ***  ROUTINE istate_init  *** 
    4958      !!  
    50       !! ** Purpose :   Initialization of the dynamics and tracers. 
    51       !! 
    52       !! ** Method  : 
    53       !! 
    54       !! History : 
    55       !!   4.0  !  91-03  ()  Original code 
    56       !!        !  91-11  (G. Madec) 
    57       !!   9.0  !  03-09  (G. Madec)  F90: Free form, modules, orthogonality 
    58       !!---------------------------------------------------------------------- 
    59       USE iom 
    60       !! * Local declarations 
    61       !CT INTEGER ::   inum 
    62       !!---------------------------------------------------------------------- 
    63  
    64  
    65       ! Initialization to zero 
    66       ! ---------------------- 
    67  
    68       !     before fields       !       now fields        !      after fields       ! 
    69       ;   ub   (:,:,:) = 0.e0   ;   un   (:,:,:) = 0.e0   ;   ua   (:,:,:) = 0.e0 
    70       ;   vb   (:,:,:) = 0.e0   ;   vn   (:,:,:) = 0.e0   ;   va   (:,:,:) = 0.e0 
    71       ;                         ;   wn   (:,:,:) = 0.e0   ; 
    72       ;   rotb (:,:,:) = 0.e0   ;   rotn (:,:,:) = 0.e0   ; 
    73       ;   hdivb(:,:,:) = 0.e0   ;   hdivn(:,:,:) = 0.e0   ; 
    74  
    75       ;   tb   (:,:,:) = 0.e0   ;   tn   (:,:,:) = 0.e0   ;   ta   (:,:,:) = 0.e0 
    76       ;   sb   (:,:,:) = 0.e0   ;   sn   (:,:,:) = 0.e0   ;   sa   (:,:,:) = 0.e0 
     59      !! ** Purpose :   Initialization of the dynamics and tracer fields. 
     60      !!---------------------------------------------------------------------- 
     61 
     62      IF(lwp) WRITE(numout,*) 
     63      IF(lwp) WRITE(numout,*) 'istate_ini : Initialization of the dynamics and tracers' 
     64      IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    7765 
    7866      rhd  (:,:,:) = 0.e0 
    7967      rhop (:,:,:) = 0.e0 
    8068      rn2  (:,:,:) = 0.e0  
    81  
    82 #if defined key_dynspg_rl 
    83       ! rigid-lid formulation 
    84       bsfb(:,:) = 0.e0      ! before barotropic stream-function 
    85       bsfn(:,:) = 0.e0      ! now    barotropic stream-function 
    86       bsfd(:,:) = 0.e0      ! barotropic stream-function trend 
    87 #endif 
    88       ! free surface formulation 
    89       sshb(:,:) = 0.e0      ! before sea-surface height 
    90       sshn(:,:) = 0.e0      ! now    sea-surface height 
    91  
    9269 
    9370      IF( ln_rstart ) THEN                    ! Restart from a file 
     
    10077         neuler = 0                              ! Set time-step indicator at nit000 (euler forward) 
    10178         adatrj = 0._wp 
     79         !                                       ! Initialization of ocean to zero 
     80         !     before fields       !       now fields           
     81         ;   ub   (:,:,:) = 0.e0   ;   un   (:,:,:) = 0.e0    
     82         ;   vb   (:,:,:) = 0.e0   ;   vn   (:,:,:) = 0.e0     
     83         ;   rotb (:,:,:) = 0.e0   ;   rotn (:,:,:) = 0.e0   
     84         ;   hdivb(:,:,:) = 0.e0   ;   hdivn(:,:,:) = 0.e0   
     85         ! 
    10286         IF( cp_cfg == 'eel' ) THEN 
    10387            CALL istate_eel                      ! EEL   configuration : start from pre-defined 
     
    10791            !                                    !                       and salinity fields  
    10892         ELSE 
    109          !                                       ! Other configurations: Initial temperature and salinity fields 
    110  
    111          !CT CALL iom_open ('ssh', inum)  
    112          !CT CALL iom_get( inum, jpdom_local, 'sshb', sshb )     ! free surface formulation (ssh) 
    113          !CT sshn(:,:) = sshb(:,:) 
    114          !CT CALL iom_close (inum) 
    115  
     93            !                                    ! Other configurations: Initial temperature and salinity fields 
    11694#if defined key_dtatem 
    11795            CALL dta_tem( nit000 )                  ! read 3D temperature data 
     
    12098#else 
    12199            IF(lwp) WRITE(numout,*)                 ! analytical temperature profile 
    122             IF(lwp) WRITE(numout,*)' Temperature initialization using an analytic profile' 
     100            IF(lwp) WRITE(numout,*)'             Temperature initialization using an analytic profile' 
    123101            CALL istate_tem 
    124102#endif 
     
    130108            ! No salinity data 
    131109            IF(lwp)WRITE(numout,*)                  ! analytical salinity profile 
    132             IF(lwp)WRITE(numout,*)' Salinity initialisation using a constant value' 
     110            IF(lwp)WRITE(numout,*)'             Salinity initialisation using a constant value' 
    133111            CALL istate_sal 
    134112#endif 
     
    139117      !                                       ! ----------------- 
    140118      CALL wzv( nit000 )                         ! from horizontal divergence 
    141  
     119      ! 
    142120   END SUBROUTINE istate_init 
    143121 
     
    153131      !! 
    154132      !! References :  Philander ??? 
    155       !! 
    156       !! History : 
    157       !!   4.0  !  89-12  (P. Andrich)  Original code 
    158       !!   6.0  !  96-01  (G. Madec)  terrain following coordinates 
    159       !!   9.0  !  02-09  (G. Madec)  F90: Free form 
    160       !!---------------------------------------------------------------------- 
    161       !! * Local declarations 
     133      !!---------------------------------------------------------------------- 
    162134      INTEGER :: ji, jj, jk 
    163135      !!---------------------------------------------------------------------- 
    164  
     136      ! 
    165137      IF(lwp) WRITE(numout,*) 
    166138      IF(lwp) WRITE(numout,*) 'istate_tem : initial temperature profile' 
     
    181153         &                 1     , jpi   , 5     , 1     , jpk   ,   & 
    182154         &                 1     , 1.    , numout                  ) 
    183  
     155      ! 
    184156   END SUBROUTINE istate_tem 
    185157 
     
    194166      !!               
    195167      !! ** Action  :   Initialize sn and sb 
    196       !! 
    197       !! History : 
    198       !!   4.0  !  89-12  (P. Andrich)  Original code 
    199       !!   8.5  !  02-09  (G. Madec)  F90: Free form 
    200       !!---------------------------------------------------------------------- 
    201       !! * Local declarations 
     168      !!---------------------------------------------------------------------- 
    202169      REAL(wp) ::   zsal = 35.50_wp 
    203170      !!---------------------------------------------------------------------- 
     
    224191      !!              - set velocity field including horizontal divergence 
    225192      !!                and relative vorticity fields 
    226       !! 
    227       !! History : 
    228       !!   8.0  !  01-09  (M. Levy, M. Ben Jelloul)  read file for EEL 2 & 6 
    229       !!   9.0  !  03-09  (G. Madec, C. Talandier)  EEL 5 
    230       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    231       !!---------------------------------------------------------------------- 
    232       !! * Modules used 
     193      !!---------------------------------------------------------------------- 
    233194      USE eosbn2     ! eq. of state, Brunt Vaisala frequency (eos     routine) 
    234195      USE divcur     ! hor. divergence & rel. vorticity      (div_cur routine) 
    235196      USE iom 
    236197  
    237       !! * Local declarations 
    238198      INTEGER  ::   inum              ! temporary logical unit 
    239199      INTEGER  ::   ji, jj, jk        ! dummy loop indices 
    240200      INTEGER  ::   ijloc 
    241       REAL(wp) ::   & 
    242          zh1, zh2, zslope, zcst,   &  ! temporary scalars 
    243          zfcor 
    244       REAL(wp) ::   & 
    245          zt1  = 12._wp,            &  ! surface temperature value (EEL R5) 
    246          zt2  =  2._wp,            &  ! bottom  temperature value (EEL R5) 
    247          zsal = 35.5_wp,           &  ! constant salinity (EEL R2, R5 and R6) 
    248          zueel = 0.1_wp               ! constant uniform zonal velocity (EEL R5) 
     201      REAL(wp) ::   zh1, zh2, zslope, zcst, zfcor   ! temporary scalars 
     202      REAL(wp) ::   zt1  = 12._wp,               &  ! surface temperature value (EEL R5) 
     203         &          zt2  =  2._wp,               &  ! bottom  temperature value (EEL R5) 
     204         &          zsal = 35.5_wp,              &  ! constant salinity (EEL R2, R5 and R6) 
     205         &          zueel = 0.1_wp                  ! constant uniform zonal velocity (EEL R5) 
    249206# if ! defined key_dynspg_rl 
    250       REAL(wp), DIMENSION(jpiglo,jpjglo) ::   & 
    251          zssh                         ! initial ssh over the global domain 
     207      REAL(wp), DIMENSION(jpiglo,jpjglo) ::   zssh  ! initial ssh over the global domain 
    252208# endif 
    253209      !!---------------------------------------------------------------------- 
     
    389345      !! ** Method  : - set temprature field 
    390346      !!              - set salinity field 
    391       !! 
    392       !! ** History :                                      
    393       !!      9.0  !  04-05  (A. Koch-Larrouy)  Original code  
    394       !!---------------------------------------------------------------------- 
    395       !! * Modules used 
    396       USE iom 
    397  
    398       !! * Local variables 
    399       INTEGER  ::   inum              ! temporary logical unit 
    400       INTEGER, PARAMETER ::   & 
    401          ntsinit = 0         ! (0/1) (analytical/input data files) T&S initialization 
    402  
     347      !!---------------------------------------------------------------------- 
    403348      INTEGER :: ji, jj, jk  ! dummy loop indices 
     349      INTEGER            ::   inum          ! temporary logical unit 
     350      INTEGER, PARAMETER ::   ntsinit = 0   ! (0/1) (analytical/input data files) T&S initialization 
    404351      !!---------------------------------------------------------------------- 
    405352 
     
    469416      ENDIF 
    470417 
    471       
    472418   END SUBROUTINE istate_gyre 
    473  
    474419 
    475420 
     
    484429      !!      surface to the bottom. 
    485430      !!                 p=integral [ rau*g dz ] 
    486       !! 
    487       !! History : 
    488       !!   8.1  !  01-09  (M. Levy, M. Ben Jelloul)  Original code 
    489       !!   8.5  !  02-09  (G. Madec)  F90: Free form 
    490       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    491       !!---------------------------------------------------------------------- 
    492       !! * Modules used 
     431      !!---------------------------------------------------------------------- 
    493432      USE eosbn2          ! eq. of state, Brunt Vaisala frequency (eos     routine) 
    494433      USE dynspg          ! surface pressure gradient             (dyn_spg routine) 
     
    496435      USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    497436 
    498       !! * Local declarations 
    499437      INTEGER ::   ji, jj, jk        ! dummy loop indices 
    500438      INTEGER ::   indic             ! ??? 
    501       REAL(wp) ::   & 
    502          zmsv, zphv, zmsu, zphu,  &  ! temporary scalars 
    503          zalfg 
    504       REAL(wp), DIMENSION (jpi,jpj,jpk) ::   & 
    505          zprn                        ! workspace 
     439      REAL(wp) ::   zmsv, zphv, zmsu, zphu, zalfg     ! temporary scalars 
     440      REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zprn     ! workspace 
    506441      !!---------------------------------------------------------------------- 
    507442 
     
    514449 
    515450      zalfg = 0.5 * grav * rau0 
    516       ! Surface value 
    517       zprn(:,:,1) = zalfg * fse3w(:,:,1) * ( 1 + rhd(:,:,1) ) 
    518  
    519       ! Vertical integration from the surface 
    520       DO jk = 2, jpkm1 
     451       
     452      zprn(:,:,1) = zalfg * fse3w(:,:,1) * ( 1 + rhd(:,:,1) )       ! Surface value 
     453 
     454      DO jk = 2, jpkm1                                              ! Vertical integration from the surface 
    521455         zprn(:,:,jk) = zprn(:,:,jk-1)   & 
    522456            &         + zalfg * fse3w(:,:,jk) * ( 2. + rhd(:,:,jk) + rhd(:,:,jk-1) ) 
     
    525459      ! Compute geostrophic balance 
    526460      ! --------------------------- 
    527  
    528461      DO jk = 1, jpkm1 
    529462         DO jj = 2, jpjm1 
     
    571504      ! after initializing u and v, we need to calculate the initial streamfunction bsf. 
    572505      ! Otherwise, only the trend will be computed and the model will blow up (inconsistency). 
    573        
    574506      ! to do that, we call dyn_spg with a special trick: 
    575       ! we fill ua and va with the velocities divided by dt, 
    576       ! and the streamfunction will be brought to the right 
    577       ! value assuming the velocities have been set up in 
    578       ! one time step. 
    579       ! we then set bsfd to zero (first guess for next step 
    580       ! is d(psi)/dt = 0.) 
    581  
    582       !  sets up s false trend to calculate the barotropic 
    583       !  streamfunction. 
     507      ! we fill ua and va with the velocities divided by dt, and the streamfunction will be brought to the 
     508      ! right value assuming the velocities have been set up in one time step. 
     509      ! we then set bsfd to zero (first guess for next step is d(psi)/dt = 0.) 
     510      !  sets up s false trend to calculate the barotropic streamfunction. 
    584511 
    585512      ua(:,:,:) = ub(:,:,:) / rdt 
     
    612539      hdivb(:,:,:) = hdivn(:,:,:)       ! set the before to the now value 
    613540      rotb (:,:,:) = rotn (:,:,:)       ! set the before to the now value 
    614  
     541      ! 
    615542   END SUBROUTINE istate_uvg 
    616543 
  • trunk/NEMO/OPA_SRC/opa.F90

    r503 r508  
    4848   USE obc_par         ! open boundary cond. parameters 
    4949   USE obcini          ! open boundary cond. initialization (obc_ini routine) 
    50    USE solver          ! solver initialization          (solver_init routine) 
    5150   USE istate          ! initial state setting          (istate_init routine) 
    5251   USE eosbn2          ! equation of state            (eos bn2 routine) 
     
    6968 
    7069   USE step            ! OPA time-stepping                  (stp     routine) 
    71    USE dynspg_oce      ! Control choice of surface pressure gradient schemes 
    7270   USE prtctl          ! Print control                 (prt_ctl_init routine) 
    7371   USE ini1d           ! re-initialization of u-v mask for the 1D configuration 
     
    245243      CALL istate_init                      ! ocean initial state (Dynamics and tracers) 
    246244 
    247       IF( lk_dynspg_flt .OR. lk_dynspg_rl ) THEN 
    248          CALL solver_init( nit000 )         ! Elliptic solver 
    249       ENDIF 
    250  
    251245!!add 
    252246                       CALL eos( tb, sb, rhd, rhop )        ! before potential and in situ densities 
  • trunk/NEMO/OPA_SRC/par_oce.F90

    r467 r508  
    2626 
    2727   INTEGER, PUBLIC, PARAMETER ::    &  !: 
    28       jpni   = 1,                   &  !: number of processors following i  
    29       jpnj   = 1,                   &  !: number of processors following j 
    30       jpnij  = 1,                   &  !: nb of local domain = nb of processors  
     28      jpni   = 2,                   &  !: number of processors following i  
     29      jpnj   = 2,                   &  !: number of processors following j 
     30      jpnij  = 4,                   &  !: nb of local domain = nb of processors  
    3131      !                                !  ( <= jpni x jpnj ) 
    3232      jpr2di = 0,                   &  !: number of columns for extra outer halo  
  • trunk/NEMO/OPA_SRC/restart.F90

    r473 r508  
    33   !!                     ***  MODULE  restart  *** 
    44   !! Ocean restart :  write the ocean restart file 
    5    !!===================================================================== 
    6  
    7    !!---------------------------------------------------------------------- 
    8    !!   rst_write  : write of the restart file 
    9    !!   rst_read   : read the restart file 
    10    !!---------------------------------------------------------------------- 
    11    !! * Modules used 
     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   !!---------------------------------------------------------------------- 
    1217   USE dom_oce         ! ocean space and time domain 
    1318   USE oce             ! ocean dynamics and tracers  
    1419   USE phycst          ! physical constants 
    15    USE in_out_manager  ! I/O manager 
    1620   USE daymod          ! calendar 
    17    USE sol_oce         ! ocean elliptic solver 
    18    USE zdf_oce         ! ??? 
    19    USE zdftke          ! turbulent kinetic energy scheme 
    2021   USE ice_oce         ! ice variables 
    2122   USE blk_oce         ! bulk variables 
    22    USE flx_oce         ! sea-ice/ocean forcings variables 
    23    USE dynspg_oce      ! free surface time splitting scheme variables 
    2423   USE cpl_oce, ONLY : lk_cpl              ! 
     24   USE in_out_manager  ! I/O manager 
     25   USE iom             ! I/O module 
    2526 
    2627   IMPLICIT NONE 
    2728   PRIVATE 
    2829 
    29    !! * Routine accessibility 
    30    PUBLIC rst_write  ! routine called by step.F90 
    31    PUBLIC rst_read   ! routine called by inidtr.F90 
    32  
    33    !! * Module variables 
    34    CHARACTER (len=48) ::   & 
    35       crestart = 'initial.nc'   ! restart file name 
    36    !!---------------------------------------------------------------------- 
    37    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     30   PUBLIC   rst_opn    ! routine called by step module 
     31   PUBLIC   rst_write  ! routine called by step module 
     32   PUBLIC   rst_read   ! routine called by opa  module 
     33 
     34   LOGICAL, PUBLIC ::   lrst_oce         !: logical to control the oce restart write  
     35   INTEGER, PUBLIC ::   nitrst           !: time step at which restart file should be written 
     36   INTEGER, PUBLIC ::   numror, numrow   !: logical unit for cean restart (read and write) 
     37 
     38   !! * Substitutions 
     39#  include "vectopt_loop_substitute.h90" 
     40   !!---------------------------------------------------------------------- 
     41   !!  OPA 9.0 , LOCEAN-IPSL (2006)  
    3842   !! $Header$  
    39    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    40    !!---------------------------------------------------------------------- 
    41  
     43   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     44   !!---------------------------------------------------------------------- 
    4245 
    4346CONTAINS 
     47 
     48   SUBROUTINE rst_opn( kt ) 
     49      !!--------------------------------------------------------------------- 
     50      !!                   ***  ROUTINE rst_opn  *** 
     51      !!                      
     52      !! ** Purpose : + initialization (should be read in the namelist) of nitrst  
     53      !!              + open the restart when we are one time step before nitrst 
     54      !!                   - restart header is defined when kt = nitrst-1 
     55      !!                   - restart data  are written when kt = nitrst 
     56      !!              + define lrst_oce to .TRUE. when we need to define or write the restart 
     57      !!---------------------------------------------------------------------- 
     58      INTEGER, INTENT(in) ::   kt     ! ocean time-step 
     59      !! 
     60      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step deine as a character 
     61      CHARACTER(LEN=50)   ::   clname   ! ice output restart file name 
     62      !!---------------------------------------------------------------------- 
     63      ! 
     64      IF( kt == nit000 ) THEN   ! default initialization, to do: should be read in the namelist... 
     65         nitrst = nitend        ! to do: should be read in the namelist in a cleaver way... 
     66         lrst_oce = .FALSE. 
     67      ENDIF 
     68       
     69      IF    ( kt == nitrst-1 .AND. lrst_oce         ) THEN 
     70         CALL ctl_stop( 'rst_opn: we cannot create an ocean restart at every time step' ) 
     71         numrow = 0 
     72      ELSEIF( kt == nitrst-1 .OR.  nitend == nit000 ) THEN   ! beware if model runs only one time step 
     73         ! beware of the format used to write kt (default is i8.8, that should be large enough) 
     74         IF( nitrst > 1.0e9 ) THEN    
     75            WRITE(clkt,*) nitrst 
     76         ELSE 
     77            WRITE(clkt,'(i8.8)') nitrst 
     78         ENDIF 
     79         ! create the file 
     80         IF(lwp) WRITE(numout,*) 
     81         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart" 
     82         IF(lwp) WRITE(numout,*) '             open ocean restart.output NetCDF file: '//clname 
     83         CALL iom_open( clname, numrow, ldwrt = .TRUE. ) 
     84         lrst_oce = .TRUE. 
     85      ENDIF 
     86      ! 
     87   END SUBROUTINE rst_opn 
     88 
    4489 
    4590#if  ( defined key_mpp_mpi   ||   defined key_mpp_shmem ) && defined key_dimgout 
     
    66111      !! ** Purpose :   Write restart fields in NetCDF format 
    67112      !! 
    68       !! ** Method  :   Write in numwrs file each nstock time step in NetCDF 
     113      !! ** Method  :   Write in numrow when kt == nitrst in NetCDF 
    69114      !!      file, save fields which are necessary for restart 
    70       !! 
    71       !! History : 
    72       !!        !  99-11  (M. Imbard)  Original code 
    73       !!   8.5  !  02-08  (G. Madec)  F90: Free form 
    74       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    75       !!---------------------------------------------------------------------- 
    76       !! * Modules used 
    77       USE ioipsl 
    78  
    79       !! * Arguments  
    80       INTEGER, INTENT( in ) ::   kt         ! ocean time-step 
    81  
    82       !! * Local declarations 
    83       LOGICAL ::   llbon 
    84       CHARACTER (len=50) ::   clname, cln 
    85       INTEGER ::   ic, jc, itime 
    86       INTEGER ::   inumwrs 
    87       REAL(wp) ::   zdate0 
    88       REAL(wp), DIMENSION( 1) ::   zfice, zfblk   ! used only in case of ice & bulk 
    89       REAL(wp), DIMENSION(10) ::   zinfo(10) 
    90       REAL(wp), DIMENSION(jpi,jpj) :: ztab  
    91 #if defined key_agrif 
    92        Integer :: knum 
    93 #endif 
    94       !!---------------------------------------------------------------------- 
    95  
    96       IF( kt == nit000 ) THEN 
    97          IF(lwp) WRITE(numout,*) 
    98          IF(lwp) WRITE(numout,*) 'rst_wri : write restart.output NetCDF file' 
    99          IF(lwp) WRITE(numout,*) '~~~~~~~' 
    100          zfice(1) = 1.e0   ;   zfblk(1) = 1.e0 
    101       ENDIF 
    102  
    103  
    104       IF( MOD( kt, nstock ) == 0 .OR. kt == nitend ) THEN 
    105           
    106          ! 0. Initializations 
    107          ! ------------------ 
    108  
    109          IF(lwp) WRITE(numout,*) ' ' 
    110          IF(lwp) WRITE(numout,*) 'rst_write : write the restart file in NetCDF format ',   & 
    111                                               'at it= ',kt,' date= ',ndastp 
    112          IF(lwp) WRITE(numout,*) '~~~~~~~~~' 
    113  
    114          ! Job informations  
    115          zinfo(:) = 0.e0  
    116          zinfo(1) = FLOAT( no        )   ! job number 
    117          zinfo(2) = FLOAT( kt        )   ! time-step 
    118          zinfo(3) = FLOAT( 2 - nsolv )   ! pcg solver 
    119          zinfo(4) = FLOAT( nsolv - 1 )   ! sor solver 
    120          IF( lk_zdftke ) THEN 
    121             zinfo(5) = 1.e0              ! TKE  
    122          ELSE 
    123             zinfo(5) = 0.e0              ! no TKE  
    124          ENDIF 
    125          zinfo(6) = FLOAT( ndastp )      ! date 
    126          zinfo(7) = adatrj               ! ??? 
    127  
    128          ! delete the restart file if it exists  
    129          INQUIRE( FILE=crestart, EXIST=llbon ) 
    130          IF(llbon) THEN 
    131 #if defined key_agrif 
    132        knum =Agrif_Get_Unit() 
    133             OPEN( UNIT=knum, FILE=crestart, STATUS='old' ) 
    134             CLOSE( knum, STATUS='delete' ) 
    135 #else             
    136             OPEN( UNIT=inumwrs, FILE=crestart, STATUS='old' ) 
    137             CLOSE( inumwrs, STATUS='delete' ) 
    138 #endif 
    139          ENDIF 
    140  
    141          ! Name of the new restart file 
    142          ic     = 1 
    143          DO jc = 1, 16 
    144             IF( cexper(jc:jc) /= ' ' )   ic = jc 
    145          END DO 
    146          WRITE(cln,'("_",i4.4,i2.2,i2.2,"_restart")') nyear, nmonth, nday 
    147          clname = cexper(1:ic)//cln 
    148          ic = 1 
    149          DO jc = 1, 48 
    150             IF( clname(jc:jc) /= ' ' ) ic = jc 
    151          END DO 
    152          crestart = clname(1:ic)//".nc" 
    153          itime = 0 
    154          CALL ymds2ju( nyear, nmonth, nday, 0.e0, zdate0 ) 
    155          CALL restini( 'NONE', jpi, jpj, glamt, gphit, jpk, gdept_0, clname,   & 
    156                         itime, zdate0, rdt*nstock ,inumwrs, domain_id=nidom ) 
    157  
    158          CALL restput( inumwrs, 'info'   , 1  , 1  , 10 , 0, zinfo   )   ! restart informations 
    159           
    160          CALL restput( inumwrs, 'ub'     , jpi, jpj, jpk, 0, ub      )   ! prognostic variables 
    161          CALL restput( inumwrs, 'vb'     , jpi, jpj, jpk, 0, vb      ) 
    162          CALL restput( inumwrs, 'tb'     , jpi, jpj, jpk, 0, tb      ) 
    163          CALL restput( inumwrs, 'sb'     , jpi, jpj, jpk, 0, sb      ) 
    164          CALL restput( inumwrs, 'rotb'   , jpi, jpj, jpk, 0, rotb    ) 
    165          CALL restput( inumwrs, 'hdivb'  , jpi, jpj, jpk, 0, hdivb   ) 
    166          CALL restput( inumwrs, 'un'     , jpi, jpj, jpk, 0, un      ) 
    167          CALL restput( inumwrs, 'vn'     , jpi, jpj, jpk, 0, vn      ) 
    168          CALL restput( inumwrs, 'tn'     , jpi, jpj, jpk, 0, tn      ) 
    169          CALL restput( inumwrs, 'sn'     , jpi, jpj, jpk, 0, sn      ) 
    170          CALL restput( inumwrs, 'rotn'   , jpi, jpj, jpk, 0, rotn    ) 
    171          CALL restput( inumwrs, 'hdivn'  , jpi, jpj, jpk, 0, hdivn   ) 
    172  
    173          ztab(:,:) = gcx(1:jpi,1:jpj) 
    174          CALL restput( inumwrs, 'gcx'    , jpi, jpj, 1  , 0, ztab    )   ! Read elliptic solver arrays 
    175          ztab(:,:) = gcxb(1:jpi,1:jpj) 
    176          CALL restput( inumwrs, 'gcxb'   , jpi, jpj, 1  , 0, ztab    ) 
    177 # if defined key_dynspg_rl 
    178          CALL restput( inumwrs, 'bsfb'   , jpi, jpj, 1  , 0, bsfb    )   ! Rigid-lid formulation (bsf) 
    179          CALL restput( inumwrs, 'bsfn'   , jpi, jpj, 1  , 0, bsfn    ) 
    180          CALL restput( inumwrs, 'bsfd'   , jpi, jpj, 1  , 0, bsfd    ) 
    181 # else 
    182          CALL restput( inumwrs, 'sshb'   , jpi, jpj, 1  , 0, sshb    )   ! free surface formulation (ssh) 
    183          CALL restput( inumwrs, 'sshn'   , jpi, jpj, 1  , 0, sshn    ) 
    184 #  if defined key_dynspg_ts 
    185          CALL restput( inumwrs, 'sshb_b' , jpi, jpj, 1  , 0, sshb_b  )   ! free surface formulation (ssh) 
    186          CALL restput( inumwrs, 'sshn_b' , jpi, jpj, 1  , 0, sshn_b  )   ! issued from barotropic loop 
    187          CALL restput( inumwrs, 'un_b'   , jpi, jpj, 1  , 0, un_b    )   ! horizontal transports 
    188          CALL restput( inumwrs, 'vn_b'   , jpi, jpj, 1  , 0, vn_b    )   ! issued from barotropic loop 
     115      !!---------------------------------------------------------------------- 
     116      INTEGER, INTENT(in) ::   kt   ! ocean time-step 
     117      !!---------------------------------------------------------------------- 
     118 
     119      IF(lwp) THEN 
     120         WRITE(numout,*) 
     121         WRITE(numout,*) 'rst_write : write ocean NetCDF restart file  kt =', kt,' date= ', ndastp 
     122         WRITE(numout,*) '~~~~~~~~~' 
     123      ENDIF 
     124       
     125      ! calendar control 
     126      CALL iom_rstput( kt, nitrst, numrow, 'kt'     , REAL( kt    , wp) )   ! time-step  
     127      CALL iom_rstput( kt, nitrst, numrow, 'ndastp' , REAL( ndastp, wp) )   ! date 
     128      CALL iom_rstput( kt, nitrst, numrow, 'adatrj' ,       adatrj      )   ! number of elapsed days since 
     129      !                                                                     ! the begining of the run [s] 
     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      CALL iom_rstput( kt, nitrst, numrow, 'tn'     , tn      ) 
     141      CALL iom_rstput( kt, nitrst, numrow, 'sn'     , sn      ) 
     142      CALL iom_rstput( kt, nitrst, numrow, 'rotn'   , rotn    ) 
     143      CALL iom_rstput( kt, nitrst, numrow, 'hdivn'  , hdivn   ) 
     144 
     145# if defined key_ice_lim         
     146      CALL iom_rstput( kt, nitrst, numrow, 'nfice'  , REAL( nfice, wp) )   !  ice computation frequency 
     147      CALL iom_rstput( kt, nitrst, numrow, 'sst_io' , sst_io  ) 
     148      CALL iom_rstput( kt, nitrst, numrow, 'sss_io' , sss_io  ) 
     149      CALL iom_rstput( kt, nitrst, numrow, 'u_io'   , u_io    ) 
     150      CALL iom_rstput( kt, nitrst, numrow, 'v_io'   , v_io    ) 
     151#  if defined key_coupled 
     152      CALL iom_rstput( kt, nitrst, numrow, 'alb_ice', alb_ice ) 
    189153#  endif 
    190154# endif 
    191 # if defined key_zdftke   ||   defined key_esopa 
    192          IF( lk_zdftke ) THEN 
    193             CALL restput( inumwrs, 'en'     , jpi, jpj, jpk, 0, en      )   ! TKE arrays 
    194          ENDIF 
    195 # endif 
    196 # if defined key_ice_lim 
    197          zfice(1) = FLOAT( nfice )                                      ! Louvain La Neuve Sea Ice Model 
    198          CALL restput( inumwrs, 'nfice'  ,   1,   1, 1  , 0, zfice   ) 
    199          CALL restput( inumwrs, 'sst_io' , jpi, jpj, 1  , 0, sst_io  ) 
    200          CALL restput( inumwrs, 'sss_io' , jpi, jpj, 1  , 0, sss_io  ) 
    201          CALL restput( inumwrs, 'u_io'   , jpi, jpj, 1  , 0, u_io    ) 
    202          CALL restput( inumwrs, 'v_io'   , jpi, jpj, 1  , 0, v_io    ) 
    203 # if defined key_coupled 
    204          CALL restput( inumwrs, 'alb_ice', jpi, jpj, 1  , 0, alb_ice ) 
    205 # endif 
    206 # endif 
    207155# if defined key_flx_bulk_monthly || defined key_flx_bulk_daily 
    208          zfblk(1) = FLOAT( nfbulk )                                 ! Bulk 
    209          CALL restput( inumwrs, 'nfbulk' ,   1,   1, 1  , 0, zfblk   ) 
    210          CALL restput( inumwrs, 'gsst'   , jpi, jpj, 1  , 0, gsst    ) 
    211 # endif 
    212  
    213          CALL restclo( inumwrs )                                         ! close the restart file 
    214           
    215       ENDIF 
    216  
     156      CALL iom_rstput( kt, nitrst, numrow, 'nfbulk' , REAL( nfbulk, wp) )   !  bulk computation frequency 
     157      CALL iom_rstput( kt, nitrst, numrow, 'gsst'   , gsst    ) 
     158# endif 
     159 
     160      IF( kt == nitrst ) THEN 
     161         CALL iom_close( numrow )     ! close the restart file (only at last time step) 
     162         lrst_oce = .FALSE. 
     163      ENDIF 
     164      ! 
    217165   END SUBROUTINE rst_write 
    218166 
     
    246194      !!       nrstdt = 2  the duration of the experiment in days (adatrj) 
    247195      !!                    has been stored in the restart file. 
    248       !! 
    249       !! History : 
    250       !!        !  99-05  (M. Imbard)  Original code 
    251       !!   8.5  !  02-09  (G. Madec)  F90: Free form 
    252       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    253       !!---------------------------------------------------------------------- 
    254       !! * Modules used 
    255       USE iom 
    256  
    257       !! * Local declarations 
    258       INTEGER  ::   & 
    259          inum                 ! temporary logical unit 
    260       REAL(wp), DIMENSION(1, 1, 10)  ::   zinfo 
    261       REAL(wp), DIMENSION(1, 1, 1)   ::   zzz  
    262       INTEGER  ::   ios 
    263 #   if defined key_ice_lim 
     196      !!---------------------------------------------------------------------- 
     197      REAL(wp) ::   zcoef, zkt, zndastp, znfice, znfbulk 
     198# if defined key_ice_lim 
    264199      INTEGER  ::   ji, jj 
    265 #   endif 
    266       !!---------------------------------------------------------------------- 
    267  
    268       IF(lwp) WRITE(numout,*) 
    269       IF(lwp) WRITE(numout,*) 'rst_read : read the NetCDF restart file' 
    270       IF(lwp) WRITE(numout,*) '~~~~~~~~' 
    271  
    272       IF(lwp) WRITE(numout,*) ' Info on the present job : ' 
    273       IF(lwp) WRITE(numout,*) '   job number          : ', no 
    274       IF(lwp) WRITE(numout,*) '   time-step           : ', nit000 
    275       IF(lwp) WRITE(numout,*) '   solver type         : ', nsolv 
    276       IF( lk_zdftke ) THEN 
    277          IF(lwp) WRITE(numout,*) '   tke option          : 1 ' 
    278       ELSE 
    279          IF(lwp) WRITE(numout,*) '   tke option          : 0 ' 
    280       ENDIF 
    281       IF(lwp) WRITE(numout,*) '   date ndastp         : ', ndastp 
    282       IF(lwp) WRITE(numout,*) 
    283  
    284       ! Time domain : restart 
    285       ! ------------------------- 
    286  
    287       IF(lwp) WRITE(numout,*) 
    288       IF(lwp) WRITE(numout,*) 
    289       IF(lwp) WRITE(numout,*) ' *** restart option' 
    290       SELECT CASE ( nrstdt ) 
    291       CASE ( 0 )  
    292          IF(lwp) WRITE(numout,*) ' nrstdt = 0 no control of nit000' 
    293       CASE ( 1 )  
    294          IF(lwp) WRITE(numout,*) ' nrstdt = 1 we control the date of nit000' 
    295       CASE ( 2 ) 
    296          IF(lwp) WRITE(numout,*) ' nrstdt = 2 the date adatrj is read in restart file' 
    297       CASE DEFAULT 
    298          IF(lwp) WRITE(numout,*) '  ===>>>> nrstdt not equal 0, 1 or 2 : no control of the date' 
    299          IF(lwp) WRITE(numout,*) ' =======                   =========' 
    300       END SELECT 
    301  
    302       CALL iom_open ( 'restart', inum ) 
    303        
    304       CALL iom_get ( inum, jpdom_unknown, 'info', zinfo ) 
    305        
    306       IF(lwp) WRITE(numout,*) 
    307       IF(lwp) WRITE(numout,*) ' Info on the restart file read : ' 
    308       IF(lwp) WRITE(numout,*) '   job number          : ', NINT( zinfo(1, 1, 1) ) 
    309       IF(lwp) WRITE(numout,*) '   time-step           : ', NINT( zinfo(1, 1, 2) ) 
    310       IF(lwp) WRITE(numout,*) '   solver type         : ', NINT( zinfo(1, 1, 4) ) + 1 
    311       IF(lwp) WRITE(numout,*) '   tke option          : ', NINT( zinfo(1, 1, 5) ) 
    312       IF(lwp) WRITE(numout,*) '   date ndastp         : ', NINT( zinfo(1, 1, 6) ) 
    313       IF(lwp) WRITE(numout,*) 
    314  
     200# endif 
     201      !!---------------------------------------------------------------------- 
     202 
     203      IF(lwp) THEN                                             ! Contol prints 
     204         WRITE(numout,*) 
     205         WRITE(numout,*) 'rst_read : read oce NetCDF restart file' 
     206         WRITE(numout,*) '~~~~~~~~' 
     207          
     208         WRITE(numout,*) ' *** Info on the present job : ' 
     209         WRITE(numout,*) '   time-step           : ', nit000 
     210!!$         WRITE(numout,*) '   solver type         : ', nsolv 
     211!!$         IF( lk_zdftke ) THEN 
     212!!$            WRITE(numout,*) '   tke option          : 1 ' 
     213!!$         ELSE 
     214!!$            WRITE(numout,*) '   tke option          : 0 ' 
     215!!$         ENDIF 
     216         WRITE(numout,*) '   date ndastp         : ', ndastp 
     217         WRITE(numout,*) 
     218         WRITE(numout,*) ' *** restart option' 
     219         SELECT CASE ( nrstdt ) 
     220         CASE ( 0 )  
     221            WRITE(numout,*) ' nrstdt = 0 no control of nit000' 
     222         CASE ( 1 )  
     223            WRITE(numout,*) ' nrstdt = 1 we control the date of nit000' 
     224         CASE ( 2 ) 
     225            WRITE(numout,*) ' nrstdt = 2 the date adatrj is read in restart file' 
     226         CASE DEFAULT 
     227            WRITE(numout,*) '  ===>>>> nrstdt not equal 0, 1 or 2 : no control of the date' 
     228            WRITE(numout,*) '  =======                  =========' 
     229         END SELECT 
     230         WRITE(numout,*) 
     231      ENDIF 
     232 
     233      CALL iom_open( 'restart', numror )                       ! Open 
     234 
     235      ! Calendar informations 
     236      CALL iom_get( numror, 'kt'    , zkt     )   ! time-step  
     237      CALL iom_get( numror, 'ndastp', zndastp )   ! date 
     238      ! Additional contol prints 
     239      IF(lwp) THEN 
     240         WRITE(numout,*) 
     241         WRITE(numout,*) ' *** Info on the restart file read : ' 
     242         WRITE(numout,*) '   time-step           : ', NINT( zkt ) 
     243!!$         WRITE(numout,*) '   solver type         : ', +++ 
     244!!$         WRITE(numout,*) '   tke option          : ', +++ 
     245         WRITE(numout,*) '   date ndastp         : ', NINT( zndastp ) 
     246         WRITE(numout,*) 
     247      ENDIF 
    315248      ! Control of date 
    316       IF( nit000 - NINT( zinfo(1, 1, 2) )  /= 1 .AND. nrstdt /= 0 ) & 
     249      IF( nit000 - NINT( zkt )  /= 1 .AND. nrstdt /= 0 ) & 
    317250           & CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart', & 
    318251           & ' verify the restart file or rerun with nrstdt = 0 (namelist)' ) 
    319  
    320252      ! re-initialisation of  adatrj0 
    321       adatrj0 =  ( FLOAT( nit000-1 ) * rdttra(1) ) / rday 
    322  
     253      adatrj0 = ( REAL( nit000-1, wp ) * rdttra(1) ) / rday 
    323254      IF ( nrstdt == 2 ) THEN 
    324255!                             by default ndatsp has been set to ndate0 in dom_nam 
    325256!                             ndate0 has been read in the namelist (standard OPA 8) 
    326257!                             here when nrstdt=2 we keep the  final date of previous run 
    327         ndastp = NINT( zinfo(1, 1, 6) ) 
    328         adatrj0 =  zinfo(1, 1, 7) 
    329       ENDIF 
    330  
    331       CALL iom_get( inum, jpdom_local, 'ub'   , ub    )   ! Read prognostic variables 
    332       CALL iom_get( inum, jpdom_local, 'vb'   , vb    ) 
    333       CALL iom_get( inum, jpdom_local, 'tb'   , tb    ) 
    334       CALL iom_get( inum, jpdom_local, 'sb'   , sb    ) 
    335       CALL iom_get( inum, jpdom_local, 'rotb' , rotb  ) 
    336       CALL iom_get( inum, jpdom_local, 'hdivb', hdivb ) 
    337       CALL iom_get( inum, jpdom_local, 'un'   , un    ) 
    338       CALL iom_get( inum, jpdom_local, 'vn'   , vn    ) 
    339       CALL iom_get( inum, jpdom_local, 'tn'   , tn    ) 
    340       CALL iom_get( inum, jpdom_local, 'sn'   , sn    ) 
    341       CALL iom_get( inum, jpdom_local, 'rotn' , rotn  ) 
    342       CALL iom_get( inum, jpdom_local, 'hdivn', hdivn ) 
    343 ! Caution : extrahallow  
    344 ! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
    345       CALL iom_get( inum, jpdom_local, 'gcx' , gcx (1:jpi,1:jpj) ) 
    346       CALL iom_get( inum, jpdom_local, 'gcxb', gcxb(1:jpi,1:jpj) )     ! Read elliptic solver arrays 
    347 # if defined key_dynspg_rl 
    348       CALL iom_get( inum, jpdom_local, 'bsfb', bsfb )     ! Rigid-lid formulation (bsf) 
    349       CALL iom_get( inum, jpdom_local, 'bsfn', bsfn ) 
    350       CALL iom_get( inum, jpdom_local, 'bsfd', bsfd ) 
    351 # else 
    352       CALL iom_get( inum, jpdom_local, 'sshb', sshb )     ! free surface formulation (ssh) 
    353       CALL iom_get( inum, jpdom_local, 'sshn', sshn ) 
    354 #  if defined key_dynspg_ts 
    355       CALL iom_get( inum, jpdom_local, 'sshb_b', sshb_b ) ! free surface formulation (ssh) 
    356       CALL iom_get( inum, jpdom_local, 'sshn_b', sshn_b ) ! issued from barotropic loop 
    357       CALL iom_get( inum, jpdom_local, 'un_b'  , un_b )   ! horizontal transports 
    358       CALL iom_get( inum, jpdom_local, 'vn_b'  , vn_b )   ! issued from barotropic loop 
    359 #  endif 
    360 # endif 
    361 # if defined key_zdftke   ||   defined key_esopa 
    362       IF( lk_zdftke ) THEN 
    363          IF( NINT( zinfo(1, 1, 5) ) == 1 ) THEN                                ! Read tke arrays 
    364             CALL iom_get( inum, jpdom_local, 'en', en ) 
    365             ln_rstke = .FALSE. 
    366          ELSE 
    367             IF(lwp) WRITE(numout,*) ' ===>>>> : the previous restart file did not used  tke scheme' 
    368             IF(lwp) WRITE(numout,*) ' =======                =======' 
    369             nrstdt = 2 
    370             ln_rstke = .TRUE. 
    371          ENDIF 
    372       ENDIF 
    373 # endif 
     258         ndastp = NINT( zndastp ) 
     259        CALL iom_get( numror, 'adatrj', adatrj )   ! number of elapsed days since the begining of last run 
     260      ENDIF 
     261 
     262      !                                                       ! Read prognostic variables 
     263      CALL iom_get( numror, jpdom_local, 'ub'   , ub    )        ! before i-component velocity 
     264      CALL iom_get( numror, jpdom_local, 'vb'   , vb    )        ! before j-component velocity 
     265      CALL iom_get( numror, jpdom_local, 'tb'   , tb    )        ! before temperature 
     266      CALL iom_get( numror, jpdom_local, 'sb'   , sb    )        ! before salinity 
     267      CALL iom_get( numror, jpdom_local, 'rotb' , rotb  )        ! before curl 
     268      CALL iom_get( numror, jpdom_local, 'hdivb', hdivb )        ! before horizontal divergence 
     269      CALL iom_get( numror, jpdom_local, 'un'   , un    )        ! now    i-component velocity 
     270      CALL iom_get( numror, jpdom_local, 'vn'   , vn    )        ! now    j-component velocity 
     271      CALL iom_get( numror, jpdom_local, 'tn'   , tn    )        ! now    temperature 
     272      CALL iom_get( numror, jpdom_local, 'sn'   , sn    )        ! now    salinity 
     273      CALL iom_get( numror, jpdom_local, 'rotn' , rotn  )        ! now    curl 
     274      CALL iom_get( numror, jpdom_local, 'hdivn', hdivn )        ! now    horizontal divergence 
     275 
     276 
     277      IF( neuler == 0 ) THEN                                  ! Euler restart (neuler=0) 
     278         tb   (:,:,:) = tn   (:,:,:)                             ! all before fields set to now field values 
     279         sb   (:,:,:) = sn   (:,:,:) 
     280         ub   (:,:,:) = un   (:,:,:) 
     281         vb   (:,:,:) = vn   (:,:,:) 
     282         rotb (:,:,:) = rotn (:,:,:) 
     283         hdivb(:,:,:) = hdivn(:,:,:) 
     284      ENDIF 
     285 
     286      !!sm: TO BE MOVED IN NEW SURFACE MODULE... 
     287 
    374288# if defined key_ice_lim 
    375289      ! Louvain La Neuve Sea Ice Model 
    376       ios = iom_varid( inum, 'nfice' ) 
    377       IF( ios > 0 ) then  
    378          CALL iom_get( inum, jpdom_unknown, 'nfice' , zzz ) 
    379          zinfo(1, 1, 8) = zzz(1, 1, 1) 
    380          CALL iom_get( inum, jpdom_local, 'sst_io', sst_io ) 
    381          CALL iom_get( inum, jpdom_local, 'sss_io', sss_io ) 
    382          CALL iom_get( inum, jpdom_local, 'u_io'  , u_io ) 
    383          CALL iom_get( inum, jpdom_local, 'v_io'  , v_io ) 
     290      IF( iom_varid( numror, 'nfice' ) > 0 ) then  
     291         CALL iom_get( numror             , 'nfice'  , znfice  )   ! ice computation frequency 
     292         CALL iom_get( numror, jpdom_local, 'sst_io' , sst_io  ) 
     293         CALL iom_get( numror, jpdom_local, 'sss_io' , sss_io  ) 
     294         CALL iom_get( numror, jpdom_local, 'u_io'   , u_io    ) 
     295         CALL iom_get( numror, jpdom_local, 'v_io'   , v_io    ) 
    384296#if defined key_coupled 
    385          CALL iom_get( inum, jpdom_local, 'alb_ice', alb_ice ) 
     297         CALL iom_get( numror, jpdom_local, 'alb_ice', alb_ice ) 
    386298#endif 
    387       ENDIF 
    388       IF( zinfo(1, 1, 8) /= FLOAT(nfice) .OR. ios == 0 ) THEN 
     299         IF( znfice /= REAL( nfice, wp ) ) THEN      ! if nfice changed between 2 runs 
     300            zcoef = REAL( nfice-1, wp ) / znfice 
     301            sst_io(:,:) = zcoef * sst_io(:,:) 
     302            sss_io(:,:) = zcoef * sss_io(:,:) 
     303            u_io  (:,:) = zcoef * u_io  (:,:) 
     304            v_io  (:,:) = zcoef * v_io  (:,:) 
     305         ENDIF 
     306      ELSE 
    389307         IF(lwp) WRITE(numout,*) 
    390308         IF(lwp) WRITE(numout,*) 'rst_read :  LLN sea Ice Model => Ice initialization' 
    391309         IF(lwp) WRITE(numout,*) 
    392          sst_io(:,:) = ( nfice-1 )*( tn(:,:,1) + rt0 )          !!bug a explanation is needed here! 
    393          sss_io(:,:) = ( nfice-1 )*  sn(:,:,1) 
     310         zcoef = REAL( nfice-1, wp ) 
     311         sst_io(:,:) = zcoef *( tn(:,:,1) + rt0 )          !!bug a explanation is needed here! 
     312         sss_io(:,:) = zcoef *  sn(:,:,1) 
     313         zcoef = 0.5 * REAL( nfice-1, wp ) 
    394314         DO jj = 2, jpj 
    395             DO ji = 2, jpi 
    396                u_io(ji,jj) = ( nfice-1 ) * 0.5 * ( un(ji-1,jj  ,1) + un(ji-1,jj-1,1) ) 
    397                v_io(ji,jj) = ( nfice-1 ) * 0.5 * ( vn(ji  ,jj-1,1) + vn(ji-1,jj-1,1) ) 
     315            DO ji = fs_2, jpi   ! vector opt. 
     316               u_io(ji,jj) = zcoef * ( un(ji-1,jj  ,1) + un(ji-1,jj-1,1) ) 
     317               v_io(ji,jj) = zcoef * ( vn(ji  ,jj-1,1) + vn(ji-1,jj-1,1) ) 
    398318            END DO 
    399319         END DO 
     
    405325# if defined key_flx_bulk_monthly || defined key_flx_bulk_daily 
    406326      ! Louvain La Neuve Sea Ice Model 
    407       ios = iom_varid( inum, 'nfbulk' ) 
    408       IF( ios > 0 ) then  
    409          CALL iom_get( inum, jpdom_unknown, 'nfbulk' , zzz ) 
    410          CALL iom_get( inum, jpdom_local, 'gsst' , gsst ) 
    411          zinfo(1, 1, 9) = zzz(1, 1, 1) 
    412       ENDIF 
    413       IF( zinfo(1, 1, 9) /= FLOAT(nfbulk) .OR. ios == 0 ) THEN 
     327      IF( iom_varid( numror, 'nfbulk' ) > 0 ) THEN  
     328         CALL iom_get( numror             , 'nfbulk', znfbulk )   ! bulk computation frequency 
     329         CALL iom_get( numror, jpdom_local, 'gsst'  , gsst    ) 
     330         IF( znfbulk /= REAL(nfbulk, wp) ) THEN      ! if you change nfbulk between 2 runs 
     331            zcoef = REAL( nfbulk-1, wp ) / znfbulk 
     332            gsst(:,:) = zcoef * gsst(:,:) 
     333         ENDIF 
     334      ELSE 
    414335         IF(lwp) WRITE(numout,*) 
    415336         IF(lwp) WRITE(numout,*) 'rst_read :  LLN sea Ice Model => Ice initialization' 
    416337         IF(lwp) WRITE(numout,*) 
    417          gsst(:,:) = 0. 
    418          gsst(:,:) = gsst(:,:) + ( nfbulk-1 )*( tn(:,:,1) + rt0 ) 
     338         gsst(:,:) = REAL( nfbulk - 1, wp )*( tn(:,:,1) + rt0 ) 
    419339      ENDIF 
    420340# endif 
    421341       
    422       CALL iom_close( inum ) 
    423  
    424   ! In case of restart with neuler = 0 then put all before fields = to now fields 
    425     IF ( neuler == 0 ) THEN 
    426       tb(:,:,:)=tn(:,:,:) 
    427       sb(:,:,:)=sn(:,:,:) 
    428       ub(:,:,:)=un(:,:,:) 
    429       vb(:,:,:)=vn(:,:,:) 
    430       rotb(:,:,:)=rotn(:,:,:) 
    431       hdivb(:,:,:)=hdivn(:,:,:) 
    432 #if defined key_dynspg_rl 
    433     ! rigid lid 
    434       bsfb(:,:)=bsfn(:,:) 
    435 #else 
    436     ! free surface formulation (eta) 
    437       sshb(:,:)=sshn(:,:) 
     342      !!sm: end of TO BE MOVED IN NEW SURFACE MODULE... 
     343      ! 
     344   END SUBROUTINE rst_read 
     345 
    438346#endif 
    439     ENDIF 
    440  
    441    END SUBROUTINE rst_read 
    442  
    443 #endif 
     347 
    444348   !!===================================================================== 
    445349END MODULE restart 
  • trunk/NEMO/OPA_SRC/step.F90

    r503 r508  
    44   !! Time-stepping    : manager of the ocean, tracer and ice time stepping 
    55   !!====================================================================== 
    6    !! History : 
    7    !!        !  91-03  ()  Original code 
    8    !!        !  91-11  (G. Madec) 
    9    !!        !  92-06  (M. Imbard)  add a first output record 
    10    !!        !  96-04  (G. Madec)  introduction of dynspg 
    11    !!        !  96-04  (M.A. Foujols)  introduction of passive tracer 
    12    !!   8.0  !  97-06  (G. Madec)  new architecture of call 
    13    !!   8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface 
    14    !!   8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit 
    15    !!   8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions 
    16    !!   9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking 
    17    !!    "   !  04-08  (C. Talandier) New trends organization 
    18    !!    "   !  05-01  (C. Ethe) Add the KPP closure scheme 
    19    !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization 
    20    !!    "   !  05-11  (G. Madec)  Reorganisation of tra and dyn calls 
     6   !! History :        !  91-03  ()  Original code 
     7   !!                  !  91-11  (G. Madec) 
     8   !!                  !  92-06  (M. Imbard)  add a first output record 
     9   !!                  !  96-04  (G. Madec)  introduction of dynspg 
     10   !!                  !  96-04  (M.A. Foujols)  introduction of passive tracer 
     11   !!             8.0  !  97-06  (G. Madec)  new architecture of call 
     12   !!             8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface 
     13   !!             8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit 
     14   !!             8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions 
     15   !!             9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking 
     16   !!             " "  !  04-08  (C. Talandier) New trends organization 
     17   !!             " "  !  05-01  (C. Ethe) Add the KPP closure scheme 
     18   !!             " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
     19   !!             " "  !  05-11  (G. Madec)  Reorganisation of tra and dyn calls 
     20   !!             " "  !  06-01  (L. Debreu, C. Mazauric)  Agrif implementation 
     21   !!             " "  !  06-07  (S. Masson)  restart using iom 
    2122   !!---------------------------------------------------------------------- 
    2223   !!   stp            : OPA system time-stepping 
    2324   !!---------------------------------------------------------------------- 
    24    !! * Modules used 
    2525   USE oce             ! ocean dynamics and tracers variables 
    2626   USE dom_oce         ! ocean space and time domain variables  
     
    3030   USE cpl_oce         ! coupled ocean-atmosphere variables 
    3131   USE in_out_manager  ! I/O manager 
     32   USE iom 
    3233   USE lbclnk 
    3334 
     
    126127   PRIVATE 
    127128 
    128    !! * Routine accessibility 
    129129   PUBLIC stp            ! called by opa.F90 
    130130 
     
    135135   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    136136   !! $Header$  
    137    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     137   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    138138   !!---------------------------------------------------------------------- 
    139139 
    140140CONTAINS 
    141141 
    142 #if !defined key_agrif 
     142#if defined key_agrif 
     143   SUBROUTINE stp( ) 
     144#else 
    143145   SUBROUTINE stp( kstp ) 
    144 #else 
    145    SUBROUTINE stp( ) 
    146146#endif 
    147147      !!---------------------------------------------------------------------- 
     
    160160      !!              -7- Compute the diagnostics variables (rd,N2, div,cur,w) 
    161161      !!              -8- Outputs and diagnostics 
    162       !! 
    163162      !!---------------------------------------------------------------------- 
    164163      !! * Arguments 
    165 #if !defined key_agrif    
     164#if defined key_agrif    
     165      INTEGER               :: kstp   ! ocean time-step index 
     166#else 
    166167      INTEGER, INTENT( in ) :: kstp   ! ocean time-step index 
    167 #else 
    168       INTEGER               :: kstp   ! ocean time-step index 
    169168#endif       
    170169 
     
    176175      kstp = nit000 + Agrif_Nb_Step() 
    177176!      IF ( Agrif_Root() .and. lwp) Write(*,*) '---' 
    178 !      IF (lwp) Write(*,*) 'Grid N°',Agrif_Fixed(),' time step ',kstp 
     177!      IF (lwp) Write(*,*) 'Grid Number',Agrif_Fixed(),' time step ',kstp 
    179178#endif    
    180179      indic = 1                    ! reset to no error condition 
     180 
    181181      adatrj = adatrj + rdt/86400._wp 
    182182 
    183183      CALL day( kstp )             ! Calendar 
     184 
     185      CALL rst_opn( kstp )         ! Open the restart file 
    184186 
    185187      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     
    326328                             CALL tra_ldf    ( kstp )       ! lateral mixing 
    327329#if defined key_agrif 
    328       IF (.NOT. Agrif_Root())  CALL Agrif_Sponge_tra( kstp )          ! tracers sponge 
     330      IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_tra( kstp )  ! tracers sponge 
    329331#endif 
    330332                             CALL tra_zdf    ( kstp )       ! vertical mixing 
     
    363365                               CALL dyn_ldf( kstp )           ! lateral mixing 
    364366#if defined key_agrif 
    365       IF (.NOT. Agrif_Root())  CALL Agrif_Sponge_dyn( kstp )         ! momemtum sponge 
     367      IF(.NOT. Agrif_Root())   CALL Agrif_Sponge_dyn( kstp )  ! momemtum sponge 
    366368#endif 
    367369                               CALL dyn_hpg( kstp )           ! horizontal gradient of Hydrostatic pressure 
     
    400402 
    401403      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    402       ! Control, diagnostics and outputs 
    403       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    404       ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
    405       !----------------------------------------------------------------------- 
    406  
    407       !                                            ! Time loop: control and print 
    408                        CALL stp_ctl( kstp, indic ) 
    409                        IF ( indic < 0 ) CALL ctl_stop( 'step: indic < 0' )  
    410  
    411       IF ( nstop == 0 ) THEN 
    412          !                                         ! Diagnostics: 
     404      ! Control, and restarts 
     405      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     406      ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
     407      !----------------------------------------------------------------------- 
     408      !                                                     ! Time loop: control and print 
     409                                 CALL stp_ctl( kstp, indic ) 
     410      IF( indic < 0          )   CALL ctl_stop( 'step: indic < 0' ) 
     411 
     412      IF( kstp == nit000     )   CALL iom_close( numror )             ! close input  ocean restart file 
     413      IF( lrst_oce           )   CALL rst_write  ( kstp )             ! write output ocean restart file 
     414      IF( lk_obc             )   CALL obc_rst_wri( kstp )             ! write open boundary restart file 
     415      IF( lk_trdmld          )   CALL trd_mld_rst_write( kstp )       ! ocean model: restart file output for trends diagnostics 
     416 
     417 
     418      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     419      ! diagnostics and outputs 
     420      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     421      ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
     422      !----------------------------------------------------------------------- 
     423 
     424      IF ( nstop == 0 ) THEN                                ! Diagnostics 
    413425         IF( lk_floats  )   CALL flo_stp( kstp )                 ! drifting Floats 
    414426         IF( lk_trddyn  )   CALL trd_dwr( kstp )                 ! trends: dynamics  
     
    423435         IF( ln_diaptr  )   CALL dia_ptr( kstp )                 ! Poleward TRansports diagnostics 
    424436 
    425          !                                         ! save and outputs 
    426                             CALL rst_write  ( kstp )             ! ocean model: restart file output 
    427          IF( lk_obc     )   CALL obc_rst_wri( kstp )             ! ocean model: open boundary restart file output 
    428          IF( lk_trdmld  )   CALL trd_mld_rst_write( kstp )       ! ocean model: restart file output for trends diagnostics 
     437         !                                                 ! Outputs 
    429438                            CALL dia_wri    ( kstp, indic )      ! ocean model: outputs 
    430  
    431439      ENDIF 
    432440 
     
    436444 
    437445      IF( lk_cpl    )   CALL cpl_stp( kstp )                 ! coupled mode : field exchanges 
    438  
     446      ! 
    439447   END SUBROUTINE stp 
    440448 
Note: See TracChangeset for help on using the changeset viewer.