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 3294 for trunk/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

Ignore:
Timestamp:
2012-01-28T17:44:18+01:00 (12 years ago)
Author:
rblod
Message:

Merge of 3.4beta into the trunk

Location:
trunk/NEMOGCM/NEMO/OPA_SRC/DOM
Files:
12 edited
1 copied

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90

    r2715 r3294  
    3434   USE restart         !  
    3535   USE trc_oce, ONLY : lk_offline ! offline flag 
     36   USE timing          ! Timing 
    3637 
    3738   IMPLICIT NONE 
     
    7071      REAL(wp) ::   zjul 
    7172      !!---------------------------------------------------------------------- 
    72  
     73      ! 
    7374      ! all calendar staff is based on the fact that MOD( rday, rdttra(1) ) == 0 
    7475      IF( MOD( rday     , rdttra(1) ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' ) 
     
    127128      ! call day to set the calendar parameters at the begining of the current simulaton. needed by iom_init 
    128129      CALL day( nit000 ) 
    129  
     130      ! 
    130131   END SUBROUTINE day_init 
    131132 
     
    204205      REAL(wp)           ::   zprec      ! fraction of day corresponding to 0.1 second 
    205206      !!---------------------------------------------------------------------- 
     207      ! 
     208      IF( nn_timing == 1 )  CALL timing_start('day') 
     209      ! 
    206210      zprec = 0.1 / rday 
    207211      !                                                 ! New time-step 
     
    255259      IF( .NOT. lk_offline ) CALL rst_opn( kt )               ! Open the restart file if needed and control lrst_oce 
    256260      IF( lrst_oce         ) CALL day_rst( kt, 'WRITE' )      ! write day restart information 
     261      ! 
     262      IF( nn_timing == 1 )  CALL timing_stop('day') 
    257263      ! 
    258264   END SUBROUTINE day 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90

    r2715 r3294  
    150150   LOGICAL, PUBLIC, PARAMETER ::   lk_vvl = .FALSE.   !: fixed grid flag 
    151151#endif 
    152    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hur  , hvr    !: inverse of u and v-points ocean depth (1/m) 
    153    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu   , hv     !: depth at u- and v-points (meters) 
    154    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu_0 , hv_0   !: refernce depth at u- and v-points (meters) 
     152   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET ::   hur  , hvr    !: inverse of u and v-points ocean depth (1/m) 
     153   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)         ::   hu   , hv     !: depth at u- and v-points (meters) 
     154   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)         ::   hu_0 , hv_0   !: refernce depth at u- and v-points (meters) 
    155155 
    156156   INTEGER, PUBLIC ::   nla10              !: deepest    W level Above  ~10m (nlb10 - 1) 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r2528 r3294  
    3535   USE c1d             ! 1D vertical configuration 
    3636   USE dyncor_c1d      ! Coriolis term (c1d case)         (cor_c1d routine) 
     37   USE timing          ! Timing 
    3738 
    3839   IMPLICIT NONE 
     
    7071      !!---------------------------------------------------------------------- 
    7172      ! 
     73      IF( nn_timing == 1 )  CALL timing_start('dom_init') 
     74      ! 
    7275      IF(lwp) THEN 
    7376         WRITE(numout,*) 
     
    102105      IF( nmsh /= 0      )   CALL dom_wri      ! Create a domain file 
    103106      IF( .NOT.ln_rstart )   CALL dom_ctl      ! Domain control 
     107      ! 
     108      IF( nn_timing == 1 )  CALL timing_stop('dom_init') 
    104109      ! 
    105110   END SUBROUTINE dom_init 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domcfg.F90

    r2715 r3294  
    1515   USE in_out_manager  ! I/O manager 
    1616   USE lib_mpp         ! distributed memory computing library 
     17   USE timing          ! Timing 
    1718 
    1819   IMPLICIT NONE 
     
    3536      !! 
    3637      !!---------------------------------------------------------------------- 
    37  
     38      ! 
     39      IF( nn_timing == 1 )  CALL timing_start('dom_cfg') 
     40      ! 
    3841      IF(lwp) THEN                   ! Control print 
    3942         WRITE(numout,*) 
     
    5659      ! 
    5760      CALL dom_glo                   ! global domain versus zoom and/or local domain 
     61      ! 
     62      IF( nn_timing == 1 )  CALL timing_stop('dom_cfg') 
    5863      ! 
    5964   END SUBROUTINE dom_cfg 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90

    r2715 r3294  
    2525   USE in_out_manager ! I/O manager 
    2626   USE lib_mpp        ! MPP library 
     27   USE timing         ! Timing 
    2728 
    2829   IMPLICIT NONE 
     
    105106      REAL(wp) ::   zphi1, zsin_alpha, zim05, zjm05 
    106107      !!---------------------------------------------------------------------- 
    107  
     108      ! 
     109      IF( nn_timing == 1 )  CALL timing_start('dom_hgr') 
     110      ! 
    108111      IF(lwp) THEN 
    109112         WRITE(numout,*) 
     
    568571         IF( znorme > 1.e-13 ) CALL ctl_stop( ' ===>>>> : symmetrical condition: rerun with good equator line' ) 
    569572      ENDIF 
    570  
     573      ! 
     574      IF( nn_timing == 1 )  CALL timing_stop('dom_hgr') 
     575      ! 
    571576   END SUBROUTINE dom_hgr 
    572577 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r2715 r3294  
    2525   USE oce             ! ocean dynamics and tracers 
    2626   USE dom_oce         ! ocean space and time domain 
    27    USE obc_oce         ! ocean open boundary conditions 
    2827   USE in_out_manager  ! I/O manager 
    2928   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    3029   USE lib_mpp 
    3130   USE dynspg_oce      ! choice/control of key cpp for surface pressure gradient 
     31   USE wrk_nemo        ! Memory allocation 
     32   USE timing          ! Timing 
    3233 
    3334   IMPLICIT NONE 
     
    3839 
    3940   !                            !!* Namelist namlbc : lateral boundary condition * 
    40    REAL(wp) ::   rn_shlat = 2.   ! type of lateral boundary condition on velocity 
     41   REAL(wp)        :: rn_shlat   = 2.   ! type of lateral boundary condition on velocity 
     42   LOGICAL, PUBLIC :: ln_vorlat  = .false.   !  consistency of vorticity boundary condition  
     43   !                                            with analytical eqs. 
     44 
    4145 
    4246   INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  icoord ! Workspace for dom_msk_nsa() 
     
    125129      !!               tmask_i  : interior ocean mask 
    126130      !!---------------------------------------------------------------------- 
    127       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released, iwrk_in_use, iwrk_not_released 
    128       USE wrk_nemo, ONLY:   zwf  =>  wrk_2d_1      ! 2D real    workspace 
    129       USE wrk_nemo, ONLY:   imsk => iwrk_2d_1      ! 2D integer workspace 
    130131      ! 
    131132      INTEGER  ::   ji, jj, jk      ! dummy loop indices 
    132133      INTEGER  ::   iif, iil, ii0, ii1, ii   ! local integers 
    133134      INTEGER  ::   ijf, ijl, ij0, ij1       !   -       - 
    134       !! 
    135       NAMELIST/namlbc/ rn_shlat 
     135      INTEGER , POINTER, DIMENSION(:,:) ::  imsk 
     136      REAL(wp), POINTER, DIMENSION(:,:) ::  zwf 
     137      !! 
     138      NAMELIST/namlbc/ rn_shlat, ln_vorlat 
    136139      !!--------------------------------------------------------------------- 
    137        
    138       IF( wrk_in_use(2, 1) .OR. iwrk_in_use(2, 1) ) THEN 
    139          CALL ctl_stop('dom_msk: requested workspace arrays unavailable')   ;   RETURN 
    140       ENDIF 
    141  
     140      ! 
     141      IF( nn_timing == 1 )  CALL timing_start('dom_msk') 
     142      ! 
     143      CALL wrk_alloc( jpi, jpj, imsk ) 
     144      CALL wrk_alloc( jpi, jpj, zwf  ) 
     145      ! 
    142146      REWIND( numnam )              ! Namelist namlbc : lateral momentum boundary condition 
    143147      READ  ( numnam, namlbc ) 
     
    148152         WRITE(numout,*) '~~~~~~' 
    149153         WRITE(numout,*) '   Namelist namlbc' 
    150          WRITE(numout,*) '      lateral momentum boundary cond.    rn_shlat = ',rn_shlat 
     154         WRITE(numout,*) '      lateral momentum boundary cond.    rn_shlat  = ',rn_shlat 
     155         WRITE(numout,*) '      consistency with analytical form   ln_vorlat = ',ln_vorlat  
    151156      ENDIF 
    152157 
     
    390395      CALL lbc_lnk( fmask, 'F', 1._wp )      ! Lateral boundary conditions on fmask 
    391396 
     397      ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) 
    392398             
    393399      IF( nprint == 1 .AND. lwp ) THEN      ! Control print 
     
    436442      ENDIF 
    437443      ! 
    438       IF( wrk_not_released(2, 1)  .OR.   & 
    439          iwrk_not_released(2, 1)  )   CALL ctl_stop('dom_msk: failed to release workspace arrays') 
     444      CALL wrk_dealloc( jpi, jpj, imsk ) 
     445      CALL wrk_dealloc( jpi, jpj, zwf  ) 
     446      ! 
     447      IF( nn_timing == 1 )  CALL timing_stop('dom_msk') 
    440448      ! 
    441449   END SUBROUTINE dom_msk 
     
    460468      REAL(wp) ::   zaa 
    461469      !!--------------------------------------------------------------------- 
    462  
     470      ! 
     471      IF( nn_timing == 1 )  CALL timing_start('dom_msk_nsa') 
     472      ! 
    463473      IF(lwp) WRITE(numout,*) 
    464474      IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' 
     
    620630      ENDIF 
    621631      ! 
     632      IF( nn_timing == 1 )  CALL timing_stop('dom_msk_nsa') 
     633      ! 
    622634   END SUBROUTINE dom_msk_nsa 
    623635 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domngb.F90

    r2715 r3294  
    1111   !!---------------------------------------------------------------------- 
    1212   USE dom_oce        ! ocean space and time domain 
     13   USE in_out_manager ! I/O manager 
    1314   USE lib_mpp        ! for mppsum 
     15   USE wrk_nemo       ! Memory allocation 
     16   USE timing         ! Timing 
    1417 
    1518   IMPLICIT NONE 
     
    3437      !!                -> not good if located at too high latitude... 
    3538      !!---------------------------------------------------------------------- 
    36       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    37       USE wrk_nemo, ONLY:   zglam => wrk_2d_2 , zgphi => wrk_2d_3 , zmask => wrk_2d_4 , zdist => wrk_2d_5 
    3839      ! 
    3940      REAL(wp)        , INTENT(in   ) ::   plon, plat   ! longitude,latitude of the point 
     
    4344      INTEGER , DIMENSION(2) ::   iloc 
    4445      REAL(wp)               ::   zlon, zmini 
     46      REAL(wp), POINTER, DIMENSION(:,:) ::  zglam, zgphi, zmask, zdist 
    4547      !!-------------------------------------------------------------------- 
    4648      ! 
    47       IF( wrk_in_use(2, 2,3,4,5) )   CALL ctl_stop('dom_ngb: Requested workspaces already in use') 
     49      IF( nn_timing == 1 )  CALL timing_start('dom_ngb') 
     50      ! 
     51      CALL wrk_alloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 
    4852      ! 
    4953      zmask(:,:) = 0._wp 
     
    7276      ENDIF 
    7377      ! 
    74       IF( wrk_not_released(2, 2,3,4,5) )   CALL ctl_stop('dom_ngb: error releasing workspaces') 
     78      CALL wrk_dealloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 
     79      ! 
     80      IF( nn_timing == 1 )  CALL timing_stop('dom_ngb') 
    7581      ! 
    7682   END SUBROUTINE dom_ngb 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r2779 r3294  
    2020   USE lib_mpp         ! distributed memory computing library 
    2121   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
     22   USE wrk_nemo        ! Memory allocation 
     23   USE timing          ! Timing 
    2224 
    2325   IMPLICIT NONE 
    2426   PRIVATE 
    2527 
    26    PUBLIC   dom_vvl       ! called by domain.F90 
    27    PUBLIC   dom_vvl_alloc ! called by nemogcm.F90 
    28  
    29    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ee_t, ee_u, ee_v, ee_f   !: ??? 
    30    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mut , muu , muv , muf    !: ???  
     28   PUBLIC   dom_vvl         ! called by domain.F90 
     29   PUBLIC   dom_vvl_2       ! called by domain.F90 
     30   PUBLIC   dom_vvl_alloc   ! called by nemogcm.F90 
     31 
     32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mut , muu , muv , muf    !: 1/H_0 at t-,u-,v-,f-points  
    3133 
    3234   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:)     ::   r2dt   ! vertical profile time-step, = 2 rdttra  
     
    4951      ! 
    5052      ALLOCATE( mut (jpi,jpj,jpk) , muu (jpi,jpj,jpk) , muv (jpi,jpj,jpk) , muf (jpi,jpj,jpk) ,     & 
    51          &      ee_t(jpi,jpj)     , ee_u(jpi,jpj)     , ee_v(jpi,jpj)     , ee_f(jpi,jpj)     ,     & 
    5253         &      r2dt        (jpk)                                                             , STAT=dom_vvl_alloc ) 
    5354         ! 
     
    6263      !!                ***  ROUTINE dom_vvl  *** 
    6364      !!                    
    64       !! ** Purpose :  compute coefficients muX at T-U-V-F points to spread 
    65       !!              ssh over the whole water column (scale factors) 
    66       !!---------------------------------------------------------------------- 
    67       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    68       USE wrk_nemo, ONLY:   zs_t   => wrk_2d_1 , zs_u_1 => wrk_2d_2 , zs_v_1 => wrk_2d_3     ! 2D workspace 
     65      !! ** Purpose :   compute mu coefficients at t-, u-, v- and f-points to  
     66      !!              spread ssh over the whole water column (scale factors) 
     67      !!                set the before and now ssh at u- and v-points  
     68      !!              (also f-point in now case) 
     69      !!---------------------------------------------------------------------- 
    6970      ! 
    7071      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    71       REAL(wp) ::   zcoefu , zcoefv   , zcoeff                   ! local scalars 
    72       REAL(wp) ::   zv_t_ij, zv_t_ip1j, zv_t_ijp1, zv_t_ip1jp1   !   -      - 
    73       !!---------------------------------------------------------------------- 
    74  
    75       IF( wrk_in_use(2, 1,2,3) ) THEN 
    76          CALL ctl_stop('dom_vvl: requested workspace arrays unavailable')   ;   RETURN 
    77       ENDIF 
    78  
     72      REAL(wp) ::   zcoefu, zcoefv , zcoeff                ! local scalars 
     73      REAL(wp) ::   zvt   , zvt_ip1, zvt_jp1, zvt_ip1jp1   !   -      - 
     74      REAL(wp), POINTER, DIMENSION(:,:) ::  zee_t, zee_u, zee_v, zee_f   ! 2D workspace 
     75      !!---------------------------------------------------------------------- 
     76      ! 
     77      IF( nn_timing == 1 )  CALL timing_start('dom_vvl') 
     78      ! 
     79      CALL wrk_alloc( jpi, jpj, zee_t, zee_u, zee_v, zee_f ) 
     80      ! 
    7981      IF(lwp) THEN 
    8082         WRITE(numout,*) 
     
    9799 
    98100      !                                 !==  mu computation  ==! 
    99       ee_t(:,:) = fse3t_0(:,:,1)                ! Lower bound : thickness of the first model level 
    100       ee_u(:,:) = fse3u_0(:,:,1) 
    101       ee_v(:,:) = fse3v_0(:,:,1) 
    102       ee_f(:,:) = fse3f_0(:,:,1) 
     101      zee_t(:,:) = fse3t_0(:,:,1)                ! Lower bound : thickness of the first model level 
     102      zee_u(:,:) = fse3u_0(:,:,1) 
     103      zee_v(:,:) = fse3v_0(:,:,1) 
     104      zee_f(:,:) = fse3f_0(:,:,1) 
    103105      DO jk = 2, jpkm1                          ! Sum of the masked vertical scale factors 
    104          ee_t(:,:) = ee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk) 
    105          ee_u(:,:) = ee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk) 
    106          ee_v(:,:) = ee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk) 
     106         zee_t(:,:) = zee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk) 
     107         zee_u(:,:) = zee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk) 
     108         zee_v(:,:) = zee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk) 
    107109         DO jj = 1, jpjm1                      ! f-point : fmask=shlat at coasts, use the product of umask 
    108             ee_f(:,jj) = ee_f(:,jj) + fse3f_0(:,jj,jk) *  umask(:,jj,jk) * umask(:,jj+1,jk) 
     110            zee_f(:,jj) = zee_f(:,jj) + fse3f_0(:,jj,jk) *  umask(:,jj,jk) * umask(:,jj+1,jk) 
    109111         END DO 
    110112      END DO   
    111113      !                                         ! Compute and mask the inverse of the local depth at T, U, V and F points 
    112       ee_t(:,:) = 1. / ee_t(:,:) * tmask(:,:,1) 
    113       ee_u(:,:) = 1. / ee_u(:,:) * umask(:,:,1) 
    114       ee_v(:,:) = 1. / ee_v(:,:) * vmask(:,:,1) 
     114      zee_t(:,:) = 1._wp / zee_t(:,:) * tmask(:,:,1) 
     115      zee_u(:,:) = 1._wp / zee_u(:,:) * umask(:,:,1) 
     116      zee_v(:,:) = 1._wp / zee_v(:,:) * vmask(:,:,1) 
    115117      DO jj = 1, jpjm1                               ! f-point case fmask cannot be used  
    116          ee_f(:,jj) = 1. / ee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1) 
    117       END DO 
    118       CALL lbc_lnk( ee_f, 'F', 1. )                  ! lateral boundary condition on ee_f 
     118         zee_f(:,jj) = 1._wp / zee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1) 
     119      END DO 
     120      CALL lbc_lnk( zee_f, 'F', 1. )                 ! lateral boundary condition on ee_f 
    119121      ! 
    120122      DO jk = 1, jpk                            ! mu coefficients 
    121          mut(:,:,jk) = ee_t(:,:) * tmask(:,:,jk)     ! T-point at T levels 
    122          muu(:,:,jk) = ee_u(:,:) * umask(:,:,jk)     ! U-point at T levels 
    123          muv(:,:,jk) = ee_v(:,:) * vmask(:,:,jk)     ! V-point at T levels 
     123         mut(:,:,jk) = zee_t(:,:) * tmask(:,:,jk)     ! T-point at T levels 
     124         muu(:,:,jk) = zee_u(:,:) * umask(:,:,jk)     ! U-point at T levels 
     125         muv(:,:,jk) = zee_v(:,:) * vmask(:,:,jk)     ! V-point at T levels 
    124126      END DO 
    125127      DO jk = 1, jpk                                 ! F-point : fmask=shlat at coasts, use the product of umask 
    126128         DO jj = 1, jpjm1 
    127                muf(:,jj,jk) = ee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk)   ! at T levels 
    128          END DO 
    129          muf(:,jpj,jk) = 0.e0 
     129               muf(:,jj,jk) = zee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk)   ! at T levels 
     130         END DO 
     131         muf(:,jpj,jk) = 0._wp 
    130132      END DO 
    131133      CALL lbc_lnk( muf, 'F', 1. )                   ! lateral boundary condition 
     
    139141      END DO 
    140142       
    141       ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations 
    142       ! for ssh and scale factors 
    143       zs_t  (:,:) =         e1t(:,:) * e2t(:,:) 
    144       zs_u_1(:,:) = 0.5 / ( e1u(:,:) * e2u(:,:) ) 
    145       zs_v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) ) 
    146  
    147143      DO jj = 1, jpjm1                          ! initialise before and now Sea Surface Height at u-, v-, f-points 
    148144         DO ji = 1, jpim1   ! NO vector opt. 
    149             zcoefu = umask(ji,jj,1) * zs_u_1(ji,jj) 
    150             zcoefv = vmask(ji,jj,1) * zs_v_1(ji,jj) 
    151             zcoeff = 0.5 * umask(ji,jj,1) * umask(ji,jj+1,1) / ( e1f(ji,jj) * e2f(ji,jj) ) 
    152             ! before fields 
    153             zv_t_ij       = zs_t(ji  ,jj  ) * sshb(ji  ,jj  ) 
    154             zv_t_ip1j     = zs_t(ji+1,jj  ) * sshb(ji+1,jj  ) 
    155             zv_t_ijp1     = zs_t(ji  ,jj+1) * sshb(ji  ,jj+1) 
    156             sshu_b(ji,jj) = zcoefu * ( zv_t_ij + zv_t_ip1j ) 
    157             sshv_b(ji,jj) = zcoefv * ( zv_t_ij + zv_t_ijp1 ) 
    158             ! now fields 
    159             zv_t_ij       = zs_t(ji  ,jj  ) * sshn(ji  ,jj  ) 
    160             zv_t_ip1j     = zs_t(ji+1,jj  ) * sshn(ji+1,jj  ) 
    161             zv_t_ijp1     = zs_t(ji  ,jj+1) * sshn(ji  ,jj+1) 
    162             zv_t_ip1jp1   = zs_t(ji  ,jj+1) * sshn(ji  ,jj+1) 
    163             sshu_n(ji,jj) = zcoefu * ( zv_t_ij + zv_t_ip1j ) 
    164             sshv_n(ji,jj) = zcoefv * ( zv_t_ij + zv_t_ijp1 ) 
    165             sshf_n(ji,jj) = zcoeff * ( zv_t_ij + zv_t_ip1j + zv_t_ijp1 + zv_t_ip1jp1 ) 
     145            zcoefu = 0.50_wp / ( e1u(ji,jj) * e2u(ji,jj) ) * umask(ji,jj,1) 
     146            zcoefv = 0.50_wp / ( e1v(ji,jj) * e2v(ji,jj) ) * vmask(ji,jj,1) 
     147            zcoeff = 0.25_wp / ( e1f(ji,jj) * e2f(ji,jj) ) * umask(ji,jj,1) * umask(ji,jj+1,1) 
     148            ! 
     149            zvt           = e1e2t(ji  ,jj  ) * sshb(ji  ,jj  )    ! before fields 
     150            zvt_ip1       = e1e2t(ji+1,jj  ) * sshb(ji+1,jj  ) 
     151            zvt_jp1       = e1e2t(ji  ,jj+1) * sshb(ji  ,jj+1) 
     152            sshu_b(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 
     153            sshv_b(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 
     154            ! 
     155            zvt           = e1e2t(ji  ,jj  ) * sshn(ji  ,jj  )    ! now fields 
     156            zvt_ip1       = e1e2t(ji+1,jj  ) * sshn(ji+1,jj  ) 
     157            zvt_jp1       = e1e2t(ji  ,jj+1) * sshn(ji  ,jj+1) 
     158            zvt_ip1jp1    = e1e2t(ji+1,jj+1) * sshn(ji+1,jj+1) 
     159            sshu_n(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 
     160            sshv_n(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 
     161            sshf_n(ji,jj) = zcoeff * ( zvt + zvt_ip1 + zvt_jp1 + zvt_ip1jp1 ) 
    166162         END DO 
    167163      END DO 
     
    169165      CALL lbc_lnk( sshv_n, 'V', 1. )   ;   CALL lbc_lnk( sshv_b, 'V', 1. ) 
    170166      CALL lbc_lnk( sshf_n, 'F', 1. ) 
    171  
    172                                                 ! initialise before scale factors at (u/v)-points 
    173       ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 
    174       DO jk = 1, jpkm1 
    175          DO jj = 1, jpjm1 
    176             DO ji = 1, jpim1 
    177                zv_t_ij           = zs_t(ji  ,jj  ) * fse3t_b(ji  ,jj  ,jk) 
    178                zv_t_ip1j         = zs_t(ji+1,jj  ) * fse3t_b(ji+1,jj  ,jk) 
    179                zv_t_ijp1         = zs_t(ji  ,jj+1) * fse3t_b(ji  ,jj+1,jk) 
    180                fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 
    181                fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 
    182             END DO 
    183          END DO 
    184       END DO 
    185       CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. )               ! lateral boundary conditions 
    186       CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. ) 
    187       ! Add initial scale factor to scale factor anomaly 
    188       fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:) 
    189       fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:) 
    190       ! 
    191       IF( wrk_not_released(2, 1,2,3) )   CALL ctl_stop('dom_vvl: failed to release workspace arrays') 
     167      ! 
     168      CALL wrk_dealloc( jpi, jpj, zee_t, zee_u, zee_v, zee_f ) 
     169      ! 
     170      IF( nn_timing == 1 )  CALL timing_stop('dom_vvl') 
    192171      ! 
    193172   END SUBROUTINE dom_vvl 
    194173 
     174 
     175   SUBROUTINE dom_vvl_2( kt, pe3u_b, pe3v_b ) 
     176      !!---------------------------------------------------------------------- 
     177      !!                ***  ROUTINE dom_vvl_2  *** 
     178      !!                    
     179      !! ** Purpose :   compute the vertical scale factors at u- and v-points 
     180      !!              in variable volume case. 
     181      !! 
     182      !! ** Method  :   In variable volume case (non linear sea surface) the  
     183      !!              the vertical scale factor at velocity points is computed 
     184      !!              as the average of the cell surface weighted e3t. 
     185      !!                It uses the sea surface heigth so it have to be initialized 
     186      !!              after ssh is read/set 
     187      !!---------------------------------------------------------------------- 
     188      INTEGER                   , INTENT(in   ) ::   kt               ! ocean time-step index 
     189      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pe3u_b, pe3v_b   ! before vertical scale factor at u- & v-pts 
     190      ! 
     191      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     192      INTEGER  ::   iku, ikv     ! local integers     
     193      INTEGER  ::   ii0, ii1, ij0, ij1   ! temporary integers 
     194      REAL(wp) ::   zvt          ! local scalars 
     195      !!---------------------------------------------------------------------- 
     196      ! 
     197      IF( nn_timing == 1 )  CALL timing_start('dom_vvl_2') 
     198      ! 
     199      IF( lwp .AND. kt == nit000 ) THEN 
     200         WRITE(numout,*) 
     201         WRITE(numout,*) 'dom_vvl_2 : Variable volume, fse3t_b initialization' 
     202         WRITE(numout,*) '~~~~~~~~~ ' 
     203         pe3u_b(:,:,jpk) = fse3u_0(:,:,jpk) 
     204         pe3v_b(:,:,jpk) = fse3u_0(:,:,jpk) 
     205      ENDIF 
     206       
     207      DO jk = 1, jpkm1           ! set the before scale factors at u- & v-points 
     208         DO jj = 2, jpjm1 
     209            DO ji = fs_2, fs_jpim1 
     210               zvt = fse3t_b(ji,jj,jk) * e1e2t(ji,jj) 
     211               pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1e2t(ji+1,jj) ) / ( e1u(ji,jj) * e2u(ji,jj) ) 
     212               pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e1e2t(ji,jj+1) ) / ( e1v(ji,jj) * e2v(ji,jj) ) 
     213            END DO 
     214         END DO 
     215      END DO 
     216 
     217      ! Correct scale factors at locations that have been individually modified in domhgr 
     218      ! Such modifications break the relationship between e1e2t and e1u*e2u etc. Recompute 
     219      ! scale factors ignoring the modified metric. 
     220      !                                                ! ===================== 
     221      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN    ! ORCA R2 configuration 
     222         !                                             ! ===================== 
     223         IF( nn_cla == 0 ) THEN 
     224            ! 
     225            ii0 = 139   ;   ii1 = 140        ! Gibraltar Strait (e2u was modified) 
     226            ij0 = 102   ;   ij1 = 102    
     227            DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     228               DO jj = mj0(ij0), mj1(ij1) 
     229                  DO ji = mi0(ii0), mi1(ii1) 
     230                     zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     231                     pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     232                  END DO 
     233               END DO 
     234            END DO 
     235            ! 
     236            ii0 = 160   ;   ii1 = 160        ! Bab el Mandeb (e2u and e1v were modified) 
     237            ij0 =  88   ;   ij1 =  88    
     238            DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     239               DO jj = mj0(ij0), mj1(ij1) 
     240                  DO ji = mi0(ii0), mi1(ii1) 
     241                     zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     242                     pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     243                  END DO 
     244               END DO 
     245            END DO 
     246            DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     247               DO jj = mj0(ij0), mj1(ij1) 
     248                  DO ji = mi0(ii0), mi1(ii1) 
     249                     zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     250                     pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     251                  END DO 
     252               END DO 
     253            END DO 
     254         ENDIF 
     255 
     256         ii0 = 145   ;   ii1 = 146        ! Danish Straits (e2u was modified) 
     257         ij0 = 116   ;   ij1 = 116    
     258         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     259            DO jj = mj0(ij0), mj1(ij1) 
     260               DO ji = mi0(ii0), mi1(ii1) 
     261                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     262                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     263               END DO 
     264            END DO 
     265         END DO 
     266         ! 
     267      ENDIF 
     268         !                                             ! ===================== 
     269      IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN    ! ORCA R1 configuration 
     270         !                                             ! ===================== 
     271 
     272         ii0 = 281   ;   ii1 = 282        ! Gibraltar Strait (e2u was modified) 
     273         ij0 = 200   ;   ij1 = 200    
     274         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     275            DO jj = mj0(ij0), mj1(ij1) 
     276               DO ji = mi0(ii0), mi1(ii1) 
     277                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     278                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     279               END DO 
     280            END DO 
     281         END DO 
     282 
     283         ii0 = 314   ;   ii1 = 315        ! Bhosporus Strait (e2u was modified) 
     284         ij0 = 208   ;   ij1 = 208    
     285         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     286            DO jj = mj0(ij0), mj1(ij1) 
     287               DO ji = mi0(ii0), mi1(ii1) 
     288                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     289                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     290               END DO 
     291            END DO 
     292         END DO 
     293 
     294         ii0 =  44   ;   ii1 =  44        ! Lombok Strait (e1v was modified) 
     295         ij0 = 124   ;   ij1 = 125    
     296         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     297            DO jj = mj0(ij0), mj1(ij1) 
     298               DO ji = mi0(ii0), mi1(ii1) 
     299                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     300                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     301               END DO 
     302            END DO 
     303         END DO 
     304 
     305         ii0 =  48   ;   ii1 =  48        ! Sumba Strait (e1v was modified) [closed from bathy_11 on] 
     306         ij0 = 124   ;   ij1 = 125    
     307         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     308            DO jj = mj0(ij0), mj1(ij1) 
     309               DO ji = mi0(ii0), mi1(ii1) 
     310                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     311                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     312               END DO 
     313            END DO 
     314         END DO 
     315 
     316         ii0 =  53   ;   ii1 =  53        ! Ombai Strait (e1v was modified) 
     317         ij0 = 124   ;   ij1 = 125    
     318         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     319            DO jj = mj0(ij0), mj1(ij1) 
     320               DO ji = mi0(ii0), mi1(ii1) 
     321                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     322                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     323               END DO 
     324            END DO 
     325         END DO 
     326 
     327         ii0 =  56   ;   ii1 =  56        ! Timor Passage (e1v was modified) 
     328         ij0 = 124   ;   ij1 = 125    
     329         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     330            DO jj = mj0(ij0), mj1(ij1) 
     331               DO ji = mi0(ii0), mi1(ii1) 
     332                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     333                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     334               END DO 
     335            END DO 
     336         END DO 
     337 
     338         ii0 =  55   ;   ii1 =  55        ! West Halmahera Strait (e1v was modified) 
     339         ij0 = 141   ;   ij1 = 142    
     340         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     341            DO jj = mj0(ij0), mj1(ij1) 
     342               DO ji = mi0(ii0), mi1(ii1) 
     343                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     344                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     345               END DO 
     346            END DO 
     347         END DO 
     348 
     349         ii0 =  58   ;   ii1 =  58        ! East Halmahera Strait (e1v was modified) 
     350         ij0 = 141   ;   ij1 = 142    
     351         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     352            DO jj = mj0(ij0), mj1(ij1) 
     353               DO ji = mi0(ii0), mi1(ii1) 
     354                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     355                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     356               END DO 
     357            END DO 
     358         END DO 
     359 
     360         ! 
     361      ENDIF 
     362      !                                                ! ====================== 
     363      IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN   ! ORCA R05 configuration 
     364         !                                             ! ====================== 
     365         ii0 = 563   ;   ii1 = 564        ! Gibraltar Strait (e2u was modified) 
     366         ij0 = 327   ;   ij1 = 327    
     367         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     368            DO jj = mj0(ij0), mj1(ij1) 
     369               DO ji = mi0(ii0), mi1(ii1) 
     370                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     371                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     372               END DO 
     373            END DO 
     374         END DO 
     375         ! 
     376         ii0 = 627   ;   ii1 = 628        ! Bosphore Strait (e2u was modified) 
     377         ij0 = 343   ;   ij1 = 343    
     378         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     379            DO jj = mj0(ij0), mj1(ij1) 
     380               DO ji = mi0(ii0), mi1(ii1) 
     381                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     382                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     383               END DO 
     384            END DO 
     385         END DO 
     386         ! 
     387         ii0 =  93   ;   ii1 =  94        ! Sumba Strait (e2u was modified) 
     388         ij0 = 232   ;   ij1 = 232    
     389         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     390            DO jj = mj0(ij0), mj1(ij1) 
     391               DO ji = mi0(ii0), mi1(ii1) 
     392                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     393                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     394               END DO 
     395            END DO 
     396         END DO 
     397         ! 
     398         ii0 = 103   ;   ii1 = 103        ! Ombai Strait (e2u was modified) 
     399         ij0 = 232   ;   ij1 = 232    
     400         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     401            DO jj = mj0(ij0), mj1(ij1) 
     402               DO ji = mi0(ii0), mi1(ii1) 
     403                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     404                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     405               END DO 
     406            END DO 
     407         END DO 
     408         ! 
     409         ii0 =  15   ;   ii1 =  15        ! Palk Strait (e2u was modified) 
     410         ij0 = 270   ;   ij1 = 270    
     411         DO jk = 1, jpkm1                 ! set the before scale factors at u-points 
     412            DO jj = mj0(ij0), mj1(ij1) 
     413               DO ji = mi0(ii0), mi1(ii1) 
     414                  zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 
     415                  pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 
     416               END DO 
     417            END DO 
     418         END DO 
     419         ! 
     420         ii0 =  87   ;   ii1 =  87        ! Lombok Strait (e1v was modified) 
     421         ij0 = 232   ;   ij1 = 233    
     422         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     423            DO jj = mj0(ij0), mj1(ij1) 
     424               DO ji = mi0(ii0), mi1(ii1) 
     425                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     426                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     427               END DO 
     428            END DO 
     429         END DO 
     430         ! 
     431         ii0 = 662   ;   ii1 = 662        ! Bab el Mandeb (e1v was modified) 
     432         ij0 = 276   ;   ij1 = 276    
     433         DO jk = 1, jpkm1                 ! set the before scale factors at v-points 
     434            DO jj = mj0(ij0), mj1(ij1) 
     435               DO ji = mi0(ii0), mi1(ii1) 
     436                  zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 
     437                  pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 
     438               END DO 
     439            END DO 
     440         END DO 
     441         ! 
     442      ENDIF 
     443      ! End of individual corrections to scale factors 
     444 
     445      IF( ln_zps ) THEN          ! minimum of the e3t at partial cell level 
     446         DO jj = 2, jpjm1 
     447            DO ji = fs_2, fs_jpim1 
     448               iku = mbku(ji,jj) 
     449               ikv = mbkv(ji,jj) 
     450               pe3u_b(ji,jj,iku) = MIN( fse3t_b(ji,jj,iku), fse3t_b(ji+1,jj  ,iku) )  
     451               pe3v_b(ji,jj,ikv) = MIN( fse3t_b(ji,jj,ikv), fse3t_b(ji  ,jj+1,ikv) )  
     452            END DO 
     453         END DO 
     454      ENDIF 
     455 
     456      pe3u_b(:,:,:) = pe3u_b(:,:,:) - fse3u_0(:,:,:)      ! anomaly to avoid zero along closed boundary/extra halos 
     457      pe3v_b(:,:,:) = pe3v_b(:,:,:) - fse3v_0(:,:,:) 
     458      CALL lbc_lnk( pe3u_b(:,:,:), 'U', 1. )               ! lateral boundary conditions 
     459      CALL lbc_lnk( pe3v_b(:,:,:), 'V', 1. ) 
     460      pe3u_b(:,:,:) = pe3u_b(:,:,:) + fse3u_0(:,:,:)      ! recover the full scale factor 
     461      pe3v_b(:,:,:) = pe3v_b(:,:,:) + fse3v_0(:,:,:) 
     462      ! 
     463      IF( nn_timing == 1 )  CALL timing_stop('dom_vvl_2') 
     464      ! 
     465   END SUBROUTINE dom_vvl_2 
     466    
    195467#else 
    196468   !!---------------------------------------------------------------------- 
     
    200472   SUBROUTINE dom_vvl 
    201473   END SUBROUTINE dom_vvl 
     474   SUBROUTINE dom_vvl_2(kdum, pudum, pvdum ) 
     475      USE par_kind 
     476      INTEGER                   , INTENT(in   ) ::   kdum        
     477      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pudum, pvdum 
     478   END SUBROUTINE dom_vvl_2 
    202479#endif 
    203480 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90

    r2715 r3294  
    2020   USE lbclnk          ! lateral boundary conditions - mpp exchanges 
    2121   USE lib_mpp         ! MPP library 
     22   USE wrk_nemo        ! Memory allocation 
     23   USE timing          ! Timing 
    2224 
    2325   IMPLICIT NONE 
     
    6365      !!                                   masks, depth and vertical scale factors 
    6466      !!---------------------------------------------------------------------- 
    65       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    66       USE wrk_nemo, ONLY:   zprt  => wrk_2d_1 , zprw  => wrk_2d_2    ! 2D workspace 
    67       USE wrk_nemo, ONLY:   zdepu => wrk_3d_1 , zdepv => wrk_3d_2    ! 3D     - 
    6867      !! 
    6968      INTEGER           ::   inum0    ! temprary units for 'mesh_mask.nc' file 
     
    7877      CHARACTER(len=21) ::   clnam4   ! filename (vertical   mesh informations) 
    7978      INTEGER           ::   ji, jj, jk   ! dummy loop indices 
    80       !!---------------------------------------------------------------------- 
    81  
    82       IF( wrk_in_use(2, 1,2) .OR. wrk_in_use(3, 1,2) )THEN 
    83          CALL ctl_stop('dom_wri: requested workspace arrays unavailable')   ;   RETURN 
    84       END IF 
    85  
     79      !                                   !  workspaces 
     80      REAL(wp), POINTER, DIMENSION(:,:  ) :: zprt, zprw  
     81      REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepu, zdepv 
     82      !!---------------------------------------------------------------------- 
     83      ! 
     84      IF( nn_timing == 1 )  CALL timing_start('dom_wri') 
     85      ! 
     86      CALL wrk_alloc( jpi, jpj, zprt, zprw ) 
     87      CALL wrk_alloc( jpi, jpj, jpk, zdepu, zdepv ) 
     88      ! 
    8689      IF(lwp) WRITE(numout,*) 
    8790      IF(lwp) WRITE(numout,*) 'dom_wri : create NetCDF mesh and mask information file(s)' 
     
    260263      END SELECT 
    261264      ! 
    262       IF( wrk_not_released(2, 1,2)  .OR.   & 
    263           wrk_not_released(3, 1,2)  )   CALL ctl_stop('dom_wri: failed to release workspace arrays') 
     265      CALL wrk_dealloc( jpi, jpj, zprt, zprw ) 
     266      CALL wrk_dealloc( jpi, jpj, jpk, zdepu, zdepv ) 
     267      ! 
     268      IF( nn_timing == 1 )  CALL timing_stop('dom_wri') 
    264269      ! 
    265270   END SUBROUTINE dom_wri 
     
    275280      !!                2) check which elements have been changed 
    276281      !!---------------------------------------------------------------------- 
    277       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    278       USE wrk_nemo, ONLY:   ztstref => wrk_2d_3      ! array with different values for each element 
    279282      ! 
    280283      CHARACTER(len=1)        , INTENT(in   ) ::   cdgrd   !  
     
    284287      INTEGER  ::  ji       ! dummy loop indices 
    285288      LOGICAL, DIMENSION(SIZE(puniq,1),SIZE(puniq,2),1) ::  lldbl  ! store whether each point is unique or not 
    286       !!---------------------------------------------------------------------- 
    287  
    288       IF( wrk_in_use(2, 3) ) THEN 
    289          CALL ctl_stop('dom_uniq: requested workspace array unavailable')   ;   RETURN 
    290       ENDIF 
    291  
     289      REAL(wp), POINTER, DIMENSION(:,:) :: ztstref 
     290      !!---------------------------------------------------------------------- 
     291      ! 
     292      IF( nn_timing == 1 )  CALL timing_start('dom_uniq') 
     293      ! 
     294      CALL wrk_alloc( jpi, jpj, ztstref ) 
     295      ! 
    292296      ! build an array with different values for each element  
    293297      ! in mpp: make sure that these values are different even between process 
     
    304308      puniq(nldi:nlei,nldj:nlej) = REAL( COUNT( lldbl(nldi:nlei,nldj:nlej,:), dim = 3 ) , wp ) 
    305309      ! 
    306       IF( wrk_not_released(2, 3) )   CALL ctl_stop('dom_uniq: failed to release workspace array') 
     310      CALL wrk_dealloc( jpi, jpj, ztstref ) 
     311      ! 
     312      IF( nn_timing == 1 )  CALL timing_stop('dom_uniq') 
    307313      ! 
    308314   END SUBROUTINE dom_uniq 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r2950 r3294  
    3838   USE lbclnk            ! ocean lateral boundary conditions (or mpp link) 
    3939   USE lib_mpp           ! distributed memory computing library 
     40   USE wrk_nemo        ! Memory allocation 
     41   USE timing          ! Timing 
    4042 
    4143   IMPLICIT NONE 
     
    8688      NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco 
    8789      !!---------------------------------------------------------------------- 
    88  
     90      ! 
     91      IF( nn_timing == 1 )  CALL timing_start('dom_zgr') 
     92      ! 
    8993      REWIND( numnam )                 ! Read Namelist namzgr : vertical coordinate' 
    9094      READ  ( numnam, namzgr ) 
     
    139143      ENDIF 
    140144      ! 
     145      IF( nn_timing == 1 )  CALL timing_stop('dom_zgr') 
     146      ! 
    141147   END SUBROUTINE dom_zgr 
    142148 
     
    170176      REAL(wp) ::   za2, zkth2, zacr2      ! Values for optional double tanh function set from parameters  
    171177      !!---------------------------------------------------------------------- 
    172  
     178      ! 
     179      IF( nn_timing == 1 )  CALL timing_start('zgr_z') 
     180      ! 
    173181      ! Set variables from parameters 
    174182      ! ------------------------------ 
     
    280288      END DO 
    281289      ! 
     290      IF( nn_timing == 1 )  CALL timing_stop('zgr_z') 
     291      ! 
    282292   END SUBROUTINE zgr_z 
    283293 
     
    319329      REAL(wp) ::   r_bump , h_bump , h_oce   ! bump characteristics  
    320330      REAL(wp) ::   zi, zj, zh, zhmin         ! local scalars 
    321       INTEGER , DIMENSION(jpidta,jpjdta) ::   idta   ! global domain integer data 
    322       REAL(wp), DIMENSION(jpidta,jpjdta) ::   zdta   ! global domain scalar data 
    323       !!---------------------------------------------------------------------- 
    324  
     331      INTEGER , POINTER, DIMENSION(:,:) ::   idta   ! global domain integer data 
     332      REAL(wp), POINTER, DIMENSION(:,:) ::   zdta   ! global domain scalar data 
     333      !!---------------------------------------------------------------------- 
     334      ! 
     335      IF( nn_timing == 1 )  CALL timing_start('zgr_bat') 
     336      ! 
     337      CALL wrk_alloc( jpidta, jpjdta, idta ) 
     338      CALL wrk_alloc( jpidta, jpjdta, zdta ) 
     339      ! 
    325340      IF(lwp) WRITE(numout,*) 
    326341      IF(lwp) WRITE(numout,*) '    zgr_bat : defines level and meter bathymetry' 
     
    440455            CALL iom_close( inum ) 
    441456            !                                                ! ===================== 
    442             IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN    ! ORCA R1 configuration 
    443                ii0 =  51   ;   ii1 =  53                      
    444                ij0 = 142   ;   ij1 = 142                     ! ===================== 
    445                DO ji = mi0(ii0), mi1(ii1)                    ! Close Halmera Strait 
    446                   DO jj = mj0(ij0), mj1(ij1) 
    447                      bathy(ji,jj) = 0._wp 
    448                   END DO 
    449                END DO 
    450                IF(lwp) WRITE(numout,*) 
    451                IF(lwp) WRITE(numout,*) '      orca_r1: Halmera strait closed at i=',ii0,' j=',ij0,'->',ij1 
    452             ENDIF 
    453             !                                                ! ===================== 
    454457            IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN    ! ORCA R2 configuration 
    455458               !                                             ! ===================== 
     
    512515      ENDIF 
    513516      ! 
     517      CALL wrk_dealloc( jpidta, jpjdta, idta ) 
     518      CALL wrk_dealloc( jpidta, jpjdta, zdta ) 
     519      ! 
     520      IF( nn_timing == 1 )  CALL timing_stop('zgr_bat') 
     521      ! 
    514522   END SUBROUTINE zgr_bat 
    515523 
     
    589597      !!              - update bathy : meter bathymetry (in meters) 
    590598      !!---------------------------------------------------------------------- 
    591       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    592       USE wrk_nemo, ONLY:   zbathy => wrk_2d_1 
    593599      !! 
    594600      INTEGER ::   ji, jj, jl                    ! dummy loop indices 
    595601      INTEGER ::   icompt, ibtest, ikmax         ! temporary integers 
    596       !!---------------------------------------------------------------------- 
    597  
    598       IF( wrk_in_use(2, 1) ) THEN 
    599          CALL ctl_stop('zgr_bat_ctl: requested workspace array unavailable')   ;   RETURN 
    600       ENDIF 
    601  
     602      REAL(wp), POINTER, DIMENSION(:,:) ::  zbathy 
     603      !!---------------------------------------------------------------------- 
     604      ! 
     605      IF( nn_timing == 1 )  CALL timing_start('zgr_bat_ctl') 
     606      ! 
     607      CALL wrk_alloc( jpi, jpj, zbathy ) 
     608      ! 
    602609      IF(lwp) WRITE(numout,*) 
    603610      IF(lwp) WRITE(numout,*) '    zgr_bat_ctl : check the bathymetry' 
     
    702709      ENDIF 
    703710      ! 
    704       IF( wrk_not_released(2, 1) )   CALL ctl_stop('zgr_bat_ctl: failed to release workspace array') 
     711      CALL wrk_dealloc( jpi, jpj, zbathy ) 
     712      ! 
     713      IF( nn_timing == 1 )  CALL timing_stop('zgr_bat_ctl') 
    705714      ! 
    706715   END SUBROUTINE zgr_bat_ctl 
     
    719728      !!                                     (min value = 1 over land) 
    720729      !!---------------------------------------------------------------------- 
    721       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    722       USE wrk_nemo, ONLY:   zmbk => wrk_2d_1 
    723730      !! 
    724731      INTEGER ::   ji, jj   ! dummy loop indices 
    725       !!---------------------------------------------------------------------- 
    726       ! 
    727       IF( wrk_in_use(2, 1) ) THEN 
    728          CALL ctl_stop('zgr_bot_level: requested 2D workspace unavailable')   ;   RETURN 
    729       ENDIF 
     732      REAL(wp), POINTER, DIMENSION(:,:) ::  zmbk 
     733      !!---------------------------------------------------------------------- 
     734      ! 
     735      IF( nn_timing == 1 )  CALL timing_start('zgr_bot_level') 
     736      ! 
     737      CALL wrk_alloc( jpi, jpj, zmbk ) 
    730738      ! 
    731739      IF(lwp) WRITE(numout,*) 
     
    745753      zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
    746754      ! 
    747       IF( wrk_not_released(2, 1) )   CALL ctl_stop('zgr_bot_level: failed to release workspace array') 
     755      CALL wrk_dealloc( jpi, jpj, zmbk ) 
     756      ! 
     757      IF( nn_timing == 1 )  CALL timing_stop('zgr_bot_level') 
    748758      ! 
    749759   END SUBROUTINE zgr_bot_level 
     
    761771      !!---------------------------------------------------------------------- 
    762772      ! 
     773      IF( nn_timing == 1 )  CALL timing_start('zgr_zco') 
     774      ! 
    763775      DO jk = 1, jpk 
    764          fsdept(:,:,jk) = gdept_0(jk) 
    765          fsdepw(:,:,jk) = gdepw_0(jk) 
    766          fsde3w(:,:,jk) = gdepw_0(jk) 
    767          fse3t (:,:,jk) = e3t_0(jk) 
    768          fse3u (:,:,jk) = e3t_0(jk) 
    769          fse3v (:,:,jk) = e3t_0(jk) 
    770          fse3f (:,:,jk) = e3t_0(jk) 
    771          fse3w (:,:,jk) = e3w_0(jk) 
    772          fse3uw(:,:,jk) = e3w_0(jk) 
    773          fse3vw(:,:,jk) = e3w_0(jk) 
    774       END DO 
     776            gdept(:,:,jk) = gdept_0(jk) 
     777            gdepw(:,:,jk) = gdepw_0(jk) 
     778            gdep3w(:,:,jk) = gdepw_0(jk) 
     779            e3t (:,:,jk) = e3t_0(jk) 
     780            e3u (:,:,jk) = e3t_0(jk) 
     781            e3v (:,:,jk) = e3t_0(jk) 
     782            e3f (:,:,jk) = e3t_0(jk) 
     783            e3w (:,:,jk) = e3w_0(jk) 
     784            e3uw(:,:,jk) = e3w_0(jk) 
     785            e3vw(:,:,jk) = e3w_0(jk) 
     786      END DO 
     787      ! 
     788      IF( nn_timing == 1 )  CALL timing_stop('zgr_zco') 
    775789      ! 
    776790   END SUBROUTINE zgr_zco 
     
    822836      !!  Reference :   Pacanowsky & Gnanadesikan 1997, Mon. Wea. Rev., 126, 3248-3270. 
    823837      !!---------------------------------------------------------------------- 
    824       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    825       USE wrk_nemo, ONLY:   zprt => wrk_3d_1 
    826838      !! 
    827839      INTEGER  ::   ji, jj, jk       ! dummy loop indices 
     
    833845      REAL(wp) ::   zdiff            ! temporary scalar 
    834846      REAL(wp) ::   zrefdep          ! temporary scalar 
     847      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zprt 
    835848      !!--------------------------------------------------------------------- 
    836       !  
    837       IF( wrk_in_use(3, 1) ) THEN 
    838          CALL ctl_stop('zgr_zps: requested workspace unavailable.')   ;   RETURN 
    839       ENDIF 
    840  
     849      ! 
     850      IF( nn_timing == 1 )  CALL timing_start('zgr_zps') 
     851      ! 
     852      CALL wrk_alloc( jpi, jpj, jpk, zprt ) 
     853      ! 
    841854      IF(lwp) WRITE(numout,*) 
    842855      IF(lwp) WRITE(numout,*) '    zgr_zps : z-coordinate with partial steps' 
     
    10171030         WRITE(numout,*) 'domzgr e3t(mbathy)'      ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10181031         WRITE(numout,*) 
    1019          WRITE(numout,*) 'domzgr e3w(mbathy)'      ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
     1032         WRITE(numout,*) 'domzgr e3w(mbathy)'      ;   CALL prihre(zprt(:,:,2),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10201033         WRITE(numout,*) 
    1021          WRITE(numout,*) 'domzgr e3u(mbathy)'      ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
     1034         WRITE(numout,*) 'domzgr e3u(mbathy)'      ;   CALL prihre(zprt(:,:,3),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10221035         WRITE(numout,*) 
    1023          WRITE(numout,*) 'domzgr e3v(mbathy)'      ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
     1036         WRITE(numout,*) 'domzgr e3v(mbathy)'      ;   CALL prihre(zprt(:,:,4),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10241037         WRITE(numout,*) 
    1025          WRITE(numout,*) 'domzgr e3f(mbathy)'      ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
     1038         WRITE(numout,*) 'domzgr e3f(mbathy)'      ;   CALL prihre(zprt(:,:,5),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10261039         WRITE(numout,*) 
    1027          WRITE(numout,*) 'domzgr gdep3w(mbathy)'   ;   CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
     1040         WRITE(numout,*) 'domzgr gdep3w(mbathy)'   ;   CALL prihre(zprt(:,:,6),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 
    10281041      ENDIF   
    10291042      ! 
    1030       IF( wrk_not_released(3, 1) )   CALL ctl_stop('zgr_zps: failed to release workspace') 
     1043      CALL wrk_dealloc( jpi, jpj, jpk, zprt ) 
     1044      ! 
     1045      IF( nn_timing == 1 )  CALL timing_stop('zgr_zps') 
    10311046      ! 
    10321047   END SUBROUTINE zgr_zps 
     
    11161131      !! Reference : Madec, Lott, Delecluse and Crepon, 1996. JPO, 26, 1393-1408. 
    11171132      !!---------------------------------------------------------------------- 
    1118       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    1119       USE wrk_nemo, ONLY:   zenv => wrk_2d_1 , ztmp => wrk_2d_2 , zmsk  => wrk_2d_3 
    1120       USE wrk_nemo, ONLY:   zri  => wrk_2d_4 , zrj  => wrk_2d_5 , zhbat => wrk_2d_6 
    1121       USE wrk_nemo, ONLY:   gsigw3  => wrk_3d_1 
    1122       USE wrk_nemo, ONLY:   gsigt3  => wrk_3d_2 
    1123       USE wrk_nemo, ONLY:   gsi3w3  => wrk_3d_3 
    1124       USE wrk_nemo, ONLY:   esigt3  => wrk_3d_4 
    1125       USE wrk_nemo, ONLY:   esigw3  => wrk_3d_5 
    1126       USE wrk_nemo, ONLY:   esigtu3 => wrk_3d_6 
    1127       USE wrk_nemo, ONLY:   esigtv3 => wrk_3d_7 
    1128       USE wrk_nemo, ONLY:   esigtf3 => wrk_3d_8 
    1129       USE wrk_nemo, ONLY:   esigwu3 => wrk_3d_9 
    1130       USE wrk_nemo, ONLY:   esigwv3 => wrk_3d_10 
    11311133      ! 
    11321134      INTEGER  ::   ji, jj, jk, jl           ! dummy loop argument 
     
    11341136      REAL(wp) ::   zcoeft, zcoefw, zrmax, ztaper   ! temporary scalars 
    11351137      ! 
     1138      REAL(wp), POINTER, DIMENSION(:,:  ) :: zenv, ztmp, zmsk, zri, zrj, zhbat 
     1139      REAL(wp), POINTER, DIMENSION(:,:,:) :: gsigw3, gsigt3, gsi3w3 
     1140      REAL(wp), POINTER, DIMENSION(:,:,:) :: esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3            
    11361141 
    11371142      NAMELIST/namzgr_sco/ rn_sbot_max, rn_sbot_min, rn_theta, rn_thetb, rn_rmax, ln_s_sigma, rn_bb, rn_hc 
    11381143      !!---------------------------------------------------------------------- 
    1139  
    1140       IF( wrk_in_use(2, 1,2,3,4,5,6) .OR. wrk_in_use(3, 1,2,3,4,5,6,7,8,9,10) ) THEN 
    1141          CALL ctl_stop('zgr_sco: ERROR - requested workspace arrays unavailable')   ;   RETURN 
    1142       ENDIF 
    1143  
     1144      ! 
     1145      IF( nn_timing == 1 )  CALL timing_start('zgr_sco') 
     1146      ! 
     1147      CALL wrk_alloc( jpi, jpj,      zenv, ztmp, zmsk, zri, zrj, zhbat                           ) 
     1148      CALL wrk_alloc( jpi, jpj, jpk, gsigw3, gsigt3, gsi3w3                                      ) 
     1149      CALL wrk_alloc( jpi, jpj, jpk, esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3 ) 
     1150      ! 
    11441151      REWIND( numnam )                       ! Read Namelist namzgr_sco : sigma-stretching parameters 
    11451152      READ  ( numnam, namzgr_sco ) 
     
    14941501 
    14951502      ! 
    1496 !!    H. Liu, POL. April 2009. Added for passing the scale check for the new released vvl code. 
     1503      where (e3t   (:,:,:).eq.0.0)  e3t(:,:,:) = 1.0 
     1504      where (e3u   (:,:,:).eq.0.0)  e3u(:,:,:) = 1.0 
     1505      where (e3v   (:,:,:).eq.0.0)  e3v(:,:,:) = 1.0 
     1506      where (e3f   (:,:,:).eq.0.0)  e3f(:,:,:) = 1.0 
     1507      where (e3w   (:,:,:).eq.0.0)  e3w(:,:,:) = 1.0 
     1508      where (e3uw  (:,:,:).eq.0.0)  e3uw(:,:,:) = 1.0 
     1509      where (e3vw  (:,:,:).eq.0.0)  e3vw(:,:,:) = 1.0 
     1510 
    14971511 
    14981512      fsdept(:,:,:) = gdept (:,:,:) 
     
    15901604!!gm bug    #endif 
    15911605      ! 
    1592       IF( wrk_not_released(2, 1,2,3,4,5,6) .OR. wrk_not_released(3, 1,2,3,4,5,6,7,8,9,10) )  & 
    1593         &  CALL ctl_stop('dom:zgr_sco: failed to release workspace arrays') 
     1606      CALL wrk_dealloc( jpi, jpj,      zenv, ztmp, zmsk, zri, zrj, zhbat                           ) 
     1607      CALL wrk_dealloc( jpi, jpj, jpk, gsigw3, gsigt3, gsi3w3                                      ) 
     1608      CALL wrk_dealloc( jpi, jpj, jpk, esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3 ) 
     1609      ! 
     1610      IF( nn_timing == 1 )  CALL timing_stop('zgr_sco') 
    15941611      ! 
    15951612   END SUBROUTINE zgr_sco 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r2777 r3294  
    1313   !!            2.0  !  2006-07  (S. Masson)  distributed restart using iom 
    1414   !!            3.3  !  2010-10  (C. Ethe) merge TRC-TRA 
     15   !!            3.4  !  2011-04  (G. Madec) Merge of dtatem and dtasal & suppression of tb,tn/sb,sn  
    1516   !!---------------------------------------------------------------------- 
    1617 
     
    3031   USE zdf_oce         ! ocean vertical physics 
    3132   USE phycst          ! physical constants 
    32    USE dtatem          ! temperature data                 (dta_tem routine) 
    33    USE dtasal          ! salinity data                    (dta_sal routine) 
     33   USE dtatsd          ! data temperature and salinity   (dta_tsd routine) 
    3434   USE restart         ! ocean restart                   (rst_read routine) 
    3535   USE in_out_manager  ! I/O manager 
     
    4242   USE dynspg_exp      ! pressure gradient schemes 
    4343   USE dynspg_ts       ! pressure gradient schemes 
    44    USE traswp          ! Swap arrays                      (tra_swp routine) 
    4544   USE lib_mpp         ! MPP library 
     45   USE wrk_nemo        ! Memory allocation 
     46   USE timing          ! Timing 
    4647 
    4748   IMPLICIT NONE 
     
    6869      ! - ML - needed for initialization of e3t_b 
    6970      INTEGER  ::  jk     ! dummy loop indice 
     71      !!---------------------------------------------------------------------- 
     72      ! 
     73      IF( nn_timing == 1 )  CALL timing_start('istate_init') 
     74      ! 
    7075 
    7176      IF(lwp) WRITE(numout,*) 
     
    7378      IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    7479 
    75       rhd  (:,:,:) = 0.e0 
    76       rhop (:,:,:) = 0.e0 
    77       rn2  (:,:,:) = 0.e0  
    78       ta   (:,:,:) = 0.e0     
    79       sa   (:,:,:) = 0.e0 
     80      CALL dta_tsd_init                       ! Initialisation of T & S input data 
     81 
     82      rhd  (:,:,:  ) = 0.e0 
     83      rhop (:,:,:  ) = 0.e0 
     84      rn2  (:,:,:  ) = 0.e0  
     85      tsa  (:,:,:,:) = 0.e0     
    8086 
    8187      IF( ln_rstart ) THEN                    ! Restart from a file 
     
    8389         neuler = 1                              ! Set time-step indicator at nit000 (leap-frog) 
    8490         CALL rst_read                           ! Read the restart file 
    85          CALL tra_swap                           ! swap 3D arrays (t,s)  in a 4D array (ts) 
     91         !                                       ! define e3u_b, e3v_b from e3t_b read in restart file 
     92         CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
    8693         CALL day_init                           ! model calendar (using both namelist and restart infos) 
    8794      ELSE 
     
    9299         CALL day_init                           ! model calendar (using both namelist and restart infos) 
    93100         !                                       ! Initialization of ocean to zero 
    94          !   before fields     !       now fields      
    95          sshb (:,:)   = 0.e0   ;   sshn (:,:)   = 0.e0 
    96          ub   (:,:,:) = 0.e0   ;   un   (:,:,:) = 0.e0 
    97          vb   (:,:,:) = 0.e0   ;   vn   (:,:,:) = 0.e0   
    98          rotb (:,:,:) = 0.e0   ;   rotn (:,:,:) = 0.e0 
    99          hdivb(:,:,:) = 0.e0   ;   hdivn(:,:,:) = 0.e0 
     101         !   before fields      !       now fields      
     102         sshb (:,:)   = 0._wp   ;   sshn (:,:)   = 0._wp 
     103         ub   (:,:,:) = 0._wp   ;   un   (:,:,:) = 0._wp 
     104         vb   (:,:,:) = 0._wp   ;   vn   (:,:,:) = 0._wp   
     105         rotb (:,:,:) = 0._wp   ;   rotn (:,:,:) = 0._wp 
     106         hdivb(:,:,:) = 0._wp   ;   hdivn(:,:,:) = 0._wp 
     107         ! 
     108         !                                       ! define e3u_b, e3v_b from e3t_b initialized in domzgr 
     109         CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
    100110         ! 
    101111         IF( cp_cfg == 'eel' ) THEN 
     
    103113         ELSEIF( cp_cfg == 'gyre' ) THEN          
    104114            CALL istate_gyre                     ! GYRE  configuration : start from pre-defined T-S fields 
    105          ELSE 
    106             !                                    ! Other configurations: Initial T-S fields 
    107 #if defined key_dtatem 
    108             CALL dta_tem( nit000 )                  ! read 3D temperature data 
    109             tb(:,:,:) = t_dta(:,:,:)   ;   tn(:,:,:) = t_dta(:,:,:) 
    110              
    111 #else 
    112             IF(lwp) WRITE(numout,*)                 ! analytical temperature profile 
    113             IF(lwp) WRITE(numout,*)'             Temperature initialization using an analytic profile' 
    114             CALL istate_tem 
    115 #endif 
    116 #if defined key_dtasal 
    117             CALL dta_sal( nit000 )                  ! read 3D salinity data 
    118             sb(:,:,:) = s_dta(:,:,:)   ;   sn(:,:,:) = s_dta(:,:,:) 
    119 #else 
    120             ! No salinity data 
    121             IF(lwp)WRITE(numout,*)                  ! analytical salinity profile 
    122             IF(lwp)WRITE(numout,*)'             Salinity initialisation using a constant value' 
    123             CALL istate_sal 
    124 #endif 
     115         ELSEIF( ln_tsd_init      ) THEN         ! Initial T-S fields read in files 
     116            CALL dta_tsd( nit000, tsb )                  ! read 3D T and S data at nit000 
     117            tsn(:,:,:,:) = tsb(:,:,:,:) 
     118            ! 
     119         ELSE                                    ! Initial T-S fields defined analytically 
     120            CALL istate_t_s 
    125121         ENDIF 
    126122         ! 
    127          CALL tra_swap                     ! swap 3D arrays (tb,sb,tn,sn)  in a 4D array 
    128123         CALL eos( tsb, rhd, rhop )        ! before potential and in situ densities 
    129124#if ! defined key_c1d 
     
    148143      ENDIF 
    149144      ! 
     145      IF( nn_timing == 1 )  CALL timing_stop('istate_init') 
     146      ! 
    150147   END SUBROUTINE istate_init 
    151148 
    152  
    153    SUBROUTINE istate_tem 
     149   SUBROUTINE istate_t_s 
    154150      !!--------------------------------------------------------------------- 
    155       !!                  ***  ROUTINE istate_tem  *** 
     151      !!                  ***  ROUTINE istate_t_s  *** 
    156152      !!    
    157153      !! ** Purpose :   Intialization of the temperature field with an  
    158154      !!      analytical profile or a file (i.e. in EEL configuration) 
    159155      !! 
    160       !! ** Method  :   Use Philander analytic profile of temperature 
     156      !! ** Method  : - temperature: use Philander analytic profile 
     157      !!              - salinity   : use to a constant value 35.5 
    161158      !! 
    162159      !! References :  Philander ??? 
    163160      !!---------------------------------------------------------------------- 
    164       INTEGER :: ji, jj, jk 
     161      INTEGER  :: ji, jj, jk 
     162      REAL(wp) ::   zsal = 35.50 
    165163      !!---------------------------------------------------------------------- 
    166164      ! 
    167165      IF(lwp) WRITE(numout,*) 
    168       IF(lwp) WRITE(numout,*) 'istate_tem : initial temperature profile' 
    169       IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     166      IF(lwp) WRITE(numout,*) 'istate_t_s : Philander s initial temperature profile' 
     167      IF(lwp) WRITE(numout,*) '~~~~~~~~~~   and constant salinity (',zsal,' psu)' 
    170168      ! 
    171169      DO jk = 1, jpk 
    172          DO jj = 1, jpj 
    173             DO ji = 1, jpi 
    174                tn(ji,jj,jk) = (  ( ( 7.5 - 0.*ABS(gphit(ji,jj))/30. )   & 
    175                   &               *( 1.-TANH((fsdept(ji,jj,jk)-80.)/30.) )   & 
    176                   &            + 10.*(5000.-fsdept(ji,jj,jk))/5000.)  ) * tmask(ji,jj,jk) 
    177                tb(ji,jj,jk) = tn(ji,jj,jk) 
    178           END DO 
    179         END DO 
     170         tsn(:,:,jk,jp_tem) = (  ( ( 7.5 - 0. * ABS( gphit(:,:) )/30. ) * ( 1.-TANH((fsdept(:,:,jk)-80.)/30.) )   & 
     171            &                + 10. * ( 5000. - fsdept(:,:,jk) ) /5000.)  ) * tmask(:,:,jk) 
     172         tsb(:,:,jk,jp_tem) = tsn(:,:,jk,jp_tem) 
    180173      END DO 
    181       ! 
    182       IF(lwp) CALL prizre( tn    , jpi   , jpj   , jpk   , jpj/2 ,   & 
    183          &                 1     , jpi   , 5     , 1     , jpk   ,   & 
    184          &                 1     , 1.    , numout                  ) 
    185       ! 
    186    END SUBROUTINE istate_tem 
    187  
    188  
    189    SUBROUTINE istate_sal 
    190       !!--------------------------------------------------------------------- 
    191       !!                  ***  ROUTINE istate_sal  *** 
    192       !! 
    193       !! ** Purpose :   Intialize the salinity field with an analytic profile 
    194       !! 
    195       !! ** Method  :   Use to a constant value 35.5 
    196       !!               
    197       !! ** Action  :   Initialize sn and sb 
    198       !!---------------------------------------------------------------------- 
    199       REAL(wp) ::   zsal = 35.50_wp 
    200       !!---------------------------------------------------------------------- 
    201       ! 
    202       IF(lwp) WRITE(numout,*) 
    203       IF(lwp) WRITE(numout,*) 'istate_sal : initial salinity : ', zsal 
    204       IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    205       ! 
    206       sn(:,:,:) = zsal * tmask(:,:,:) 
    207       sb(:,:,:) = sn(:,:,:) 
    208       ! 
    209    END SUBROUTINE istate_sal 
     174      tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 
     175      tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
     176      ! 
     177   END SUBROUTINE istate_t_s 
    210178 
    211179 
     
    254222            ! 
    255223            DO jk = 1, jpk 
    256                tn(:,:,jk) = ( zt2 + zt1 * exp( - fsdept(:,:,jk) / 1000 ) ) * tmask(:,:,jk) 
    257                tb(:,:,jk) = tn(:,:,jk) 
     224               tsn(:,:,jk,jp_tem) = ( zt2 + zt1 * exp( - fsdept(:,:,jk) / 1000 ) ) * tmask(:,:,jk) 
     225               tsb(:,:,jk,jp_tem) = tsn(:,:,jk,jp_tem) 
    258226            END DO 
    259227            ! 
    260             IF(lwp) CALL prizre( tn    , jpi   , jpj   , jpk   , jpj/2 ,   & 
    261                &                 1     , jpi   , 5     , 1     , jpk   ,   & 
    262                &                 1     , 1.    , numout                  ) 
     228            IF(lwp) CALL prizre( tsn(:,:,:,jp_tem), jpi   , jpj   , jpk   , jpj/2 ,   & 
     229               &                             1     , jpi   , 5     , 1     , jpk   ,   & 
     230               &                             1     , 1.    , numout                  ) 
    263231            ! 
    264232            ! set salinity field to a constant value 
     
    268236            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    269237            ! 
    270             sn(:,:,:) = zsal * tmask(:,:,:) 
    271             sb(:,:,:) = sn(:,:,:) 
     238            tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 
     239            tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    272240            ! 
    273241            ! set the dynamics: U,V, hdiv, rot (and ssh if necessary) 
     
    323291            ! 
    324292            CALL iom_open ( 'eel.initemp', inum ) 
    325             CALL iom_get ( inum, jpdom_data, 'initemp', tb ) ! read before temprature (tb) 
     293            CALL iom_get ( inum, jpdom_data, 'initemp', tsb(:,:,:,jp_tem) ) ! read before temprature (tb) 
    326294            CALL iom_close( inum ) 
    327295            ! 
    328             tn(:,:,:) = tb(:,:,:)                            ! set nox temperature to tb 
    329             ! 
    330             IF(lwp) CALL prizre( tn    , jpi   , jpj   , jpk   , jpj/2 ,   & 
    331                &                 1     , jpi   , 5     , 1     , jpk   ,   & 
    332                &                 1     , 1.    , numout                  ) 
     296            tsn(:,:,:,jp_tem) = tsb(:,:,:,jp_tem)                            ! set nox temperature to tb 
     297            ! 
     298            IF(lwp) CALL prizre( tsn(:,:,:,jp_tem), jpi   , jpj   , jpk   , jpj/2 ,   & 
     299               &                            1     , jpi   , 5     , 1     , jpk   ,   & 
     300               &                            1     , 1.    , numout                  ) 
    333301            ! 
    334302            ! set salinity field to a constant value 
     
    338306            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    339307            ! 
    340             sn(:,:,:) = zsal * tmask(:,:,:) 
    341             sb(:,:,:) = sn(:,:,:) 
     308            tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 
     309            tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    342310            ! 
    343311            !                                    ! =========================== 
     
    377345            DO jj = 1, jpj 
    378346               DO ji = 1, jpi 
    379                   tn(ji,jj,jk) = (  16. - 12. * TANH( (fsdept(ji,jj,jk) - 400) / 700 )         )   & 
     347                  tsn(ji,jj,jk,jp_tem) = (  16. - 12. * TANH( (fsdept(ji,jj,jk) - 400) / 700 )         )   & 
    380348                       &           * (-TANH( (500-fsdept(ji,jj,jk)) / 150 ) + 1) / 2               & 
    381349                       &       + (      15. * ( 1. - TANH( (fsdept(ji,jj,jk)-50.) / 1500.) )       & 
     
    383351                       &                + 7.  * (1500. - fsdept(ji,jj,jk)) / 1500.             )   &  
    384352                       &           * (-TANH( (fsdept(ji,jj,jk) - 500) / 150) + 1) / 2 
    385                   tn(ji,jj,jk) = tn(ji,jj,jk) * tmask(ji,jj,jk) 
    386                   tb(ji,jj,jk) = tn(ji,jj,jk) 
    387  
    388                   sn(ji,jj,jk) =  (  36.25 - 1.13 * TANH( (fsdept(ji,jj,jk) - 305) / 460 )  )  & 
     353                  tsn(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) 
     354                  tsb(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) 
     355 
     356                  tsn(ji,jj,jk,jp_sal) =  (  36.25 - 1.13 * TANH( (fsdept(ji,jj,jk) - 305) / 460 )  )  & 
    389357                     &              * (-TANH((500 - fsdept(ji,jj,jk)) / 150) + 1) / 2          & 
    390358                     &          + (  35.55 + 1.25 * (5000. - fsdept(ji,jj,jk)) / 5000.         & 
     
    393361                     &                + 0.2  * TANH( (fsdept(ji,jj,jk) - 1000.) / 5000.)    )  & 
    394362                     &              * (-TANH((fsdept(ji,jj,jk) - 500) / 150) + 1) / 2  
    395                   sn(ji,jj,jk) = sn(ji,jj,jk) * tmask(ji,jj,jk) 
    396                   sb(ji,jj,jk) = sn(ji,jj,jk) 
     363                  tsn(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 
     364                  tsb(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) 
    397365               END DO 
    398366            END DO 
     
    408376         ! ---------------------- 
    409377         CALL iom_open ( 'data_tem', inum ) 
    410          CALL iom_get ( inum, jpdom_data, 'votemper', tn )  
     378         CALL iom_get ( inum, jpdom_data, 'votemper', tsn(:,:,:,jp_tem) )  
    411379         CALL iom_close( inum ) 
    412380 
    413          tn(:,:,:) = tn(:,:,:) * tmask(:,:,:)  
    414          tb(:,:,:) = tn(:,:,:) 
     381         tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * tmask(:,:,:)  
     382         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) 
    415383 
    416384         ! Read salinity field 
    417385         ! ------------------- 
    418386         CALL iom_open ( 'data_sal', inum ) 
    419          CALL iom_get ( inum, jpdom_data, 'vosaline', sn )  
     387         CALL iom_get ( inum, jpdom_data, 'vosaline', tsn(:,:,:,jp_sal) )  
    420388         CALL iom_close( inum ) 
    421389 
    422          sn(:,:,:)  = sn(:,:,:) * tmask(:,:,:)  
    423          sb(:,:,:)  = sn(:,:,:) 
     390         tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) * tmask(:,:,:)  
     391         tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    424392 
    425393      END SELECT 
     
    429397         WRITE(numout,*) '              Initial temperature and salinity profiles:' 
    430398         WRITE(numout, "(9x,' level   gdept_0   temperature   salinity   ')" ) 
    431          WRITE(numout, "(10x, i4, 3f10.2)" ) ( jk, gdept_0(jk), tn(2,2,jk), sn(2,2,jk), jk = 1, jpk ) 
     399         WRITE(numout, "(10x, i4, 3f10.2)" ) ( jk, gdept_0(jk), tsn(2,2,jk,jp_tem), tsn(2,2,jk,jp_sal), jk = 1, jpk ) 
    432400      ENDIF 
    433401 
     
    446414      !!                 p=integral [ rau*g dz ] 
    447415      !!---------------------------------------------------------------------- 
    448       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    449       USE wrk_nemo, ONLY:   zprn => wrk_3d_1    ! 3D workspace 
    450  
    451416      USE dynspg          ! surface pressure gradient             (dyn_spg routine) 
    452417      USE divcur          ! hor. divergence & rel. vorticity      (div_cur routine) 
     
    456421      INTEGER ::   indic             ! ??? 
    457422      REAL(wp) ::   zmsv, zphv, zmsu, zphu, zalfg     ! temporary scalars 
    458       !!---------------------------------------------------------------------- 
    459  
    460       IF(wrk_in_use(3, 1) ) THEN 
    461          CALL ctl_stop('istate_uvg: requested workspace array unavailable')   ;   RETURN 
    462       ENDIF 
    463  
     423      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprn 
     424      !!---------------------------------------------------------------------- 
     425      ! 
     426      CALL wrk_alloc( jpi, jpj, jpk, zprn) 
     427      ! 
    464428      IF(lwp) WRITE(numout,*)  
    465429      IF(lwp) WRITE(numout,*) 'istate_uvg : Start from Geostrophy' 
     
    557521      rotb (:,:,:) = rotn (:,:,:)       ! set the before to the now value 
    558522      ! 
    559       IF( wrk_not_released(3, 1) ) THEN 
    560          CALL ctl_stop('istate_uvg: failed to release workspace array') 
    561       ENDIF 
     523      CALL wrk_dealloc( jpi, jpj, jpk, zprn) 
    562524      ! 
    563525   END SUBROUTINE istate_uvg 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90

    r2528 r3294  
    99   !!              -   !  2006-08  (G. Madec)  style  
    1010   !!             3.2  !  2006-08  (S. Masson, G. Madec)  suppress useless variables + style  
     11   !!             3.4  !  2011-11  (C. Harris)  minor changes for CICE constants  
    1112   !!---------------------------------------------------------------------- 
    1213 
     
    4849#endif 
    4950 
     51#if defined key_cice 
     52   REAL(wp), PUBLIC ::   rau0     = 1026._wp      !: reference volumic mass (density)  (kg/m3) 
     53#else 
    5054   REAL(wp), PUBLIC ::   rau0     = 1035._wp      !: reference volumic mass (density)  (kg/m3) 
     55#endif 
    5156   REAL(wp), PUBLIC ::   rau0r                    !: reference specific volume         (m3/kg) 
    5257   REAL(wp), PUBLIC ::   rcp      =    4.e+3_wp   !: ocean specific heat 
    5358   REAL(wp), PUBLIC ::   ro0cpr                   !: = 1. / ( rau0 * rcp ) 
    5459 
    55 #if defined key_lim3 
     60#if defined key_lim3 || defined key_cice 
    5661   REAL(wp), PUBLIC ::   rcdsn   =   0.31_wp      !: thermal conductivity of snow 
    5762   REAL(wp), PUBLIC ::   rcdic   =   2.034396_wp  !: thermal conductivity of fresh ice 
     
    100105      rsiyea = 365.25 * rday * 2. * rpi / 6.283076 
    101106      rsiday = rday / ( 1. + rday / rsiyea ) 
     107#if defined key_cice 
     108      omega =  7.292116e-05 
     109#else 
    102110      omega  = 2. * rpi / rsiday  
     111#endif 
    103112 
    104113      rau0r  = 1. /   rau0   
     
    155164         WRITE(numout,*) '          latent heat of fusion of fresh ice / snow = ', lfus    , ' J/kg' 
    156165         WRITE(numout,*) '          latent heat of subl.  of fresh ice / snow = ', lsub    , ' J/kg' 
     166#elif defined key_cice 
     167         WRITE(numout,*) '          latent heat of fusion of fresh ice / snow = ', lfus    , ' J/kg' 
    157168#else 
    158169         WRITE(numout,*) '          density times specific heat for snow      = ', rcpsn   , ' J/m^3/K'  
Note: See TracChangeset for help on using the changeset viewer.