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 3790 – NEMO

Changeset 3790


Ignore:
Timestamp:
2013-02-10T14:16:06+01:00 (11 years ago)
Author:
cetlod
Message:

2013/dev_r3411_CNRS4_IOCRS : major bugs corrections

Location:
branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM/CONFIG/ORCA2_LIM_CRS/EXP00/iodef.xml

    r3778 r3790  
    291291         <field id="soce_crs"   description="salinity"                unit="psu"  axis_ref="deptht" /> 
    292292   <field id="ssh_crs"    description="sea surface height"      unit="m"                      /> 
    293         <field id="ssh2_crs"    description="sea surface height"      unit="m"                      /> 
    294293   <field id="sst_crs"    description="sea surface temperature" unit="degC"                   /> 
    295294   <field id="sss_crs"    description="sea surface salinity"    unit="psu"                    /> 
     
    312311      <group id="gcrs_W"  axis_ref="depthw" grid_ref="grid_W_crs"> 
    313312   <field id="woce_crs"   description="ocean vertical velocity"   unit="m/s" /> 
    314         <field id="woce2_crs"   description="ocean vertical velocity"   unit="m/s" /> 
    315    <field id="wocet_crs"  description="ocean vertical velocity times temperature" unit="degC.m/s" /> 
    316    <field id="woces_crs"  description="ocean vertical velocity times salinity"    unit="psu.m/s" /> 
    317313      </group> 
    318314 
     
    348344             <field ref="toce_crs"    name="votemper" /> 
    349345             <field ref="soce_crs"    name="vosaline" /> 
    350              <field ref="ssh_crs"     name="sossheig" /> 
    351              <field ref="ssh2_crs"    name="sossheig2" /> 
    352              <field ref="hdiv_crs"    name="vohdiver" /> 
    353346              <field ref="sst_crs"     name="sosstsst" /> 
    354347        <field ref="sss_crs"     name="sosaline" /> 
     348             <field ref="ssh_crs"     name="sossheig" /> 
     349             <field ref="hdiv_crs"    name="vohdiver" /> 
    355350          </file> 
    356351<!--  
     
    433428           <field ref="toce_crs"    name="votemper" /> 
    434429           <field ref="soce_crs"    name="vosaline" /> 
    435            <field ref="ssh_crs"     name="sossheig" /> 
    436            <field ref="ssh2_crs"     name="sossheig2" /> 
    437430      <field ref="sst_crs"     name="sosstsst" /> 
    438431      <field ref="sss_crs"     name="sosaline" /> 
     432           <field ref="ssh_crs"     name="sossheig" /> 
    439433           <field ref="hdiv_crs"    name="vohdiver" /> 
    440434        </file>    
     
    478472        <file id="5d_gcrs_W" name="auto" description="ocean U grid coarsened variables" > 
    479473           <field ref="woce_crs"    name="vovecrtz"  /> 
    480            <field ref="woce2_crs"   name="vovecrtz2"  /> 
    481            <field ref="wocet_crs"   name="vovewoct" /> 
    482            <field ref="woces_crs"   name="vovewocs" /> 
    483474        </file>    
    484475    
  • branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM/CONFIG/ORCA2_LIM_CRS/MY_SRC/oce_trc.F90

    r3622 r3790  
    1212   !!---------------------------------------------------------------------- 
    1313 
    14    IF ( ln_crs ) THEN 
    1514 
    1615   !* Domain size * 
     
    544543# endif 
    545544 
    546    ENDIF 
    547545 
    548546#else 
  • branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM/NEMO/OPA_SRC/CRS/crs.F90

    r3779 r3790  
    5454CONTAINS 
    5555 
    56    SUBROUTINE crsfun_mask( p_ptmask, cd_type, p_pmask2d, p_cmask, p_cmask2d ) 
     56   SUBROUTINE crsfun_mask( p_ptmask, cd_type, p_cmask ) 
    5757      !!---------------------------------------------------------------- 
    5858      !!               *** SUBROUTINE crsfun_mask *** 
     
    7777      !!---------------------------------------------------------------- 
    7878      ! Arguments 
    79       REAL, DIMENSION(jpi,jpj,jpk)                 , INTENT(in)  :: p_ptmask  ! Parent grid tmask 
    80       CHARACTER(len=1)                             , INTENT(in)  :: cd_type   ! grid type (T,U,V,F) 
    81       REAL, DIMENSION(:,:),   OPTIONAL, ALLOCATABLE, INTENT(in)  :: p_pmask2d ! 2D mask, such as rnfmsk        
    82       REAL, DIMENSION(:,:,:), OPTIONAL, ALLOCATABLE, INTENT(out) :: p_cmask   ! Coarse grid [t,u,v,f]mask 
    83       REAL, DIMENSION(:,:),   OPTIONAL, ALLOCATABLE, INTENT(out) :: p_cmask2d ! Coarse grid [t,u,v,f]mask 
     79      REAL, DIMENSION(jpi,jpj,jpk)        , INTENT(in)  :: p_ptmask  ! Parent grid tmask 
     80      CHARACTER(len=1)                    , INTENT(in)  :: cd_type   ! grid type (T,U,V,F) 
     81      REAL, DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out) :: p_cmask   ! Coarse grid [t,u,v,f]mask 
    8482 
    8583 
     
    8987 
    9088      ! Initialize 
    91       IF ( PRESENT( p_cmask ) ) THEN 
    92          ALLOCATE( p_cmask(jpi_crs,jpj_crs,jpk) ) 
    93          p_cmask(:,:,:) = 0.0 
    94          ijpk = jpk 
    95       ENDIF 
    96  
    97       IF ( PRESENT( p_cmask2d ) ) THEN 
    98          ALLOCATE( p_cmask2d(jpi_crs,jpj_crs) ) 
    99          p_cmask2d(:,:) = 0.0 
    100          ijpk = 1 
    101       ENDIF 
    102   
    103  
    104       DO jk = 1, ijpk 
     89      DO jk = 1, jpk 
    10590         DO ji = 2, jpi_crsm1 
    10691            ijie = mie_crs(ji) 
     
    11196               ijjs = mjs_crs(jj)                    
    11297 
    113                IF ( PRESENT( p_cmask ) ) THEN 
    114                 
    11598                  IF ( cd_type == 'T' ) THEN   
    11699                     p_cmask(ji,jj,jk) = SUM( p_ptmask(ijis:ijie,ijjs:ijje,jk) )  
     
    128111                  ENDIF 
    129112 
    130                ENDIF 
    131  
    132                IF ( PRESENT( p_cmask2d ) ) THEN 
    133                   IF ( cd_type == 'T' ) THEN   
    134                      p_cmask2d(ji,jj) = SUM( p_pmask2d(ijis:ijie,ijjs:ijje) )  
    135                      IF ( p_cmask2d(ji,jj) > 0 )  p_cmask2d(ji,jj) = 1 
    136                   ELSEIF ( cd_type == 'V' ) THEN                
    137                      p_cmask2d(ji,jj) = SUM( p_pmask2d(ijis:ijie,ijje) * p_pmask2d(ijis:ijie,ijje+1) ) 
    138                      IF ( p_cmask2d(ji,jj) > 0 ) p_cmask2d(ji,jj) = 1 
    139                   ELSEIF ( cd_type == 'U' ) THEN 
    140                      p_cmask2d(ji,jj) = SUM( p_pmask2d(ijie,ijjs:ijje) * p_pmask2d(ijie+1,ijjs:ijje) )  
    141                      IF ( p_cmask2d(ji,jj) > 0 ) p_cmask2d(ji,jj) = 1 
    142                   ELSE   ! fmask 
    143                      p_cmask2d(ji,jj) =  p_pmask2d(ijie,ijje) + p_pmask2d(ijie+1,ijje)        & 
    144                          &            + p_pmask2d(ijie,ijje+1) + p_pmask2d(ijie+1,ijje+1)  
    145                      IF ( p_cmask2d(ji,jj) > 0 ) p_cmask2d(ji,jj) = 1 
    146                   ENDIF 
    147  
    148                ENDIF 
    149  
    150113            ENDDO 
    151114         ENDDO 
     
    155118! Retroactively add back the halo cells j=1, j=jpj_crs and i=1, i=jpi_crs. 
    156119 
    157       IF ( PRESENT( p_cmask ) ) THEN 
    158          IF( nperio /= 0 ) THEN 
    159             CALL crs_lbc_lnk( cd_type,1.0,pt3d1=p_cmask ) 
    160          ELSE 
    161             p_cmask( 1  ,:   ,:) = p_cmask(     2,:   ,:)          ! all points 
    162             p_cmask(jpi_crs,:,:) = p_cmask(jpi_crsm1,:,:) 
    163             p_cmask(   :,1   ,:) = p_cmask(     :,2   ,:) 
    164             p_cmask(:,jpj_crs,:) = p_cmask(:,jpj_crsm1,:) 
    165          ENDIF 
     120      IF( nperio /= 0 ) THEN 
     121          CALL crs_lbc_lnk( cd_type,1.0,pt3d1=p_cmask ) 
     122      ELSE 
     123         p_cmask( 1  ,:   ,:) = p_cmask(     2,:   ,:)          ! all points 
     124         p_cmask(jpi_crs,:,:) = p_cmask(jpi_crsm1,:,:) 
     125         p_cmask(   :,1   ,:) = p_cmask(     :,2   ,:) 
     126         p_cmask(:,jpj_crs,:) = p_cmask(:,jpj_crsm1,:) 
    166127      ENDIF 
    167  
    168       IF ( PRESENT( p_cmask2d ) ) THEN 
    169          IF( nperio /= 0 ) THEN 
    170             CALL crs_lbc_lnk( cd_type,1.0,pt2d=p_cmask2d ) 
    171          ELSE 
    172             p_cmask2d( 1  ,:   ) = p_cmask2d(     2,:   )          ! all points 
    173             p_cmask2d(jpi_crs,:) = p_cmask2d(jpi_crsm1,:) 
    174             p_cmask2d(   :,1   ) = p_cmask2d(     :,2   ) 
    175             p_cmask2d(:,jpj_crs) = p_cmask2d(:,jpj_crsm1) 
    176          ENDIF 
    177       ENDIF 
    178  
    179128 
    180129 
     
    207156      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  :: p_pglam  ! Parent grid longitude 
    208157      CHARACTER(len=1),             INTENT(in)  :: cd_type   ! grid type (T,U,V,F)  
    209       REAL(wp), ALLOCATABLE, DIMENSION(:,:), INTENT(out) :: p_cgphi  ! Coarse grid latitude 
    210       REAL(wp), ALLOCATABLE, DIMENSION(:,:), INTENT(out) :: p_cglam  ! Coarse grid longitude 
     158      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_cgphi  ! Coarse grid latitude 
     159      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_cglam  ! Coarse grid longitude 
    211160 
    212161      !! Local variables 
     
    217166 
    218167      !! Initialize output fields 
    219       IF ( .NOT. ALLOCATED( p_cgphi ) ) THEN 
    220         ALLOCATE( p_cgphi(jpi_crs,jpj_crs), p_cglam(jpi_crs,jpj_crs) ) 
    221       ENDIF 
    222168      p_cgphi(:,:) = 0.e0 
    223169      p_cglam(:,:) = 0.e0 
     
    343289      REAL(wp), DIMENSION(jpi,jpj),     INTENT(in) :: p_e2     ! Parent grid U,V scale factors (e2) 
    344290 
    345       REAL(wp), DIMENSION(:,:),   INTENT(out), ALLOCATABLE, OPTIONAL :: p_cfield2d_1 ! Coarse grid box 2D quantity 
    346       REAL(wp), DIMENSION(:,:),   INTENT(out), ALLOCATABLE, OPTIONAL :: p_cfield2d_2 ! Coarse grid box 2D quantity 
    347       REAL(wp), DIMENSION(:,:,:), INTENT(out), ALLOCATABLE, OPTIONAL :: p_cfield3d_1 ! Coarse grid box 3D quantity  
    348       REAL(wp), DIMENSION(:,:,:), INTENT(out), ALLOCATABLE, OPTIONAL :: p_cfield3d_2 ! Coarse grid box 3D quantity  
     291      REAL(wp), DIMENSION(:,:),   INTENT(out), OPTIONAL :: p_cfield2d_1 ! Coarse grid box 2D quantity 
     292      REAL(wp), DIMENSION(:,:),   INTENT(out), OPTIONAL :: p_cfield2d_2 ! Coarse grid box 2D quantity 
     293      REAL(wp), DIMENSION(:,:,:), INTENT(out), OPTIONAL :: p_cfield3d_1 ! Coarse grid box 3D quantity  
     294      REAL(wp), DIMENSION(:,:,:), INTENT(out), OPTIONAL :: p_cfield3d_2 ! Coarse grid box 3D quantity  
    349295 
    350296      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL :: p_fse3     ! Parent grid vertical level thickness (fse3u, fse3v) 
     
    381327      ze1(:,:) = p_e1(:,:) 
    382328      ze2(:,:) = p_e2(:,:) 
    383  
    384       IF ( PRESENT(p_cfield2d_1) ) ALLOCATE( p_cfield2d_1(jpi_crs,jpj_crs) ) 
    385       IF ( PRESENT(p_cfield2d_2) ) ALLOCATE( p_cfield2d_2(jpi_crs,jpj_crs) )  
    386       IF ( PRESENT(p_cfield3d_1) ) ALLOCATE( p_cfield3d_1(jpi_crs,jpj_crs,jpk) )  
    387       IF ( PRESENT(p_cfield3d_2) ) ALLOCATE( p_cfield3d_2(jpi_crs,jpj_crs,jpk) )  
    388329 
    389330      IF ( PRESENT(p_cfield2d_1) ) p_cfield2d_1(:,:) = 0.0 
     
    1012953         ze1e2(:,:) = p_e1e2t(:,:) 
    1013954         ze3(:,:,:) = 1.0 
    1014          IF ( PRESENT(p_pfield3d_1) ) THEN 
     955         IF ( PRESENT(p_pfield3d_2) ) THEN 
    1015956            ! Prep to do the area-weighted average of (3D) W 
    1016             zpfield3d(:,:,:) = p_pfield3d_1(:,:,:) 
     957            zpfield3d(:,:,:) = p_pfield3d_2(:,:,:) 
    1017958         ELSEIF ( PRESENT(p_pfield2d) ) THEN 
    1018959            ! Prep to do the area-weighted average of some 2D quantity  
     
    11251066            ENDDO 
    11261067 
    1127  
    11281068            IF ( ijpk == 1 ) THEN 
    11291069               IF ( PRESENT(p_cfield2d) ) p_cfield2d(:,:) = zcfield2d(:,:) * zcmask(:,:,jk) 
     
    11711111   END SUBROUTINE crsfun_TW 
    11721112 
    1173 SUBROUTINE crs_e3_max( p_e3, cd_type, p_mask, p_e3_crs) 
     1113   SUBROUTINE crs_e3_max( p_e3, cd_type, p_mask, p_e3_crs) 
    11741114      !!---------------------------------------------------------------- 
    11751115      !!               *** SUBROUTINE crsfun_TW *** 
     
    13261266 
    13271267 
    1328 SUBROUTINE crs_surf(p_e1, p_e2 ,p_e3, cd_type, p_mask, surf_crs, surf_msk_crs) 
     1268   SUBROUTINE crs_surf(p_e1, p_e2 ,p_e3, cd_type, p_mask, surf_crs, surf_msk_crs) 
    13291269      !!---------------------------------------------------------------- 
    13301270      !!               *** SUBROUTINE crsfun_TW *** 
     
    13861326      INTEGER ::  ji, jj, jk                   ! dummy loop indices 
    13871327      INTEGER :: ijie,ijis,ijje,ijjs,ijpk,jii,jjj 
    1388       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze1, ze2, ze3, zsurf_crs, zsurf_msk_crs, zpmask   
     1328      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze1, ze2, ze3   
     1329      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zsurf_crs, zsurf_msk_crs, zpmask   
    13891330      !!----------------------------------------------------------------   
    13901331      ! Initialize 
  • branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM/NEMO/OPA_SRC/CRS/crs_dom.F90

    r3778 r3790  
    108108      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE      :: tsn_crs 
    109109      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: un_crs, vn_crs, wn_crs 
    110       REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: ut_crs, vt_crs, wt_crs, us_crs, vs_crs, ws_crs  
    111       REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: rhd_crs, rhop_crs, hdivn_crs     
    112       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: sshb_crs, sshn_crs, ssha_crs     
    113       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: sshun_crs, sshvn_crs, sshfn_crs  
    114       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: sshub_crs, sshvb_crs, sshua_crs, sshva_crs  
    115       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: hu_crs, hv_crs  
    116       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: hdivbt_crs 
    117       REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: ssh_b_crs, ssh_a_crs, ssh_un_crs, ssh_vn_crs        ! instantaneous fields 
    118  
     110      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: hdivn_crs     
     111      REAL(wp), DIMENSION(:,:)    , ALLOCATABLE      :: sshn_crs     
    119112 
    120113      !  
    121114      ! Surface fluxes to pass to TOP 
    122       REAL(wp), DIMENSION(:,:)  , ALLOCATABLE      ::  utau_crs, vtau_crs, wndm_crs, qsr_crs 
    123       REAL(wp), DIMENSION(:,:)  , ALLOCATABLE      ::  del_emp_crs, sum_emp_crs 
    124       REAL(wp), DIMENSION(:,:)  , ALLOCATABLE      ::  emp_crs, emp_b_crs, emps_crs  
    125       REAL(wp), DIMENSION(:,:)  , ALLOCATABLE      ::  rnf_crs, fr_i_crs, h_rnf_crs 
    126      
    127       ! 
    128       ! Lateral diffusivity (tracers) to pass to TOP 
    129       REAL(wp)                                     ::  rldf_crs, rn_aht_0_crs, aht0_crs, ahtb0_crs  
    130  
    131 #if defined key_traldf_c3d 
    132       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) ::   ahtt_crs, ahtu_crs, ahtv_crs, ahtw_crs   !: 3D coefficients at T-,U-,V-,W-points 
    133       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) ::   aeiu_crs, aeiv_crs, aeiw_crs    
    134 #elif defined key_traldf_c2d 
    135       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   ahtt_crs, ahtu_crs, ahtv_crs, ahtw_crs   !: 2D coefficients at T-,U-,V-,W-points 
    136       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   aeiu_crs, aeiv_crs, aeiw_crs    
    137 #elif defined key_traldf_c1d 
    138       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:)     ::   ahtt_crs, ahtu_crs, ahtv_crs, ahtw_crs   !: 1D coefficients at T-,U-,V-,W-points 
    139       REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:)     ::   aeiu_crs, aeiv_crs, aeiw_crs    
    140 #else 
    141       REAL(wp), PUBLIC                                ::   ahtt_crs, ahtu_crs, ahtv_crs, ahtw_crs   !: scalar coefficients at T-,U-,V-,W-points 
    142 #endif 
    143  
    144   
     115      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE        ::  wndm_crs, qsr_crs, emp_crs, emps_crs 
     116 
    145117      ! Vertical diffusion 
    146118      REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)  ::  avt_crs           !: vert. diffusivity coef. [m2/s] at w-point for temp   
     
    153125 
    154126      ! Direction of lateral diffusion 
    155  
    156127 
    157128 
     
    164135      !!------------------------------------------------------------------- 
    165136      !! Local variables 
    166       INTEGER, DIMENSION(14) :: ierr 
     137      INTEGER, DIMENSION(15) :: ierr 
    167138 
    168139      ierr(:) = 0 
     
    219190 
    220191 
    221       ALLOCATE( un_crs(jpi_crs,jpj_crs,jpk) ,   vn_crs(jpi_crs,jpj_crs,jpk) , & 
    222          &      wn_crs(jpi_crs,jpj_crs,jpk) ,                                 & 
    223          &      ut_crs(jpi_crs,jpj_crs,jpk) ,   vt_crs(jpi_crs,jpj_crs,jpk) , & 
    224          &      us_crs(jpi_crs,jpj_crs,jpk) ,   vs_crs(jpi_crs,jpj_crs,jpk) , & 
    225          &      wt_crs(jpi_crs,jpj_crs,jpk) ,   ws_crs(jpi_crs,jpj_crs,jpk) , & 
    226          &      rhd_crs(jpi_crs,jpj_crs,jpk) , rhop_crs(jpi_crs,jpj_crs,jpk) , & 
    227          &      hdivn_crs(jpi_crs,jpj_crs,jpk), & 
    228          &      STAT=ierr(11)) 
    229  
    230       ALLOCATE( sshb_crs(jpi_crs,jpj_crs)   , sshn_crs(jpi_crs,jpj_crs)   , & 
    231          &      sshun_crs(jpi_crs,jpj_crs)  , sshvn_crs(jpi_crs,jpj_crs)  , & 
    232          &      sshfn_crs(jpi_crs,jpj_crs)  , emp_crs(jpi_crs,jpj_crs)    , & 
    233          &      del_emp_crs(jpi_crs,jpj_crs), sum_emp_crs(jpi_crs,jpj_crs), & 
    234          &      emp_b_crs(jpi_crs,jpj_crs)  , emps_crs(jpi_crs,jpj_crs)   , & 
    235          &      ssh_b_crs(jpi_crs,jpj_crs)  , ssh_a_crs(jpi_crs,jpj_crs)  , & 
    236          &      ssh_un_crs(jpi_crs,jpj_crs) , ssh_vn_crs(jpi_crs,jpj_crs)  , & 
    237          &      STAT=ierr(12)  ) 
     192      ALLOCATE( un_crs(jpi_crs,jpj_crs,jpk) , vn_crs(jpi_crs,jpj_crs,jpk) , & 
     193         &      wn_crs(jpi_crs,jpj_crs,jpk) , hdivn_crs(jpi_crs,jpj_crs,jpk), STAT=ierr(11)) 
     194 
     195      ALLOCATE( sshn_crs(jpi_crs,jpj_crs),  emp_crs(jpi_crs,jpj_crs)    , & 
     196         &      qsr_crs(jpi_crs,jpj_crs) ,  wndm_crs(jpi_crs,jpj_crs)    , & 
     197         &      emps_crs(jpi_crs,jpj_crs),        STAT=ierr(12)  ) 
    238198  
    239199      ALLOCATE( tsn_crs(jpi_crs,jpj_crs,jpk,jpts), avt_crs(jpi_crs,jpj_crs,jpk),    & 
     
    245205      ALLOCATE( nmln_crs(jpi_crs,jpj_crs) , hmld_crs(jpi_crs,jpj_crs) , & 
    246206         &      hmlp_crs(jpi_crs,jpj_crs) , hmlpt_crs(jpi_crs,jpj_crs) , STAT=ierr(14) ) 
     207 
    247208 
    248209      crs_dom_alloc = MAXVAL(ierr) 
  • branches/2013/dev_r3411_CNRS4_IOCRS/NEMOGCM/NEMO/OPA_SRC/CRS/crsdiawri.F90

    r3778 r3790  
    9292      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    9393      !! 
    94       INTEGER                      ::   ji, jj, jk              ! dummy loop indices 
    95       !! 
    96       REAL(wp), POINTER, DIMENSION(:,:)   :: z2d    ! 2D workspace 
    97       REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d    ! 3D workspace 
    98       REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d1   ! 3D workspace 
    99       REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3t ! 3D workspace for e3 
    100       REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3u ! 3D workspace for e3 
    101       REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3v ! 3D workspace for e3 
    102       REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3w ! 3D workspace for e3 
    103       REAL(wp), POINTER, DIMENSION(:,:)   :: ze1e2u ! 2D workspace 
    104       REAL(wp), POINTER, DIMENSION(:,:)   :: ze1e2v ! 2D workspace 
    105       REAL(wp), POINTER, DIMENSION(:,:,:) :: z3dcrs, zw ! 3D workspace for coarse grid !cc 
    106       REAL(wp), POINTER, DIMENSION(:,:)   :: z2dcrs, ssh_crs2 ! 2D workspace for coarse grid 
    107       INTEGER        :: ialloc, iiki, iikn 
    108       INTEGER, SAVE  :: itsct 
    109       REAL(wp)       :: zdtj, zrtsct 
    110 !!cc1 
     94      INTEGER               ::   ji, jj, jk              ! dummy loop indices 
     95      !! 
     96      REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3t, zfse3u, zfse3v, zfse3w ! 3D workspace for e3 
     97      REAL(wp), POINTER, DIMENSION(:,:,:) :: zt, zs  
     98      REAL(wp), POINTER, DIMENSION(:,:,:) :: zt_crs, zs_crs ! 
    11199      REAL(wp)       :: z2dcrsu, z2dcrsv 
    112 !!cc1 
    113       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:), SAVE :: ztsnm 
    114  
    115       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: zum, zvm, zwm  ! 3d work arrays for 
    116       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: zutm, zusm, zvtm, zvsm, zwtm, zwsm              ! accumulating sums for averaging 
    117       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: zvolw, zsuru, zsurv, zsurw, zhdivnm 
    118       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: zrn_absm, zrhopm, zrhdm, zrotbm 
    119       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: zavtm, zavsm, zuslpm,zvslpm, zwslpim, zwslpjm 
    120  
    121       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ), SAVE :: zsshbm, zsshm, zssham, zsshfm                    ! 2d work arrays  
    122       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ), SAVE :: zsshubm, zsshum, zsshuam, zsurtaux, zutaum       ! for accumulating sums for averaging 
    123       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ), SAVE :: zsshvbm, zsshvm, zsshvam, zsurtauy, zvtaum        
    124       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ), SAVE :: zhdivbtm, zempm, zemp_bm, zempsm 
    125       REAL(wp), ALLOCATABLE, DIMENSION(:,:  ), SAVE :: zwndmm, zqsrm, zrnfm, zfr_im, zh_rnfm 
    126  
    127       !! 
    128       REAL(wp)                                      :: zcoefu, zcoefv, zcoeff, z2dt, z1_2dt, z1_rau0, zraur 
    129       REAL(wp)                                      :: zij, zip1j, zijp1 
    130 !!cc1      REAL(wp), ALLOCATABLE, DIMENSION(:,:)         :: z2dcrsu, z2dcrsv, z2dcrsf, zhdivbt 
    131       REAL(wp), ALLOCATABLE, DIMENSION(:,:)         :: z2dcrsf, zhdivbt 
    132       REAL(wp), ALLOCATABLE, DIMENSION(:,:)         :: zsshub, zsshua, zsshvb, zsshva   ! temp work arrays for instantaneous fields 
    133       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)       :: zhdiv, zvolt_wgt, zrhd, zrhop, zavt 
    134       REAL(wp), ALLOCATABLE, DIMENSION(:,:)         :: zee_t, zee_f, zee_u, zee_v 
    135       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)       :: zmut, zmuf 
    136       REAL(wp), ALLOCATABLE, DIMENSION(:,:)         :: zcrs_surtaux_wgt , zcrs_surtauy_wgt 
    137       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:)     :: ztsn       
    138       REAL(wp) ::   zrho_c = 0.01_wp    ! density criterion for mixed layer depth 
    139       REAL(wp) ::   zavt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth 
    140  
    141       INTEGER , ALLOCATABLE, DIMENSION(:,:  )       :: imld 
    142  
    143       CHARACTER(LEN=80) :: clmode 
    144       !!---------------------------------------------------------------------- 
     100      !! 
     101       !!---------------------------------------------------------------------- 
    145102      !  
    146103 
    147104      IF( nn_timing == 1 )   CALL timing_start('crs_dia_wri') 
    148       !  
    149       ! 1. Initialize arrays 
    150       ! -------------------------------------------------------- 
    151       CALL wrk_alloc( jpi , jpj , z2d, ze1e2u, ze1e2v ) 
    152       CALL wrk_alloc( jpi , jpj , jpk , z3d, z3d1 ) 
    153       CALL wrk_alloc( jpi , jpj , jpk , zfse3u, zfse3v ) 
    154       CALL wrk_alloc( jpi , jpj , jpk , zfse3t, zfse3w ) 
    155       CALL wrk_alloc( jpi_crs , jpj_crs , jpk , z3dcrs, zw ) !cc 
    156       CALL wrk_alloc( jpi_crs , jpj_crs , z2dcrs, ssh_crs2 ) 
    157        
    158  
    159       IF ( .NOT. ALLOCATED(ztsnm) ) THEN 
    160          ALLOCATE( ztsnm(jpi_crs,jpj_crs,jpk,jpts), STAT=ialloc ) 
    161          IF ( ialloc /= 0 ) CALL ctl_warn('crsdiawri : failed to allocate cumulative 4d array') 
    162       ENDIF 
    163  
    164       IF ( .NOT. ALLOCATED(zum) ) THEN 
    165          ALLOCATE( zum(jpi_crs,jpj_crs,jpk) ,   zvm(jpi_crs,jpj_crs,jpk) ,  & 
    166             &    zavtm(jpi_crs,jpj_crs,jpk) , zvolw(jpi_crs,jpj_crs,jpk) ,  & 
    167             &    zsuru(jpi_crs,jpj_crs,jpk) , zsurv(jpi_crs,jpj_crs,jpk) ,  & 
    168             &    zsurw(jpi_crs,jpj_crs,jpk) ,   zwm(jpi_crs,jpj_crs,jpk) ,  & 
    169             &     zutm(jpi_crs,jpj_crs,jpk) ,  zusm(jpi_crs,jpj_crs,jpk) ,  & 
    170             &     zvtm(jpi_crs,jpj_crs,jpk) ,  zvsm(jpi_crs,jpj_crs,jpk) ,  & 
    171             &     zwtm(jpi_crs,jpj_crs,jpk) ,  zwsm(jpi_crs,jpj_crs,jpk) ,  & 
    172             &     zrhopm(jpi_crs,jpj_crs,jpk), zrhdm(jpi_crs,jpj_crs,jpk) ,  & 
    173             &   zhdivnm(jpi_crs,jpj_crs,jpk) ,  STAT=ialloc ) 
    174          IF( ialloc /= 0 )   CALL ctl_warn('crsdiawri : failed to allocate cumulative 3d arrays') 
    175       ENDIF 
    176  
    177       IF ( .NOT. ALLOCATED(zsshm) ) THEN 
    178          ALLOCATE( zsshm(jpi_crs,jpj_crs) ,   zsshfm(jpi_crs,jpj_crs) , & 
    179             &     zsshbm(jpi_crs,jpj_crs) ,   zssham(jpi_crs,jpj_crs) , & 
    180             &     zsshum(jpi_crs,jpj_crs) ,   zsshvm(jpi_crs,jpj_crs) , & 
    181             &   zhdivbtm(jpi_crs,jpj_crs) ,                             & 
    182             &      zempm(jpi_crs,jpj_crs) ,  zemp_bm(jpi_crs,jpj_crs) , & 
    183             &     zempsm(jpi_crs,jpj_crs) ,                             & 
    184             &     zvtaum(jpi_crs,jpj_crs) ,   zutaum(jpi_crs,jpj_crs) , & 
    185             &   zsurtaux(jpi_crs,jpj_crs) , zsurtauy(jpi_crs,jpj_crs) , & 
    186             &     zwndmm(jpi_crs,jpj_crs) ,    zqsrm(jpi_crs,jpj_crs) , & 
    187             &      zrnfm(jpi_crs,jpj_crs) ,   zfr_im(jpi_crs,jpj_crs) , & 
    188             &    zh_rnfm(jpi_crs,jpj_crs) ,                             & 
    189             &     STAT=ialloc ) 
    190          IF( ialloc /= 0 )   CALL ctl_warn('crsdiawri : failed to allocate cumulative 2d arrays') 
    191       ENDIF 
    192  
    193       IF ( .NOT. ALLOCATED(z2dcrsf) ) THEN 
    194 !!cc1      IF ( .NOT. ALLOCATED(z2dcrsu) ) THEN 
    195 !!cc1         ALLOCATE( z2dcrsu(jpi_crs,jpj_crs) , z2dcrsv(jpi_crs,jpj_crs) , & 
    196          ALLOCATE(   & 
    197             &      z2dcrsf(jpi_crs,jpj_crs) , zhdivbt(jpi_crs,jpj_crs) , & 
    198             &      zsshub(jpi_crs,jpj_crs)  , zsshvb(jpi_crs,jpj_crs)  , & 
    199             &      zsshua(jpi_crs,jpj_crs)  , zsshva(jpi_crs,jpj_crs)  , & 
    200             &      zee_t(jpi_crs,jpj_crs)   , zee_u(jpi_crs,jpj_crs)    , & 
    201             &      zee_v(jpi_crs,jpj_crs)   , zee_f(jpi_crs,jpj_crs)    , & 
    202             &      zmut(jpi_crs,jpj_crs,jpk), zmuf(jpi_crs,jpj_crs,jpk)    , & 
    203             &      zcrs_surtaux_wgt(jpi_crs,jpj_crs) , zcrs_surtauy_wgt(jpi_crs,jpj_crs) , & 
    204             &        imld(jpi_crs,jpj_crs)  ,                            & 
    205             &      STAT=ialloc ) 
    206          IF( ialloc /= 0 )   CALL ctl_warn('crsdiawri : failed to allocate temp 2d arrays') 
    207       ENDIF 
    208   
    209       IF ( .NOT. ALLOCATED(zhdiv) ) THEN 
    210          ALLOCATE( zhdiv(jpi_crs,jpj_crs,jpk) , zvolt_wgt(jpi_crs,jpj_crs,jpk) , & 
    211             &       zrhd(jpi_crs,jpj_crs,jpk) ,     zrhop(jpi_crs,jpj_crs,jpk) , & 
    212             &       zavt(jpi_crs,jpj_crs,jpk) , STAT=ialloc ) 
    213          IF( ialloc /= 0 )   CALL ctl_warn('crsdiawri : failed to allocate temp 3d arrays') 
    214       ENDIF 
    215  
    216       IF ( .NOT. ALLOCATED(ztsn) ) THEN 
    217          ALLOCATE( ztsn(jpi_crs,jpj_crs,jpk,jpts),    STAT=ialloc ) 
    218          IF( ialloc /= 0 )   CALL ctl_warn('crsdiawri : failed to allocate temp 4d arrays') 
    219       ENDIF 
    220  
    221       zw(:,:,:)=0.0 
    222       ! generic work arrays 
    223       z2d(:,:)      = 0.0 
    224       z3d(:,:,:)    = 0.0 
    225       z3d1(:,:,:)   = 0.0 
     105       
     106      !  Initialize arrays 
     107      CALL wrk_alloc( jpi, jpj, jpk, zfse3t, zfse3w ) 
     108      CALL wrk_alloc( jpi, jpj, jpk, zfse3u, zfse3v ) 
     109      CALL wrk_alloc( jpi, jpj, jpk, zt, zs       ) 
     110      ! 
     111      CALL wrk_alloc( jpi_crs, jpj_crs, jpk, zt_crs, zs_crs ) 
    226112 
    227113      ! Depth work arrrays 
     
    231117      zfse3w(:,:,:) = fse3w(:,:,:) 
    232118 
    233       ! Horizontal divergence work array 
    234       zhdiv(:,:,:) = 0._wp 
    235       zhdivbt(:,:) = 0._wp 
    236  
    237       zvolt_wgt(:,:,:) = 1._wp 
    238  
    239       ! Sea-surface height work arrays 
    240       zee_t(:,:) = 0._wp    ;   zee_f(:,:) = 0._wp 
    241       zee_u(:,:) = 0._wp    ;   zee_v(:,:) = 0._wp 
    242       zmut(:,:,:) = 0._wp   ;   zmuf(:,:,:) = 0._wp 
    243  
    244       ! work arrays 
    245       ze1e2u(:,:) = e1u(:,:) * e2u(:,:) 
    246       ze1e2v(:,:) = e1v(:,:) * e2v(:,:) 
    247   
    248119      IF( kt == nit000  ) THEN 
    249          ! reset work arrays for average computation 
    250          itsct = 0 
    251          ! 
    252          zum(:,:,:)   = 0._wp 
    253          zvm(:,:,:)   = 0._wp 
    254          ztsnm(:,:,:,:) = 0._wp 
    255          zwm(:,:,:)   = 0._wp 
    256          zutm(:,:,:)  = 0._wp 
    257          zusm(:,:,:)  = 0._wp 
    258          zvtm(:,:,:)  = 0._wp 
    259          zvsm(:,:,:)  = 0._wp 
    260          zwtm(:,:,:)  = 0._wp 
    261          zwsm(:,:,:)  = 0._wp 
    262          zavtm(:,:,:)  = 0._wp 
    263  
    264          zvolw(:,:,:) = 0._wp   ! t-centered grid box volume, masked 
    265          zsuru(:,:,:) = 0._wp   ! east face area, masked 
    266          zsurv(:,:,:) = 0._wp   ! north face area, masked 
    267          zsurw(:,:,:) = 0._wp   ! top face area, masked on gridT 
    268          zsurtaux(:,:) = 0._wp  ! top face area, masked on gridU 
    269          zsurtauy(:,:) = 0._wp  ! top face area, masked on gridV 
    270  
    271          zsshm(:,:)   = 0._wp 
    272          zsshum(:,:)  = 0._wp 
    273          zsshvm(:,:)  = 0._wp 
    274          zsshfm(:,:)  = 0._wp 
    275          zsshbm(:,:) = 0._wp 
    276          zsshubm(:,:) = 0._wp 
    277          zsshvbm(:,:) = 0._wp 
    278          zssham(:,:) = 0._wp 
    279          zsshuam(:,:) = 0._wp 
    280          zsshvam(:,:) = 0._wp 
    281          zwndmm(:,:) = 0._wp 
    282          zutaum(:,:) = 0._wp 
    283          zvtaum(:,:) = 0._wp 
    284          zqsrm(:,:) = 0._wp 
    285          zrnfm(:,:) = 0._wp 
    286          zh_rnfm(:,:) = 0._wp 
    287          zfr_im(:,:) = 0._wp 
    288  
    289          tsn_crs(:,:,:,:) = 0._wp               ! passive tracer array, now  
    290          un_crs(:,:,:) = 0._wp                  ! u-velocity 
    291          vn_crs(:,:,:) = 0._wp                  ! v-velocity 
    292          wn_crs(:,:,:) = 0._wp                  ! bottom boundary condition: w(:,:,jpk)=0 (set once for all)  
    293          sshn_crs(:,:) = 0._wp                  ! sea-surface height now 
    294          ut_crs(:,:,:) = 0._wp                  ! U*T 
    295          us_crs(:,:,:) = 0._wp                  ! U*S 
    296          vt_crs(:,:,:) = 0._wp                  ! V*T 
    297          vs_crs(:,:,:) = 0._wp                  ! V*S 
    298          avt_crs(:,:,:) = 0._wp                  ! avt 
    299  
    300          zcrs_surtaux_wgt(:,:) = ze1e2u(:,:) * umask(:,:,1) 
    301          zcrs_surtauy_wgt(:,:) = ze1e2v(:,:) * vmask(:,:,1) 
    302  
     120         tsn_crs(:,:,:,:) = 0._wp    ! temp/sal  array, now  
     121         un_crs (:,:,:  ) = 0._wp    ! u-velocity 
     122         vn_crs (:,:,:  ) = 0._wp    ! v-velocity 
     123         wn_crs (:,:,:  ) = 0._wp    ! w 
     124         avt_crs(:,:,:  ) = 0._wp    ! avt 
     125         sshn_crs(:,:   ) = 0._wp    ! ssh 
    303126      ENDIF 
    304127 
    305       ! 2. Coarsen fields and accumulate sums at each time step 
     128 
     129      ! 2. Coarsen fields at each time step 
    306130      ! -------------------------------------------------------- 
    307 ! jes, TODO. 31 aug 2012 
    308 ! need to pass the  
    309 ! 2D. utau, vtau, wndm, qsr, rnf, fi_i(:,:), h_rnf , nmln, hmld, hmlp, hmlpt 
    310 ! 3D. rn_abs(:,:,:), rhop, rhd=(rho-rau0)/rau0,  rotb, avt, avs, slopes,  
    311 ! DIM?. aht0, ahtu, ahtv, ahtw, ahtt, aeiu, aeiv, aeiw 
    312  
    313  
    314       itsct = itsct+1     ! counter for the number of time steps  
    315  
    316       ! 2.0  Weights 
    317       zvolw(:,:,:) = zvolw(:,:,:) + crs_volt_wgt(:,:,:)   ! masked ocean volume 
    318       zsurv(:,:,:) = zsurv(:,:,:) + crs_surfv_wgt(:,:,:)  ! masked north face surface area 
    319       zsuru(:,:,:) = zsuru(:,:,:) + crs_surfu_wgt(:,:,:)  ! masked east face surface area 
    320       zsurw(:,:,:) = zsurw(:,:,:) + crs_surfw_wgt(:,:,:)  ! masked top face surface area of gridT 
    321       zsurtaux(:,:) = zsurtaux(:,:) + zcrs_surtaux_wgt(:,:) ! masked top face surface area of gridU 
    322       zsurtauy(:,:) = zsurtauy(:,:) + zcrs_surtauy_wgt(:,:) ! masked top face surface area of gridV 
    323  
    324       ! 2.1 Temperature 
    325       !    2.1.1. now-temperature 
    326       z3d(:,:,:) = tsn(:,:,:,jp_tem) 
     131 
     132      !  Temperature 
     133      zt(:,:,:) = tsn(:,:,:,jp_tem)  ;      zt_crs(:,:,:) = 0._wp 
    327134      CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='VOL', p_cmask=tmask_crs, p_ptmask=tmask, & 
    328          &         p_pfield3d_1=zfse3t, p_pfield3d_2=z3d, p_cfield3d=z3dcrs )  
    329       ztsn(:,:,:,1) = z3dcrs(:,:,:) 
    330       ztsnm(:,:,:,1) = ztsnm(:,:,:,1) + ztsn(:,:,:,1)  
    331  
    332       ! 2.2 Salinity 
    333       !    2.2.1. now-salinity 
    334       z3d(:,:,:) = tsn(:,:,:,jp_sal)   ; z3dcrs(:,:,:) = 0.0  
     135         &         p_pfield3d_1=zfse3t, p_pfield3d_2=zt, p_cfield3d=zt_crs ) 
     136      tsn_crs(:,:,:,jp_tem) = zt_crs(:,:,:) 
     137      CALL crs_iom_put( "toce_crs", pv_r3d=tsn_crs(:,:,:,jp_tem) )    ! temp 
     138      CALL crs_iom_put( "sst_crs" , pv_r2d=tsn_crs(:,:,1,jp_tem) )    ! sst 
     139 
     140      !  Salinity 
     141      zt(:,:,:) = tsn(:,:,:,jp_sal)  ;      zt_crs(:,:,:) = 0._wp 
    335142      CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='VOL', p_cmask=tmask_crs, p_ptmask=tmask, & 
    336          &         p_pfield3d_1=zfse3t, p_pfield3d_2=z3d, p_cfield3d=z3dcrs ) 
    337       ztsn(:,:,:,2) = z3dcrs(:,:,:) 
    338       ztsnm(:,:,:,2) = ztsnm(:,:,:,2) + ztsn(:,:,:,2)  
    339  
    340       !    2.5.9 Vertical eddy diffusion coefficient (avt) on gridW 
    341       z3d(:,:,:) = avt(:,:,:)   ;  z3dcrs(:,:,:) = 0.0  
    342       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    343          &         p_pfield3d_1=zfse3t, p_pfield3d_2=z3d, p_cfield3d=z3dcrs ) 
    344       zavtm(:,:,:) = zavtm(:,:,:) + z3dcrs(:,:,:)        
    345  
    346       ! 2.3. U-Velocity 
    347       !    2.3.1. U 
    348       z3dcrs(:,:,:) = 0.0 
     143         &         p_pfield3d_1=zfse3t, p_pfield3d_2=zt, p_cfield3d=zt_crs ) 
     144      tsn_crs(:,:,:,jp_sal) = zt_crs(:,:,:) 
     145      CALL crs_iom_put( "soce_crs" , pv_r3d=tsn_crs(:,:,:,jp_sal) )    ! sal 
     146      CALL crs_iom_put( "sss_crs"  , pv_r2d=tsn_crs(:,:,1,jp_sal) )    ! sss 
     147 
     148      !  U-velocity 
    349149      CALL crsfun( p_e1_e2=e2u, cd_type='U', psgn=-1.0, p_pmask=umask, & 
    350          &         p_fse3=zfse3u, p_pfield=un, p_surf_crs=crs_surfu_wgt, p_cfield3d=z3dcrs ) 
    351       zum(:,:,:) = z3dcrs(:,:,:) 
    352       !    2.3.2.  UT, US 
    353       z3d(:,:,:) = 0.e0 ; z3d1(:,:,:) = 0.0 
     150         &         p_fse3=zfse3u, p_pfield=un, p_surf_crs=crs_surfu_wgt, p_cfield3d=un_crs ) 
     151      CALL crs_iom_put( "uoce_crs"  , pv_r3d=un_crs  )   ! i-current  
     152      ! 
     153      zt(:,:,:) = 0._wp     ;    zs(:,:,:) = 0._wp  ;   zt_crs(:,:,:) = 0._wp   ;    zs_crs(:,:,:) = 0._wp 
    354154      DO jk = 1, jpkm1 
    355155         DO jj = 2, jpjm1 
    356156            DO ji = 2, jpim1    
    357                z3d(ji,jj,jk)  = un(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )  
    358                z3d1(ji,jj,jk) = un(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) )  
     157               zt(ji,jj,jk)  = un(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_tem) + tsn(ji+1,jj,jk,jp_tem) )  
     158               zs(ji,jj,jk) = un(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_sal) + tsn(ji+1,jj,jk,jp_sal) )  
    359159            END DO 
    360160         END DO 
    361161      END DO 
    362       !   UT 
    363       z3dcrs(:,:,:) = 0.0 
    364       CALL crsfun(  p_e1_e2=e2u, cd_type='U', psgn=-1.0, p_pmask=umask, & 
    365          &         p_fse3=zfse3u, p_pfield=z3d, p_cfield3d=z3dcrs ) 
    366       zutm(:,:,:) = zutm(:,:,:) + z3dcrs(:,:,:) 
    367       !   US 
    368       z3dcrs(:,:,:) = 0.0 
    369       CALL crsfun(  p_e1_e2=e2u, cd_type='U', psgn=-1.0, p_pmask=umask, & 
    370          &         p_fse3=zfse3u, p_pfield=z3d1, p_cfield3d=z3dcrs ) 
    371       zusm(:,:,:) = zusm(:,:,:) + z3dcrs(:,:,:) 
    372  
    373  
    374       ! 2.4. V-Velocity 
    375       !    2.4.1. V 
    376       z3dcrs(:,:,:) = 0.0 
     162      CALL crsfun( p_e1_e2=e2u, cd_type='U', psgn=-1.0, p_pmask=umask, & 
     163         &         p_fse3=zfse3u, p_pfield=zt,  p_cfield3d=zt_crs ) 
     164      CALL crs_iom_put( "uocet_crs" , pv_r3d=zt_crs )   ! uT 
     165      CALL crsfun( p_e1_e2=e2u, cd_type='U', psgn=-1.0, p_pmask=umask, & 
     166         &         p_fse3=zfse3u, p_pfield=zs,  p_cfield3d=zs_crs ) 
     167      CALL crs_iom_put( "uoces_crs" , pv_r3d=zs_crs )   ! uS 
     168 
     169 
     170      !  V-velocity 
    377171      CALL crsfun( p_e1_e2=e1v, cd_type='V', psgn=-1.0, p_pmask=vmask, & 
    378          &         p_fse3=zfse3v, p_pfield=vn, p_surf_crs=crs_surfv_wgt, p_cfield3d=z3dcrs ) 
    379       zvm(:,:,:) = z3dcrs(:,:,:) 
    380       !    2.4.2  VT, VS 
    381       z3d(:,:,:) = 0.e0 ; z3d1(:,:,:) = 0.0 
     172         &         p_fse3=zfse3v, p_pfield=vn, p_surf_crs=crs_surfv_wgt, p_cfield3d=vn_crs ) 
     173      CALL crs_iom_put( "voce_crs"  , pv_r3d=vn_crs  )   ! v-current  
     174      !                                                                                  
     175      zt(:,:,:) = 0._wp     ;    zs(:,:,:) = 0._wp  ;   zt_crs(:,:,:) = 0._wp   ;    zs_crs(:,:,:) = 0._wp 
    382176      DO jk = 1, jpkm1 
    383177         DO jj = 2, jpjm1 
    384178            DO ji = 2, jpim1    
    385                z3d(ji,jj,jk)  = vn(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )  
    386                z3d1(ji,jj,jk) = vn(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) )  
     179               zt(ji,jj,jk)  = vn(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_tem) + tsn(ji,jj+1,jk,jp_tem) )  
     180               zs(ji,jj,jk) = vn(ji,jj,jk) * 0.5 * ( tsn(ji,jj,jk,jp_sal) + tsn(ji,jj+1,jk,jp_sal) )  
    387181            END DO 
    388182         END DO 
    389183      END DO 
    390       !   VT 
    391       z3dcrs(:,:,:) = 0.0 
    392184      CALL crsfun( p_e1_e2=e1v, cd_type='V', psgn=-1.0, p_pmask=vmask, & 
    393          &         p_fse3=zfse3v, p_pfield=z3d, p_cfield3d=z3dcrs ) 
    394       zvtm(:,:,:) = zvtm(:,:,:) + z3dcrs(:,:,:) 
    395       !   VS 
    396       z3dcrs(:,:,:) = 0.0 
     185         &         p_fse3=zfse3v, p_pfield=zt, p_cfield3d=zt_crs ) 
     186      CALL crs_iom_put( "vocet_crs" , pv_r3d=zt_crs )   ! vT 
     187 
    397188      CALL crsfun( p_e1_e2=e1v, cd_type='V', psgn=-1.0, p_pmask=vmask, & 
    398          &         p_fse3=zfse3v, p_pfield=z3d1, p_cfield3d=z3dcrs ) 
    399       zvsm(:,:,:) = zvsm(:,:,:) + z3dcrs(:,:,:) 
    400        
    401        
    402       ! Vitesse vertical !cc 
    403       z3dcrs(:,:,:) = 0.0 
    404       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    405         &           p_pfield3d_1=wn, p_cfield3d=z3dcrs) 
    406       zw(:,:,:) = z3dcrs(:,:,:) 
    407  
    408  
    409       ! 2.5. Surface boundary conditions 
    410       !    2.5.1 Evaporation minus Precipitation (emp) 
    411       z2d(:,:) = emp_b(:,:)   ; z2dcrs(:,:) = 0._wp 
    412       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    413          &         p_pfield2d=z2d, p_cfield2d=z2dcrs )    
    414       zemp_bm(:,:) = zemp_bm(:,:) + z2dcrs     ! summation of emp_b_crs 
    415       ! 
    416       z2d(:,:) = emp(:,:)   ; z2dcrs(:,:) = 0._wp 
    417       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    418          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    419       zempm(:,:) = zempm(:,:) + z2dcrs     ! summation of emp_crs 
    420       ! 
    421       z2d(:,:) = emps(:,:)   ; z2dcrs(:,:) = 0._wp 
    422       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    423          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    424       zempsm(:,:) = zempsm(:,:) + z2dcrs     ! summation of emps_crs 
    425       !         terms needed for ssh calculation 
    426       z2d(:,:) = emp_b(:,:)-emp(:,:) 
    427       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    428          &         p_pfield2d=z2d, p_cfield2d=del_emp_crs ) 
    429       z2d(:,:) = emp_b(:,:)+emp(:,:) 
    430       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
    431          &         p_pfield2d=z2d, p_cfield2d=sum_emp_crs ) 
    432       ! 
    433       !    2.5.2 Wind Stress (utau, vtau)  
    434       !          Taux 
    435       z2d(:,:) = utau(:,:)   ; z2dcrs(:,:) = 0._wp 
    436       CALL crsfun( p_e1e2t=ze1e2u, cd_type='U', cd_op='SUM', p_cmask=umask_crs, p_ptmask=umask, & 
    437          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    438       zutaum(:,:) = zutaum(:,:) + z2dcrs(:,:)       
    439       !          Tauy 
    440       z2d(:,:) = vtau(:,:)   ; z2dcrs(:,:) = 0._wp 
    441       CALL crsfun( p_e1e2t=ze1e2v, cd_type='V', cd_op='SUM', p_cmask=vmask_crs, p_ptmask=vmask, & 
    442          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    443       zvtaum(:,:) = zvtaum(:,:) + z2dcrs(:,:) 
    444   
    445       !    2.5.3 Scalar Wind (wndm) 
    446       z2d(:,:) = wndm(:,:)   ; z2dcrs(:,:) = 0._wp 
    447       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    448          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    449       zwndmm(:,:) = zwndmm(:,:) + z2dcrs(:,:) 
    450  
    451       !    2.5.4 Surface shortwave radiative flux (qsr) 
    452       z2d(:,:) = qsr(:,:)   ; z2dcrs(:,:) = 0._wp 
    453       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    454          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    455       zqsrm(:,:) = zqsrm(:,:) + z2dcrs(:,:) 
    456  
    457       !    2.5.5 Runoff freshwater flux 
    458       z2d(:,:) = rnf(:,:)   ; z2dcrs(:,:) = 0._wp 
    459       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    460          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    461       zrnfm(:,:) = zrnfm(:,:) + z2dcrs(:,:) 
    462  
    463       !    2.5.6 Depth of runoff in meters h_rnf(:,:)) 
    464       z2d(:,:) = h_rnf(:,:) ; z2dcrs(:,:) = 0._wp 
    465       CALL crsfun(p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    466          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    467       zh_rnfm(:,:) = zh_rnfm(:,:) + z2dcrs(:,:) 
    468  
    469       !    2.5.7 Sea-ice fraction (fr_i) 
    470       z2d(:,:) = fr_i(:,:)   ; z2dcrs(:,:) = 0._wp 
    471       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='SUM', p_cmask=tmask_crs, p_ptmask=tmask, & 
    472          &         p_pfield2d=z2d, p_cfield2d=z2dcrs ) 
    473       zfr_im(:,:) = zfr_im(:,:) + z2dcrs(:,:) 
    474  
    475  
    476  
    477       !    2.5.8 Potential density (rhop)  
    478 ! jes need to declare rho cumulative arrays 
    479   !    CALL crs_eos( ztsn, zrhd, zrhop ) 
    480  
    481       CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain 
    482  
    483       CALL eos( ztsn, zrhd, zrhop ) 
    484  
    485       CALL dom_grid_glo   ! Return to parent grid domain 
    486  
    487       zrhdm(:,:,:) = zrhdm(:,:,:) + zrhd 
    488       zrhopm(:,:,:) = zrhopm(:,:,:) + zrhop 
    489  
    490       ! 2.6. Diagnostic variables: horizontal divergence and vertical velocity 
    491       !    2.6.1 ocean depth at u- v- points 
    492       ! hu_crs(:,:) = ssh_un_crs(:,:)                ! now ocean depth (at u- and v-points) 
    493       ! hv_crs(:,:) = ssh_vn_crs(:,:) 
    494       !                                            ! now masked inverse of the ocean depth (at u- and v-points) 
    495       ! jes. not hur_crs, huv_crs are needed. remove? 
    496       !hur_crs(:,:) = umask_crs(:,:,1) / ( hu_crs(:,:) + 1._wp - umask_crs(:,:,1) ) 
    497       !hvr_crs(:,:) = vmask_crs(:,:,1) / ( hv_crs(:,:) + 1._wp - vmask_crs(:,:,1) ) 
    498       ! 
    499       !    2.6.2. Horizontal divergence 
    500       z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog) 
    501       IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt 
    502       z1_2dt = 0.e0 / (z2dt) 
    503       ! 
    504       WHERE ( crs_volt_wgt /= 0 ) zvolt_wgt(:,:,:) = 1.0/crs_volt_wgt(:,:,:) 
    505       ! 
    506       fse3t_crs(:,:,:) = e3t_crs(:,:,:) 
    507       IF ( lk_vvl ) fse3t_crs(:,:,:) = fse3t_n_crs(:,:,:) 
    508  
    509       zhdivbt(:,:) = 0.0 
     189         &         p_fse3=zfse3v, p_pfield=zs, p_cfield3d=zs_crs ) 
     190      CALL crs_iom_put( "voces_crs" , pv_r3d=zs_crs )   ! vS 
     191      ! 
     192 
     193 
     194      !  Horizontal divergence ( following OPA_SRC/DYN/divcur.F90 )  
    510195      DO jk = 1, jpkm1 
    511 !!cc1         z2dcrsu(:,:) = 0.0; z2dcrsv(:,:) = 0.0    
    512196         DO ji = 2, jpi_crsm1 
    513197            DO jj = 2, jpj_crsm1 
    514  
    515                ! Horizontal divergence ( following OPA_SRC/DYN/divcur.F90 )  
    516                ! with partial steps and/or variable layer thicknesses for W 
    517 !!cc1               z2dcrsu(ji,jj) =  ( zum(ji,  jj,jk) * crs_surfu_wgt(ji,  jj,jk)  ) & 
    518                z2dcrsu =  ( zum(ji,  jj,jk) * crs_surfu_wgt(ji,  jj,jk)  ) & 
    519                   &            - ( zum(ji-1,jj,jk) * crs_surfu_wgt(ji-1,jj,jk)  ) 
    520 !!cc1               z2dcrsv(ji,jj) =  ( zvm(ji,jj,  jk) * crs_surfv_wgt(ji,jj  ,jk)  ) & 
    521                z2dcrsv =  ( zvm(ji,jj,  jk) * crs_surfv_wgt(ji,jj  ,jk)  ) & 
    522                   &            - ( zvm(ji,jj-1,jk) * crs_surfv_wgt(ji,jj-1,jk)  ) 
    523  
    524 !!cc1               zhdiv(ji,jj,jk) = ( z2dcrsu(ji,jj) + z2dcrsv(ji,jj) ) * zvolt_wgt(ji,jj,jk)  
    525                zhdiv(ji,jj,jk) = ( z2dcrsu + z2dcrsv ) * zvolt_wgt(ji,jj,jk)  
    526   
     198               IF( tmask_crs(ji,jj,jk ) > 0 ) THEN 
     199                   z2dcrsu =  ( un_crs(ji  ,jj  ,jk) * crs_surfu_wgt(ji  ,jj  ,jk) ) & 
     200                      &     - ( un_crs(ji-1,jj  ,jk) * crs_surfu_wgt(ji-1,jj  ,jk) ) 
     201                   z2dcrsv =  ( vn_crs(ji  ,jj  ,jk) * crs_surfv_wgt(ji  ,jj  ,jk) ) & 
     202                      &     - ( vn_crs(ji  ,jj-1,jk) * crs_surfv_wgt(ji  ,jj-1,jk) ) 
     203                   ! 
     204                   hdivn_crs(ji,jj,jk) = ( z2dcrsu + z2dcrsv ) / crs_volt_wgt(ji,jj,jk)  
     205               ENDIF 
    527206            ENDDO 
    528207         ENDDO 
    529  
    530          ! Horizontal barotropic divergence for ssh ( following OPA_SRC/DYN/sshwzv.F90 ) 
    531  
    532          zhdivbt(:,:) =  zhdivbt(:,:) + fse3t_crs(:,:,jk) * zhdiv(:,:,jk)  
    533  
    534208      ENDDO 
    535  
    536       CALL crs_lbc_lnk( 'T', 1.0, pt3d1=zhdiv ) 
    537       CALL crs_lbc_lnk( 'T', 1.0, pt3d1=zhdivbt ) 
    538       ! accumulate  
    539       zhdivnm(:,:,:) = zhdivnm(:,:,:) + zhdiv 
    540       zhdivbtm(:,:) = zhdivbtm(:,:) + zhdivbt 
    541  
    542 !      !!! Calcul de l'energie cinétique   !cc  !! DECLARE LES VARIABLES 
    543 !       
    544 !      zun2(:,:,:) = un(:,:,:) * un(:,:,:)    ! U^2 
    545 !      zvn2(:,:,:) = vn(:,:,:) * vn(:,:,:)    ! V^2 
    546 !       
    547 !!      ! moyenne sur i de U^2 
    548 !       
    549 !      DO ji = 1, jpiglo-1 
    550 !         zun2(ji,:,:) = 0.5 * (zun2(ji,:,:) + zun2(ji+1,:,:)) 
    551 !      END 
    552 !      ji = jpiglo 
    553 !      zun2(ji,:,:) = 0.5 * zun2(ji,:,:)  
    554 !      uun2(:,:,:) = zun2(:,:,:) 
    555 !       
    556 !       
    557 !         CALL crs_lbc_lnk( cd_type='T', 1.0, pt3d1=uun2) 
    558 !         ! lbcnlk met la ligne jpj = 1 a 0 donc il faut la remettre en ne pas oubliant le cyclique est-ouest  
    559 !         ! a faire un case pour cyclique est-ouest ou non.       
    560 !         uun2(   :   ,1,:) = zun2(    :    ,1,:)    
    561 !          
    562 !         IF ((jperio==4) .OR. (jperio==6)) THEN 
    563 !            uun2(   1   ,1,:) = zun2(jpi_crsm1,1,:) 
    564 !            uun2(jpi_crs,1,:) = zun2(    2    ,1,:) 
    565 !         ENDIF   
    566 ! 
    567 ! 
    568  !     DO jj = 1, jpjglo-1 
    569  !        zvn2(:,jj,:) = 0.5 * (zvn2(:,jj,:) + zvn2(:,jj+1,:)) 
    570  !     END 
    571  !     jj = jpjglo  
    572  !     zvn2(:,jj,:) = 0.5 * zvn2(:,jj,:) 
    573  !     vvn2(:,:,:) = zvn2(:,:,:) 
    574  !      
    575  !     CALL crs_lbc_lnk( cd_type='T', 1.0, pt3d1=vvn2) 
    576  !        ! lbcnlk met la ligne jpj = 1 a 0 donc il faut la remettre en ne pas oubliant le cyclique est-ouest  
    577  !        ! a faire un case pour cyclique est-ouest ou non.       
    578  !        vvn2(   :   ,1,:) = zvn2(    :    ,1,:)    
    579  !         
    580  !        IF ((jperio==4) .OR. (jperio==6)) THEN 
    581  !           vvn2(   1   ,1,:) = zvn2(jpi_crsm1,1,:) 
    582  !           vvn2(jpi_crs,1,:) = zvn2(    2    ,1,:) 
    583  !        ENDIF   
    584  
    585  
    586  
    587        
    588       !    2.6.3. Sea-surface Height  ( See DOM/istate.F90: ssh init; OPA_SRC/DYN/sshwzv.F90: ssh update ) 
    589       !cc 
    590       z2dcrs(:,:) = 0.0 
    591       CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, p_pfield2d=sshn, & 
    592       &                  p_cfield2d=z2dcrs ) 
    593       ssh_crs2(:,:) = z2dcrs(:,:) 
    594   !    WRITE(numout,*) 'test', sshn(:,:) 
    595   !    WRITE(numout,*) 'test', ssh_crs2(:,:) 
    596        
    597      !cc     
    598        
    599       !       set some temp variables 
    600       z2dcrs(:,:) = 0.0; z2dcrsf(:,:) = 0.0 
    601       zcoefu = 0.0; zcoefv = 0.0; zcoeff = 0.0;  
    602       zraur = 1._wp / rau0 
    603       z1_rau0 = 0.5 / rau0 
    604       !       sshn, sshb to accumulate 
    605       DO jj = 2, jpj_crsm1                                      ! leap-frog on ssh_e 
    606          DO ji = 2, jpi_crsm1 
    607             z2dcrs(ji,jj) = ( ssh_b_crs(ji,jj) - z2dt * ( zraur * ( del_emp_crs(ji,jj) ) + zhdivbt(ji,jj) ) ) * tmask_crs(ji,jj,1)  
    608          END DO 
    609       END DO 
    610       CALL crs_lbc_lnk('T',1.0,pt2d=z2dcrs)    
    611       zsshm(:,:)  = zsshm(:,:)  + z2dcrs(:,:) 
    612       zsshbm(:,:) = zsshbm(:,:) + ssh_b_crs(:,:) 
    613  
    614       ! 
    615       !       sshu, sshv, sshf, for vvl case. Calculate then accumulate the fields for time-mean 
    616       ! 
    617             DO jj = 1, jpj_crsm1 
    618                DO ji = 1, jpi_crsm1                     
    619                   zcoefu = 0.5  * umask_crs(ji,jj,1) / ( e1u_crs(ji,jj) * e2u_crs(ji,jj) ) 
    620                   zcoefv = 0.5  * vmask_crs(ji,jj,1) / ( e1v_crs(ji,jj) * e2v_crs(ji,jj) ) 
    621                   zij   =  e1e2t_crs(ji  ,jj  ) * z2dcrs(ji  ,jj  )  
    622                   zip1j =  e1e2t_crs(ji+1,jj  ) * z2dcrs(ji+1,jj  )  
    623                   zijp1 =  e1e2t_crs(ji  ,jj+1) * z2dcrs(ji  ,jj+1)  
    624                
    625                   ssh_un_crs(ji,jj) = zcoefu * ( zij + zip1j ) 
    626                   ssh_vn_crs(ji,jj) = zcoefv * ( zij + zijp1 ) 
    627          IF( lk_vvl ) THEN 
    628                   zsshub(ji,jj) = zcoefu * ( e1e2t_crs(ji,  jj) * ssh_b_crs(ji  ,jj)   & 
    629                      &                        +  e1e2t_crs(ji+1,jj) * ssh_b_crs(ji+1,jj) ) 
    630                   zsshvb(ji,jj) = zcoefv * ( e1e2t_crs(ji,  jj) * ssh_b_crs(ji,  jj)   & 
    631                      &                        +  e1e2t_crs(ji,jj+1) * ssh_b_crs(ji,jj+1) ) 
    632          ENDIF 
    633                END DO 
    634             END DO 
    635  
    636             DO jj = 1, jpj_crsm1 
    637                DO ji = 1, jpi_crsm1      
    638                   z2dcrsf(ji,jj) = 0.5  * umask_crs(ji,jj,1) * umask_crs(ji,jj+1,1)                   & 
    639                        &               / ( e1f_crs(ji,jj  ) * e2f_crs(ji,jj  ) )                     & 
    640                        &               * ( e1u_crs(ji,jj  ) * e2u_crs(ji,jj  ) * ssh_un_crs(ji,jj  )     & 
    641                        &                 + e1u_crs(ji,jj+1) * e2u_crs(ji,jj+1) * ssh_un_crs(ji,jj+1) ) 
    642                END DO 
    643             END DO 
    644       !     accumulate fields in time 
    645             CALL crs_lbc_lnk( 'U', 1., pt2d=ssh_un_crs )   ! sshun_crs 
    646             CALL crs_lbc_lnk( 'V', 1., pt2d=ssh_vn_crs )   ! sshvn_crs 
    647             CALL crs_lbc_lnk('F', 1., pt2d=z2dcrsf ) 
    648             zsshum(:,:)  = zsshum(:,:)  + ssh_un_crs(:,:) 
    649             zsshvm(:,:)  = zsshvm(:,:)  + ssh_vn_crs(:,:) 
    650             zsshfm(:,:)  = zsshfm(:,:)  + z2dcrsf(:,:) 
    651          IF( lk_vvl ) THEN 
    652             CALL crs_lbc_lnk( 'U', 1., pt2d=zsshub ) 
    653             CALL crs_lbc_lnk( 'V', 1., pt2d=zsshvb ) 
    654             zsshubm(:,:) = zsshubm(:,:) + zsshub(:,:)  
    655             zsshvbm(:,:) = zsshvbm(:,:) + zsshvb(:,:)  
    656          ENDIF 
    657       !       update sshb to sshn 
    658       ssh_b_crs(:,:) = z2dcrs(:,:)  
    659       !       update ssha from sshn 
    660       ssh_a_crs(:,:) = (  z2dcrs(:,:) - z2dt * ( z1_rau0 * ( sum_emp_crs(:,:) ) + zhdivbt(:,:) )  ) * tmask(:,:,1) 
    661       CALL crs_lbc_lnk('T', 1.0, pt2d=ssh_a_crs) 
    662       zssham(:,:) = zssham(:,:) + ssh_a_crs(:,:) 
    663       !       accumulate fields in time and update sshua from sshun 
    664       IF ( lk_vvl ) THEN 
    665          DO jj = 1, jpj_crsm1 
    666             DO ji = 1, jpi_crsm1      
    667                zsshua(ji,jj) = 0.5  * umask_crs(ji,jj,1) / ( e1u_crs(ji  ,jj) * e2u_crs(ji  ,jj) )                   & 
    668                   &                                  * ( e1t_crs(ji  ,jj) * e2t_crs(ji  ,jj) * ssh_a_crs(ji  ,jj)     & 
    669                   &                                    + e1t_crs(ji+1,jj) * e2t_crs(ji+1,jj) * ssh_a_crs(ji+1,jj) ) 
    670                zsshva(ji,jj) = 0.5  * vmask_crs(ji,jj,1) / ( e1v_crs(ji,jj  ) * e2v_crs(ji,jj  ) )                   & 
    671                   &                                  * ( e1t_crs(ji,jj  ) * e2t_crs(ji,jj  ) * ssh_a_crs(ji,jj  )     & 
    672                   &                                    + e1t_crs(ji,jj+1) * e2t_crs(ji,jj+1) * ssh_a_crs(ji,jj+1) ) 
    673             END DO 
    674          END DO 
    675  
    676          CALL crs_lbc_lnk( 'U', 1., pt2d=zsshua ) 
    677          CALL crs_lbc_lnk( 'V', 1., pt2d=zsshva ) 
    678  
    679          zsshuam(:,:) = zsshuam(:,:) + zsshua(:,:)  
    680          zsshvam(:,:) = zsshvam(:,:) + zsshva(:,:)      
    681       ENDIF 
    682       ! 
    683       !    2.5.2.  W-velocity ( See OPA_SRC/DYN/sshwzv.F90 ) 
    684       IF ( lk_vvl ) THEN 
    685          !                                 !==  mu computation  ==! 
    686          zee_t(:,:) = e3t_crs(:,:,1)       ! Lower bound : thickness of the first model level 
    687          zee_f(:,:) = e3f_crs(:,:,1)       
    688          zee_u(:,:) = e3u_crs(:,:,1)       
    689          zee_v(:,:) = e3v_crs(:,:,1)       
    690          ! 
    691          DO jk = 2, jpkm1                  ! Sum of the masked vertical scale factors 
    692             zee_t(:,:) = zee_t(:,:) + e3t_crs(:,:,jk) * tmask_crs(:,:,jk) 
    693             zee_u(:,:) = zee_u(:,:) + e3u_crs(:,:,jk) * umask_crs(:,:,jk) 
    694             zee_v(:,:) = zee_v(:,:) + e3v_crs(:,:,jk) * vmask_crs(:,:,jk) 
    695             DO jj = 1, jpj_crsm1                      ! f-point : fmask=shlat at coasts, use the product of umask 
    696                zee_f(:,jj) = zee_f(:,jj) + e3f_crs(:,jj,jk) *  umask_crs(:,jj,jk) * umask_crs(:,jj+1,jk) 
    697             END DO 
    698          END DO 
    699          !                                 ! Compute and mask the inverse of the local depth at T, U, V and F points 
    700          zee_t(:,:) = 1._wp / zee_t(:,:) * tmask_crs(:,:,1) 
    701          zee_u(:,:) = 1._wp / zee_u(:,:) * umask_crs(:,:,1) 
    702          zee_v(:,:) = 1._wp / zee_v(:,:) * vmask_crs(:,:,1) 
    703          DO jj = 1, jpj_crsm1                               ! f-point case fmask cannot be used 
    704             zee_f(:,jj) = 1._wp / zee_f(:,jj) * umask_crs(:,jj,1) * umask_crs(:,jj+1,1) 
    705          END DO 
    706          CALL crs_lbc_lnk( 'F', 1., pt2d=zee_f )                 ! lateral boundary condition on ee_f 
    707          ! 
    708          DO jk = 1, jpk                    ! mu coefficients 
    709             zmut(:,:,jk) = zee_t(:,:) * tmask_crs(:,:,jk)     ! T-point at T levels 
    710          END DO 
    711          DO jk = 1, jpk                                 ! F-point : fmask=shlat at coasts, use the product of umask 
    712             DO jj = 1, jpj_crsm1 
    713                zmuf(:,jj,jk) = zee_f(:,jj) * umask_crs(:,jj,jk) * umask_crs(:,jj+1,jk)   ! at T levels 
    714             END DO 
    715             zmuf(:,jpj,jk) = 0._wp 
    716          END DO 
    717          CALL crs_lbc_lnk( 'F', 1., pt3d1=zmuf )                   ! lateral boundary condition 
    718  
    719          DO jk = 1, jpk                    ! mu coefficients 
    720             fse3t_crs(:,:,jk)   = e3t_crs(:,:,jk) 
    721             fse3t_n_crs(:,:,jk) = e3t_crs(:,:,jk) + z2dcrs(:,:)*zmut(:,:,jk) 
    722             fse3t_b_crs(:,:,jk) = fse3t_b_crs(:,:,jk) 
    723             fse3t_a_crs(:,:,jk) = e3t_crs(:,:,jk) + ssh_a_crs(:,:)*zmut(:,:,jk) 
    724          ENDDO 
    725 ! jes. 30 aug 2012.  Need to add here fse3u, fse3v, fse3f - n,b,a 
    726       ELSE 
    727          fse3t_crs = e3t_crs 
    728 !         fse3w_crs = e3w_crs 
    729 !         fse3u_crs = e3u_crs 
    730 !         fse3v_crs = e3v_crs 
    731       ENDIF 
    732       ! 
    733       z3dcrs(:,:,:) = 0.0       
    734       IF ( lk_vvl ) THEN 
    735          DO jk = jpkm1, 1, -1                             ! integrate from the bottom the hor. divergence 
    736             z3dcrs(:,:,jk) = z3dcrs(:,:,jk+1) -   fse3t_n_crs(:,:,jk) * zhdiv(:,:,jk)        & 
    737                &                      - ( fse3t_a_crs(:,:,jk) - fse3t_b_crs(:,:,jk) )    & 
    738                &                         * tmask_crs(:,:,jk) * z1_2dt 
    739          END DO 
    740       ELSE 
    741          DO jk = jpkm1, 1, -1                             ! calculate from the bottom-up 
    742             z3dcrs(:,:,jk) = z3dcrs(:,:,jk+1) - fse3t_crs(:,:,jk) * zhdiv(:,:,jk)         
    743          ENDDO 
    744       ENDIF 
    745       CALL crs_lbc_lnk( 'T', 1.0, pt3d1=z3dcrs ) 
    746       zwm(:,:,:) = zwm(:,:,:) + z3dcrs(:,:,:)     
    747  
    748       ! 2.5.3 Slopes  for slopes look at OPA_SRC/DYN/dynldf_iso.F90 and OPA_SRC/LDF/ldfslp.F90 
    749  
    750  
    751       ! 3. Output of some dynamics and tracer fields and selected fields 
    752       !   * Within the nn_fcrs frequence for iom  
    753       ! -------------------------------------------------------- 
    754       ! 
    755       CALL iom_setkt( kt + nn_fcrs - 1 )    ! in crs, iom_put is called every nn_fcrs time step(s) 
    756   
    757       IF( MOD( kt-1, nn_fcrs ) == 0 ) THEN 
    758  
    759          ! It is time to dump some fields 
    760          ! compute average by scaling the sum by spatial dimension then time 
    761          ! OR by time only 
    762  
    763          ! 3.1  Ocean fields 
    764          !    3.1.1 Weight for time average 
    765          zdtj=rdt/86400.   ! time step in days 
    766          WRITE(clmode,'(f5.1,a)' ) zdtj,' days average' 
    767          WRITE(numout,*) 'crsdiawri.', clmode 
    768          zrtsct = 1.0/REAL(itsct, wp) 
    769  
    770          ! 
    771          !    3.1.2 Weights for spatial averages 
    772          WHERE ( zvolw /= 0 ) zvolw(:,:,:) = 1.0/zvolw(:,:,:) 
    773          WHERE ( zsurv /= 0 ) zsurv(:,:,:) = 1.0/zsurv(:,:,:) 
    774          WHERE ( zsuru /= 0 ) zsuru(:,:,:) = 1.0/zsuru(:,:,:) 
    775          WHERE ( zsurw /= 0 ) zsurw(:,:,:) = 1.0/zsurw(:,:,:) 
    776          WHERE ( zsurtaux /= 0 ) zsurtaux(:,:) = 1.0 / zsurtaux(:,:) 
    777          WHERE ( zsurtauy /= 0 ) zsurtauy(:,:) = 1.0 / zsurtauy(:,:) 
    778  
    779          !    3.1.3 Temperature 
    780          tsn_crs(:,:,:,1) = ztsnm(:,:,:,1) * zrtsct ! volume-weighted- , time- mean   
    781          z2dcrs(:,:) = tsn_crs(:,:,1,1) 
    782          CALL crs_iom_put( "toce_crs"   , pv_r3d=tsn_crs(:,:,:,1)  )    ! temperature 
    783          CALL crs_iom_put( "sst_crs"    , pv_r2d=z2dcrs            )    ! sst 
    784       
    785          ! 
    786          !    3.1.4 Salinity 
    787          tsn_crs(:,:,:,2) = ztsnm(:,:,:,2) * zrtsct ! volume-weighted- , time- mean  
    788          z2dcrs(:,:) = tsn_crs(:,:,1,2) 
    789          CALL crs_iom_put( "soce_crs"   , pv_r3d=tsn_crs(:,:,:,2)   )   !  salinity 
    790          CALL crs_iom_put( "sss_crs"    , pv_r2d=z2dcrs             )   !  sss 
    791          ! 
    792  
    793          !    3.1.5 U-velocity 
    794          un_crs(:,:,:) = zum(:,:,:) * zrtsct ! area-weighted- , time- mean  
    795          CALL crs_iom_put( "uoce_crs"   , pv_r3d=un_crs             )   ! i-current  
    796          ut_crs(:,:,:) = zutm(:,:,:) * zsuru(:,:,:) ! area-weighted- , time- mean       
    797          CALL crs_iom_put( "uocet_crs"     , pv_r3d=ut_crs             )   ! uT 
    798          us_crs(:,:,:) = zusm(:,:,:) * zsuru(:,:,:) ! area-weighted- , time- mean       
    799          CALL crs_iom_put( "uoces_crs"     , pv_r3d=us_crs             )   ! uS 
    800                
    801          !    3.1.6 V-velocity 
    802          vn_crs(:,:,:) = zvm(:,:,:) * zrtsct ! area-weighted- , time- mean  
    803          CALL crs_iom_put( "voce_crs"   , pv_r3d=vn_crs                   )    ! j-current 
    804          vt_crs(:,:,:) = zvtm(:,:,:) * zsurv(:,:,:) ! area-weighted- , time- mean       
    805          CALL crs_iom_put( "vocet_crs"     , pv_r3d=vt_crs             )   ! vT 
    806          vs_crs(:,:,:) = zvsm(:,:,:) * zsurv(:,:,:) ! area-weighted- , time- mean       
    807          CALL crs_iom_put( "voces_crs"     , pv_r3d=vs_crs             )   ! vS 
    808  
    809          !    3.1.7 W-velocity 
    810          wn_crs(:,:,:) = zwm(:,:,:) * zrtsct ! area-weighted- , time- mean  
    811          CALL crs_iom_put( "woce_crs"   , pv_r3d=wn_crs                   )    ! W-velocity 
    812          CALL crs_iom_put( "woce2_crs"   , pv_r3d=zw                      )    ! cc 
    813  
    814          !    3.1.8 Horizontal divergence 
    815          hdivn_crs(:,:,:) = zhdivnm(:,:,:) * zrtsct 
    816       !   CALL crs_iom_put( "hdivn_crs"   , pv_r3d=hdivn_crs            )   ! hdiv 
    817          CALL crs_iom_put( "hdiv_crs"   , pv_r3d=hdivn_crs            )   
    818  
    819          hdivbt_crs(:,:) = zhdivbtm(:,:) * zrtsct 
    820  
    821          !    3.1.9 Sea-surface Height 
    822          sshn_crs(:,:) = zsshm(:,:)  * zrtsct 
    823          sshb_crs(:,:) = zsshbm(:,:) * zrtsct 
    824          ssha_crs(:,:) = zssham(:,:) * zrtsct 
    825  
    826          IF( lk_vvl ) THEN                    ! now Sea SSH at u-, v-, f-points (vvl case only) 
    827             sshun_crs(:,:) = zsshum(:,:)  * zrtsct                              ! sshun  
    828             sshub_crs(:,:) = zsshubm(:,:) * zrtsct                              ! sshub  
    829             sshua_crs(:,:) = zsshuam(:,:) * zrtsct                              ! sshua  
    830             sshvn_crs(:,:) = zsshvm(:,:)  * zrtsct                              ! sshvn  
    831             sshvb_crs(:,:) = zsshvbm(:,:) * zrtsct                              ! sshvb  
    832             sshva_crs(:,:) = zsshvam(:,:) * zrtsct                              ! sshva  
    833             sshfn_crs(:,:) = zsshfm(:,:) * zrtsct                               ! sshfn  
    834          ENDIF 
    835          CALL crs_iom_put( "ssh_crs"   , pv_r2d=sshn_crs                  )   ! ssh output  
    836          CALL crs_iom_put( "ssh2_crs"  , pv_r2d=ssh_crs2                  )   !cc 
    837  
    838          !    3.1.10 Potential density 
    839          rhop_crs(:,:,:) = zrhopm(:,:,:) * zrtsct 
    840          rhd_crs(:,:,:)  =  zrhdm(:,:,:) * zrtsct  
    841  
    842          !    3.1.10 Vertical eddy diffusivity 
    843          avt_crs(:,:,:) = zavtm(:,:,:) * zsurw(:,:,:)  
    844  
    845          !    3.1.11  Mixed-layer diagnostics 
    846          !            Recalculated following OPA_SRC/ZDF/zdfmxl.F90 
    847          !    
    848          ! w-level of the mixing and mixed layers 
    849          nmln_crs(:,:) = mbkt_crs(:,:) + 1        ! Initialization to the number of w ocean point 
    850          imld(:,:) = mbkt_crs(:,:) + 1 
    851          DO jk = jpkm1, nlb10, -1         ! from the bottom to nlb10 
    852             DO jj = 1, jpj_crs 
    853                DO ji = 1, jpi_crs 
    854                   IF( rhop_crs(ji,jj,jk) > rhop_crs(ji,jj,nla10) + zrho_c )   nmln_crs(ji,jj) = jk      ! Mixed layer 
    855                   IF(  avt_crs(ji,jj,jk) < zavt_c                     )           imld(ji,jj) = jk      ! Turbocline 
    856                END DO 
    857             END DO 
    858          END DO 
    859          ! depth of the mixing and mixed layers 
    860          DO jj = 1, jpj_crs 
    861             DO ji = 1, jpi_crs 
    862                iiki = imld(ji,jj) 
    863                iikn = nmln_crs(ji,jj) 
    864                ! jes. TODO.  update here gdep for vvl case 
    865                hmld_crs (ji,jj) = gdepw_crs(ji,jj,iiki  ) * tmask_crs(ji,jj,1)    ! Turbocline depth 
    866                hmlp_crs (ji,jj) = gdepw_crs(ji,jj,iikn  ) * tmask_crs(ji,jj,1)    ! Mixed layer depth 
    867                hmlpt_crs(ji,jj) = gdept_crs(ji,jj,iikn-1)                         ! depth of the last T-point inside the mixed layer 
    868             END DO 
    869          END DO 
    870  
    871  
    872          ! 3.2. Surface fluxes 
    873          !    3.2.1. Evaporation minus Precipitation 
    874          emp_crs(:,:)   = zempm(:,:)   * zrtsct                      
    875          emp_b_crs(:,:) = zemp_bm(:,:) * zrtsct                      
    876          emps_crs(:,:)  = zempsm(:,:)  * zrtsct                      
    877          !    3.2.2. Wind stresses, Scalar Wind 
    878          utau_crs(:,:) = zutaum(:,:) * zsurtaux(:,:) 
    879          vtau_crs(:,:) = zvtaum(:,:) * zsurtauy(:,:) 
    880          wndm_crs(:,:) = zwndmm(:,:) * zsurw(:,:,1) 
    881          !    3.2.3. Shortwave radiative flux 
    882          qsr_crs(:,:) = zqsrm(:,:) * zsurw(:,:,1) 
    883          !    3.2.4. Runoff freshwater flux 
    884          rnf_crs(:,:) = zrnfm(:,:) * zsurw(:,:,1) 
    885          !    3.2.4. Depth of freshwater flux input 
    886          h_rnf_crs(:,:) = zh_rnfm(:,:) * zrtsct 
    887          !    3.2.5. Ice fraction 
    888          fr_i_crs(:,:) = zfr_im(:,:) * zsurw(:,:,1) 
    889  
    890  
    891          ! Reset work arrays and counter to 0 after writing 
    892          ! 
    893          itsct = 0 
    894          ! 
    895          zum(:,:,:)   = 0._wp 
    896          zvm(:,:,:)   = 0._wp 
    897          ztsnm(:,:,:,:) = 0._wp 
    898          zwm(:,:,:)   = 0._wp 
    899          zutm(:,:,:)  = 0._wp 
    900          zusm(:,:,:)  = 0._wp 
    901          zvtm(:,:,:)  = 0._wp 
    902          zvsm(:,:,:)  = 0._wp 
    903          zwtm(:,:,:)  = 0._wp 
    904          zwsm(:,:,:)  = 0._wp 
    905  
    906          zvolw(:,:,:) = 0._wp 
    907          zsuru(:,:,:) = 0._wp 
    908          zsurv(:,:,:) = 0._wp 
    909          zsurw(:,:,:) = 0._wp    
    910          zsurtaux(:,:) = 0._wp 
    911          zsurtauy(:,:) = 0._wp 
    912  
    913          zsshm(:,:) = 0._wp 
    914          zsshum(:,:) = 0._wp 
    915          zsshvm(:,:) = 0._wp 
    916          zsshfm(:,:) = 0._wp 
    917          zsshbm(:,:) = 0._wp 
    918          zsshubm(:,:) = 0._wp 
    919          zsshvbm(:,:) = 0._wp 
    920          zssham(:,:) = 0._wp 
    921          zsshuam(:,:) = 0._wp 
    922          zsshvam(:,:) = 0._wp 
    923          zwndmm(:,:) = 0._wp 
    924          zutaum(:,:) = 0._wp 
    925          zvtaum(:,:) = 0._wp 
    926          zqsrm(:,:) = 0._wp 
    927          zrnfm(:,:) = 0._wp 
    928          zh_rnfm(:,:) = 0._wp 
    929          zfr_im(:,:) = 0._wp 
    930  
    931       ENDIF 
    932  
    933       CALL iom_setkt( kt )                  ! iom_put outside of crs is called at every time step 
    934  
    935  
    936  
    937       ! 5. Clean-up 
    938       ! -------------------------------------------------------- 
    939       ! 
    940       CALL wrk_dealloc( jpi , jpj , z2d, ze1e2u, ze1e2v ) 
    941       CALL wrk_dealloc( jpi , jpj , jpk , z3d, z3d1 ) 
    942       CALL wrk_dealloc( jpi , jpj , jpk , zfse3u, zfse3v ) 
    943       CALL wrk_dealloc( jpi , jpj , jpk , zfse3t, zfse3w ) 
    944       CALL wrk_dealloc( jpi_crs , jpj_crs , jpk , z3dcrs, zw ) !cc 
    945       CALL wrk_dealloc( jpi_crs , jpj_crs , z2dcrs, ssh_crs2 ) 
    946  
    947 !!cc1      DEALLOCATE( z2dcrsu, z2dcrsv, z2dcrsf, zhdivbt ) ! probleme de malloc au 65 eme pas de temps 
    948       DEALLOCATE( z2dcrsf, zhdivbt ) ! probleme de malloc au 65 eme pas de temps 
    949       DEALLOCATE( zsshub, zsshua, zsshvb, zsshva ) 
    950       DEALLOCATE( zee_t, zee_f, zee_u, zee_v ) 
    951       DEALLOCATE( zmut, zmuf ) 
    952       DEALLOCATE( zcrs_surtaux_wgt, zcrs_surtauy_wgt ) 
    953       DEALLOCATE( zhdiv, zvolt_wgt)  
    954       ! 
    955                          
     209      CALL crs_lbc_lnk( 'T', 1.0,    pt3d1=hdivn_crs ) 
     210      CALL crs_iom_put( "hdiv_crs", pv_r3d=hdivn_crs )   
     211 
     212       
     213      !  Sea-surface Height  
     214      CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
     215           &       p_pfield2d=sshn, p_cfield2d=sshn_crs ) 
     216      CALL crs_iom_put( "ssh_crs" , pv_r2d=sshn_crs  )   ! ssh output  
     217 
     218      !  W-velocity 
     219      CALL crsfun( p_e1e2t=e1e2t, cd_type='T', cd_op='ARE', p_cmask=tmask_crs, p_ptmask=tmask, & 
     220         &         p_pfield3d_2=wn, p_cfield3d=wn_crs ) 
     221      CALL crs_iom_put( "woce_crs"  , pv_r3d=wn_crs  )   ! i-current  
     222       
     223 
     224 
     225      !  Clean-up 
     226       
     227      CALL wrk_dealloc( jpi, jpj, jpk, zfse3t, zfse3w ) 
     228      CALL wrk_dealloc( jpi, jpj, jpk, zfse3u, zfse3v ) 
     229      CALL wrk_dealloc( jpi, jpj, jpk, zt, zs       ) 
     230      ! 
     231      CALL wrk_dealloc( jpi_crs, jpj_crs, jpk, zt_crs, zs_crs ) 
     232      ! 
    956233      IF( nn_timing == 1 )   CALL timing_stop('crs_dia_wri') 
    957234      ! 
Note: See TracChangeset for help on using the changeset viewer.