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 7309 for branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2016-11-22T18:43:11+01:00 (8 years ago)
Author:
clem
Message:

first implementations

Location:
branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC
Files:
20 edited

Legend:

Unmodified
Added
Removed

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90

    r5836 r7309  
    2727#elif defined key_lim3 
    2828   USE ice             ! LIM_3 ice variables 
    29    USE dom_ice         ! sea-ice domain 
    3029   USE limvar 
     30   USE limctl 
    3131#endif  
    3232   USE par_oce         ! ocean parameters 
    … …  
    8282      ! 
    8383#if defined key_lim3 
    84       CALL lim_var_zapsmall 
    85       CALL lim_var_agg(1) 
     84                        CALL lim_var_zapsmall 
     85                        CALL lim_var_agg(1) 
     86      IF( ln_limctl )   CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 
    8687#endif 
    8788      ! 
    … …  
    121122      ! 
    122123#if defined key_lim2 
    123       DO jb = 1, idx%nblen(jgrd) 
     124      DO jb = 1, idx%nblenrim(jgrd) 
    124125         ji    = idx%nbi(jb,jgrd) 
    125126         jj    = idx%nbj(jb,jgrd) 
    … …  
    141142 
    142143      DO jl = 1, jpl 
    143          DO jb = 1, idx%nblen(jgrd) 
     144         DO jb = 1, idx%nblenrim(jgrd) 
    144145            ji    = idx%nbi(jb,jgrd) 
    145146            jj    = idx%nbj(jb,jgrd) 
    … …  
    177178 
    178179      DO jl = 1, jpl 
    179          DO jb = 1, idx%nblen(jgrd) 
     180         DO jb = 1, idx%nblenrim(jgrd) 
    180181            ji    = idx%nbi(jb,jgrd) 
    181182            jj    = idx%nbj(jb,jgrd) 
    … …  
    236237            END SELECT 
    237238            ! 
    238             IF( nn_icesal == 1 ) THEN     ! constant salinity : overwrite rn_ice_sal 
     239            IF( nn_icesal == 1 ) THEN     ! constant salinity : overwrite rn_icesal 
    239240               sm_i(ji,jj  ,jl) = rn_icesal 
    240241               s_i (ji,jj,:,jl) = rn_icesal 
    … …  
    325326            CASE ( 'U' )   
    326327               jgrd = 2      ! u velocity 
    327                DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 
     328               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    328329                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    329330                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 
    … …  
    352353            CASE ( 'V' ) 
    353354               jgrd = 3      ! v velocity 
    354                DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 
     355               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    355356                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    356357                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r6140 r7309  
    3838   PUBLIC   dia_hsb        ! routine called by step.F90 
    3939   PUBLIC   dia_hsb_init   ! routine called by nemogcm.F90 
    40    PUBLIC   dia_hsb_rst    ! routine called by step.F90 
    4140 
    4241   LOGICAL, PUBLIC ::   ln_diahsb   !: check the heat and salt budgets 
    … …  
    8685      !!--------------------------------------------------------------------------- 
    8786      IF( nn_timing == 1 )   CALL timing_start('dia_hsb')       
     87      ! 
    8888      CALL wrk_alloc( jpi,jpj,   z2d0, z2d1 ) 
    8989      ! 
    … …  
    171171      ENDDO 
    172172 
    173       ! Substract forcing from heat content, salt content and volume variations 
     173      ! ------------------------ ! 
     174      ! 3 -  Drifts              ! 
     175      ! ------------------------ ! 
    174176      zdiff_v1 = zdiff_v1 - frc_v 
    175177      IF( .NOT.ln_linssh )   zdiff_v2 = zdiff_v2 - frc_v 
    … …  
    184186 
    185187      ! ----------------------- ! 
    186       ! 3 - Diagnostics writing ! 
     188      ! 4 - Diagnostics writing ! 
    187189      ! ----------------------- ! 
    188190      zvol_tot = 0._wp                    ! total ocean volume (calculated with scale factors) 
    … …  
    197199!!gm end 
    198200 
    199       IF( ln_linssh ) THEN 
    200         CALL iom_put( 'bgtemper' , zdiff_hc1 / zvol_tot)              ! Heat content variation (C)  
    201         CALL iom_put( 'bgsaline' , zdiff_sc1 / zvol_tot)              ! Salt content variation (psu) 
    202         CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * rau0 * rcp )  ! Heat content variation (1.e20 J)  
    203         CALL iom_put( 'bgsaltco' , zdiff_sc1 * 1.e-9   )              ! Salt content variation (psu*km3) 
    204         CALL iom_put( 'bgvolssh' , zdiff_v1  * 1.e-9   )              ! volume ssh variation (km3)   
    205         CALL iom_put( 'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
    206         CALL iom_put( 'bgfrctem' , frc_t / zvol_tot    )              ! hc  - surface forcing (C)  
    207         CALL iom_put( 'bgfrcsal' , frc_s / zvol_tot    )              ! sc  - surface forcing (psu)  
     201      CALL iom_put(   'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
     202      CALL iom_put(   'bgfrctem' , frc_t    * rau0 * rcp * 1.e-20 )   ! hc  - surface forcing (1.e20 J)  
     203      CALL iom_put(   'bgfrchfx' , frc_t    * rau0 * rcp /  &         ! hc  - surface forcing (W/m2)  
     204         &                       ( surf_tot * kt * rdt )        ) 
     205      CALL iom_put(   'bgfrcsal' , frc_s    * 1.e-9    )              ! sc  - surface forcing (psu*km3)  
     206 
     207      IF( .NOT. ln_linssh ) THEN 
     208        CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot )              ! Temperature drift     (C)  
     209        CALL iom_put( 'bgsaline' , zdiff_sc / zvol_tot )              ! Salinity    drift     (pss) 
     210        CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * rau0 * rcp )   ! Heat content drift    (1.e20 J)  
     211        CALL iom_put( 'bgheatfx' , zdiff_hc * rau0 * rcp /  &         ! Heat flux drift       (W/m2)  
     212           &                       ( surf_tot * kt * rdt )        ) 
     213        CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9    )              ! Salt content drift    (psu*km3) 
     214        CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh drift      (km3)   
     215        CALL iom_put( 'bgvole3t' , zdiff_v2 * 1.e-9    )              ! volume e3t drift      (km3)   
     216      ELSE 
     217        CALL iom_put( 'bgtemper' , zdiff_hc1 / zvol_tot)              ! Heat content drift    (C)  
     218        CALL iom_put( 'bgsaline' , zdiff_sc1 / zvol_tot)              ! Salt content drift    (pss) 
     219        CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * rau0 * rcp )  ! Heat content drift    (1.e20 J)  
     220        CALL iom_put( 'bgheatfx' , zdiff_hc1 * rau0 * rcp /  &        ! Heat flux drift       (W/m2)  
     221           &                       ( surf_tot * kt * rdt )         ) 
     222        CALL iom_put( 'bgsaltco' , zdiff_sc1 * 1.e-9    )             ! Salt content drift    (psu*km3) 
     223        CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh drift      (km3)   
    208224        CALL iom_put( 'bgmistem' , zerr_hc1 / zvol_tot )              ! hc  - error due to free surface (C) 
    209225        CALL iom_put( 'bgmissal' , zerr_sc1 / zvol_tot )              ! sc  - error due to free surface (psu) 
    210       ELSE 
    211         CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot )              ! Temperature variation (C)  
    212         CALL iom_put( 'bgsaline' , zdiff_sc / zvol_tot )              ! Salinity    variation (psu) 
    213         CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * rau0 * rcp )   ! Heat content variation (1.e20 J)  
    214         CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9    )              ! Salt content variation (psu*km3) 
    215         CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh variation (km3)   
    216         CALL iom_put( 'bgvole3t' , zdiff_v2 * 1.e-9    )              ! volume e3t variation (km3)   
    217         CALL iom_put( 'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
    218         CALL iom_put( 'bgfrctem' , frc_t / zvol_tot    )              ! hc  - surface forcing (C)  
    219         CALL iom_put( 'bgfrcsal' , frc_s / zvol_tot    )              ! sc  - surface forcing (psu)  
    220226      ENDIF 
    221227      ! 
    … …  
    231237   SUBROUTINE dia_hsb_rst( kt, cdrw ) 
    232238     !!--------------------------------------------------------------------- 
    233      !!                   ***  ROUTINE limdia_rst  *** 
     239     !!                   ***  ROUTINE dia_hsb_rst  *** 
    234240     !!                      
    235241     !! ** Purpose :   Read or write DIA file in restart file 
    … …  
    241247     ! 
    242248     INTEGER ::   ji, jj, jk   ! dummy loop indices 
    243      INTEGER ::   id1          ! local integers 
    244249     !!---------------------------------------------------------------------- 
    245250     ! 
    246251     IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise  
    247252        IF( ln_rstart ) THEN                   !* Read the restart file 
    248            !id1 = iom_varid( numror, 'frc_vol'  , ldstop = .FALSE. ) 
    249253           ! 
    250254           IF(lwp) WRITE(numout,*) '~~~~~~~' 
    … …  
    259263           ENDIF 
    260264           CALL iom_get( numror, jpdom_autoglo, 'surf_ini', surf_ini ) ! ice sheet coupling 
    261            CALL iom_get( numror, jpdom_autoglo, 'ssh_ini', ssh_ini ) 
    262            CALL iom_get( numror, jpdom_autoglo, 'e3t_ini', e3t_ini ) 
    263            CALL iom_get( numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini ) 
    264            CALL iom_get( numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini ) 
     265           CALL iom_get( numror, jpdom_autoglo, 'ssh_ini', ssh_ini(:,:) ) 
     266           CALL iom_get( numror, jpdom_autoglo, 'e3t_ini', e3t_ini(:,:,:) ) 
     267           CALL iom_get( numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini(:,:,:) ) 
     268           CALL iom_get( numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini(:,:,:) ) 
    265269           IF( ln_linssh ) THEN 
    266               CALL iom_get( numror, jpdom_autoglo, 'ssh_hc_loc_ini', ssh_hc_loc_ini ) 
    267               CALL iom_get( numror, jpdom_autoglo, 'ssh_sc_loc_ini', ssh_sc_loc_ini ) 
     270              CALL iom_get( numror, jpdom_autoglo, 'ssh_hc_loc_ini', ssh_hc_loc_ini(:,:) ) 
     271              CALL iom_get( numror, jpdom_autoglo, 'ssh_sc_loc_ini', ssh_sc_loc_ini(:,:) ) 
    268272           ENDIF 
    269273       ELSE 
    … …  
    313317        ENDIF 
    314318        CALL iom_rstput( kt, nitrst, numrow, 'surf_ini', surf_ini )      ! ice sheet coupling 
    315         CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini', ssh_ini ) 
    316         CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini', e3t_ini ) 
    317         CALL iom_rstput( kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini ) 
    318         CALL iom_rstput( kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini ) 
     319        CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini', ssh_ini(:,:) ) 
     320        CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini', e3t_ini(:,:,:) ) 
     321        CALL iom_rstput( kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini(:,:,:) ) 
     322        CALL iom_rstput( kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini(:,:,:) ) 
    319323        IF( ln_linssh ) THEN 
    320            CALL iom_rstput( kt, nitrst, numrow, 'ssh_hc_loc_ini', ssh_hc_loc_ini ) 
    321            CALL iom_rstput( kt, nitrst, numrow, 'ssh_sc_loc_ini', ssh_sc_loc_ini ) 
     324           CALL iom_rstput( kt, nitrst, numrow, 'ssh_hc_loc_ini', ssh_hc_loc_ini(:,:) ) 
     325           CALL iom_rstput( kt, nitrst, numrow, 'ssh_sc_loc_ini', ssh_sc_loc_ini(:,:) ) 
    322326        ENDIF 
    323327        ! 
    … …  
    339343      !!             - Compute coefficients for conversion 
    340344      !!--------------------------------------------------------------------------- 
    341       INTEGER ::   jk       ! dummy loop indice 
    342345      INTEGER ::   ierror   ! local integer 
    343346      INTEGER ::   ios 
    … …  
    345348      NAMELIST/namhsb/ ln_diahsb 
    346349      !!---------------------------------------------------------------------- 
    347  
    348       IF(lwp) THEN 
    349          WRITE(numout,*) 
    350          WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 
    351          WRITE(numout,*) '~~~~~~~~ ' 
    352       ENDIF 
    353350 
    354351      REWIND( numnam_ref )              ! Namelist namhsb in reference namelist 
    … …  
    361358      IF(lwm) WRITE ( numond, namhsb ) 
    362359 
    363       ! 
    364       IF(lwp) THEN                   ! Control print 
     360      IF(lwp) THEN 
    365361         WRITE(numout,*) 
    366          WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 
    367          WRITE(numout,*) '~~~~~~~~~~~~' 
    368          WRITE(numout,*) '   Namelist namhsb : set hsb parameters' 
    369          WRITE(numout,*) '      Switch for hsb diagnostic (T) or not (F)  ln_diahsb  = ', ln_diahsb 
    370          WRITE(numout,*) 
    371       ENDIF 
    372  
     362         WRITE(numout,*) 'dia_hsb_init' 
     363         WRITE(numout,*) '~~~~~~~~ ' 
     364         WRITE(numout,*) '  check the heat and salt budgets (T) or not (F)       ln_diahsb = ', ln_diahsb 
     365      ENDIF 
     366      ! 
    373367      IF( .NOT. ln_diahsb )   RETURN 
    374368         !      IF( .NOT. lk_mpp_rep ) & 
    … …  
    388382      IF( ln_linssh )   ALLOCATE( ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj),STAT=ierror ) 
    389383      IF( ierror > 0 ) THEN 
    390          CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' )   ;   RETURN 
     384         CALL ctl_stop( 'dia_hsb: unable to allocate ssh_hc_loc_ini' )   ;   RETURN 
    391385      ENDIF 
    392386 
    … …  
    394388      ! 2 - Time independant variables and file opening ! 
    395389      ! ----------------------------------------------- ! 
    396       IF(lwp) WRITE(numout,*) "dia_hsb: heat salt volume budgets activated" 
    397       IF(lwp) WRITE(numout,*) '~~~~~~~' 
    398390      surf(:,:) = e1t(:,:) * e2t(:,:) * tmask_i(:,:)      ! masked surface grid cell area 
    399       surf_tot  = glob_sum( surf(:,:) )                                       ! total ocean surface area 
     391      surf_tot  = glob_sum( surf(:,:) )                   ! total ocean surface area 
    400392 
    401393      IF( lk_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' )          

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90

    r6140 r7309  
    166166   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e1e2f , r1_e1e2f                !: associated metrics at f-point 
    167167   ! 
    168    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ff                              !: coriolis factor                   [1/s] 
    169  
     168   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ff                              !: coriolis factor at F_point [1/s] 
     169   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ff_t                            !: coriolis factor at T-point [1/s] 
    170170   !!---------------------------------------------------------------------- 
    171171   !! vertical coordinate and scale factors 
    … …  
    289289   !!---------------------------------------------------------------------- 
    290290   !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
    291    !! $Id$ 
     291   !! $Id$  
    292292   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    293293   !!---------------------------------------------------------------------- 
    … …  
    309309   INTEGER FUNCTION dom_oce_alloc() 
    310310      !!---------------------------------------------------------------------- 
    311       INTEGER, DIMENSION(13) :: ierr 
     311      INTEGER, DIMENSION(12) :: ierr 
    312312      !!---------------------------------------------------------------------- 
    313313      ierr(:) = 0 
    … …  
    332332         &      e1e2v(jpi,jpj) , r1_e1e2v(jpi,jpj) , e1_e2v(jpi,jpj)                   ,     & 
    333333         &      e1e2f(jpi,jpj) , r1_e1e2f(jpi,jpj)                                     ,     & 
    334          &        ff (jpi,jpj)                                                         , STAT=ierr(3) ) 
     334         &        ff (jpi,jpj) , ff_t    (jpi,jpj)                                     , STAT=ierr(3) ) 
    335335         ! 
    336336      ALLOCATE( gdept_0(jpi,jpj,jpk) , gdepw_0(jpi,jpj,jpk) , gde3w_0(jpi,jpj,jpk) ,     & 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r6981 r7309  
    9292      IF( ln_sco )   CALL dom_stiff             ! Maximum stiffness ratio/hydrostatic consistency 
    9393      ! 
    94       ht_0(:,:) = e3t_0(:,:,1) * tmask(:,:,1)   ! Reference ocean thickness 
    95       hu_0(:,:) = e3u_0(:,:,1) * umask(:,:,1) 
    96       hv_0(:,:) = e3v_0(:,:,1) * vmask(:,:,1) 
    97       DO jk = 2, jpk 
     94      ht_0(:,:) = 0._wp  ! Reference ocean thickness 
     95      hu_0(:,:) = 0._wp 
     96      hv_0(:,:) = 0._wp 
     97      DO jk = 1, jpk 
    9898         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk) 
    9999         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk) 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domcfg.F90

    r6140 r7309  
    5656         IF( jperio == 5 )   WRITE(numout,*) '      jperio= 5, north fold with F-point pivot' 
    5757         IF( jperio == 6 )   WRITE(numout,*) '      jperio= 6, cyclic east-west and north fold with F-point pivot' 
    58       ENDIF 
    59       ! 
    60       IF( jperio <  0 .OR. jperio > 6 )   CALL ctl_stop( 'jperio is out of range' ) 
     58         IF( jperio == 7 )   WRITE(numout,*) '      jperio= 7, cyclic east-west and north-south' 
     59      ENDIF 
     60      ! 
     61      IF( jperio <  0 .OR. jperio > 7 )   CALL ctl_stop( 'jperio is out of range' ) 
    6162      ! 
    6263      CALL dom_glo                   ! global domain versus zoom and/or local domain 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90

    r6140 r7309  
    3838 
    3939   !!---------------------------------------------------------------------- 
    40    !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
     40   !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
    4141   !! $Id$  
    4242   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    321321         e1v(:,:) =  ze1     ;      e2v(:,:) = ze1 
    322322         e1f(:,:) =  ze1     ;      e2f(:,:) = ze1 
    323          ! 
     323 
     324      CASE ( 6 )                   ! clem: f-plane with irregular grid-spacing 
     325 
     326         IF(lwp) WRITE(numout,*) 
     327         IF(lwp) WRITE(numout,*) '          f-plane with irregular grid-spacing (+- 10%)' 
     328         IF(lwp) WRITE(numout,*) '          the max is given by ppe1_m and ppe2_m'  
     329 
     330         ! Position coordinates (in kilometers) 
     331         !                          ========== 
     332         glam0 = 0._wp 
     333         gphi0 = 0._wp 
     334          
     335#if defined key_agrif  
     336         IF( .NOT. Agrif_Root() ) THEN 
     337            glam0  = Agrif_Parent(glam0) + (Agrif_ix())*Agrif_Parent(ppe1_m) * 1.e-5 
     338            gphi0  = Agrif_Parent(gphi0) + (Agrif_iy())*Agrif_Parent(ppe2_m) * 1.e-5 
     339            ppe1_m = Agrif_Parent(ppe1_m)/Agrif_Rhox() 
     340            ppe2_m = Agrif_Parent(ppe2_m)/Agrif_Rhoy()           
     341         ENDIF 
     342#endif          
     343 
     344         DO jj = 1, jpj 
     345            DO ji = 1, jpi 
     346               zti = FLOAT( ji - 1 + nimpp - 1 )         ;   ztj = FLOAT( jj - 1 + njmpp - 1 ) 
     347               zui = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zuj = FLOAT( jj - 1 + njmpp - 1 ) 
     348               zvi = FLOAT( ji - 1 + nimpp - 1 )         ;   zvj = FLOAT( jj - 1 + njmpp - 1 ) + 0.5 
     349               zfi = FLOAT( ji - 1 + nimpp - 1 ) + 0.5   ;   zfj = FLOAT( jj - 1 + njmpp - 1 ) + 0.5 
     350 
     351               glamt(ji,jj) = glam0 + ppe1_m * 1.e-5 * zti 
     352               glamu(ji,jj) = glam0 + ppe1_m * 1.e-5 * zui 
     353               glamv(ji,jj) = glam0 + ppe1_m * 1.e-5 * zvi 
     354               glamf(ji,jj) = glam0 + ppe1_m * 1.e-5 * zfi 
     355    
     356               gphit(ji,jj) = gphi0 + ppe2_m * 1.e-5 * ztj 
     357               gphiu(ji,jj) = gphi0 + ppe2_m * 1.e-5 * zuj 
     358               gphiv(ji,jj) = gphi0 + ppe2_m * 1.e-5 * zvj 
     359               gphif(ji,jj) = gphi0 + ppe2_m * 1.e-5 * zfj 
     360            END DO 
     361         END DO 
     362          
     363         ! Horizontal scale factors (in meters) 
     364         !                              ====== 
     365!! ==> EITHER 1) variable scale factors 
     366         DO jj = 1, jpj 
     367            DO ji = 1, jpi 
     368               !!e1t(ji,jj) = ppe1_m * EXP( -0.8/REAL(jpiglo**2) * (mi0(ji)-REAL(jpiglo+1)*0.5)**2 )  ! gaussian shape 
     369               !!e2t(ji,jj) = ppe2_m * EXP( -0.8/REAL(jpjglo**2) * (mj0(jj)-REAL(jpjglo+1)*0.5)**2 )  ! gaussian shape 
     370               e1t(ji,jj) = ppe1_m * ( 1. -0.1 * ABS(REAL(mi0(ji))-REAL(jpiglo+1)*0.5) / (1.-REAL(jpiglo+1)*0.5) ) ! linear shape 
     371               e2t(ji,jj) = ppe2_m * ( 1. -0.1 * ABS(REAL(mj0(jj))-REAL(jpjglo+1)*0.5) / (1.-REAL(jpjglo+1)*0.5) ) ! linear shape 
     372            END DO 
     373         END DO 
     374#if defined key_agrif  
     375         IF( .NOT. Agrif_Root() ) THEN ! only works if the zoom is positioned at the center of the parent grid 
     376            DO jj = 1, jpj 
     377               DO ji = 1, jpi 
     378                  e1t(ji,jj) = ppe1_m * ( 1. -0.1 * ABS(REAL(mi0(ji))-REAL(jpiglo+1)*0.5) / (1.-REAL(jpiglo+1)*0.5)  & 
     379                     &                            * REAL(jpiglo) / REAL(Agrif_Parent(jpiglo) * Agrif_Rhox()) )       ! factor to match parent grid 
     380                  e2t(ji,jj) = ppe2_m * ( 1. -0.1 * ABS(REAL(mj0(jj))-REAL(jpjglo+1)*0.5) / (1.-REAL(jpjglo+1)*0.5)  & 
     381                     &                            * REAL(jpjglo) / REAL(Agrif_Parent(jpjglo) * Agrif_Rhoy()) )       ! factor to match parent grid 
     382               END DO 
     383            END DO 
     384         ENDIF 
     385#endif 
     386!! ==> OR 2) constant scale factors 
     387!!         e1t(:,:) = ppe1_m 
     388!!         e2t(:,:) = ppe2_m 
     389          
     390         e1u(:,:) = e1t(:,:)      ;      e2u(:,:) = e2t(:,:) 
     391         e1v(:,:) = e1t(:,:)      ;      e2v(:,:) = e2t(:,:) 
     392         e1f(:,:) = e1t(:,:)      ;      e2f(:,:) = e2t(:,:) 
     393 
    324394      CASE DEFAULT 
    325395         WRITE(ctmp1,*) '          bad flag value for jphgr_msh = ', jphgr_msh 
    … …  
    377447      CASE ( 0, 1, 4 )               ! mesh on the sphere 
    378448         ! 
    379          ff(:,:) = 2. * omega * SIN( rad * gphif(:,:) )  
     449         ff  (:,:) = 2. * omega * SIN( rad * gphif(:,:) )  
     450         ff_t(:,:) = 2. * omega * SIN( rad * gphit(:,:) )  
    380451         ! 
    381452      CASE ( 2 )                     ! f-plane at ppgphi0  
    382453         ! 
    383          ff(:,:) = 2. * omega * SIN( rad * ppgphi0 ) 
     454         ff  (:,:) = 2. * omega * SIN( rad * ppgphi0 ) 
     455         ff_t(:,:) = 2. * omega * SIN( rad * ppgphi0 )  ! clem: coriolis at T-point 
    384456         ! 
    385457         IF(lwp) WRITE(numout,*) '          f-plane: Coriolis parameter = constant = ', ff(1,1) 
    … …  
    399471         zf0     = 2. * omega * SIN( rad * zphi0 )                              ! compute f0 1st point south 
    400472         ! 
    401          ff(:,:) = ( zf0  + zbeta * gphif(:,:) * 1.e+3 )                        ! f = f0 +beta* y ( y=0 at south) 
     473         ff  (:,:) = ( zf0  + zbeta * gphif(:,:) * 1.e+3 )                        ! f = f0 +beta* y ( y=0 at south) 
     474         ff_t(:,:) = ( zf0  + zbeta * gphit(:,:) * 1.e+3 )                        ! clem: coriolis at T-point 
    402475         ! 
    403476         IF(lwp) THEN 
    … …  
    420493         zf0   = 2. * omega * SIN( rad * zphi0 )                            ! compute f0 1st point south 
    421494         ! 
    422          ff(:,:) = ( zf0 + zbeta * ABS( gphif(:,:) - zphi0 ) * rad * ra )   ! f = f0 +beta* y ( y=0 at south) 
     495         ff  (:,:) = ( zf0 + zbeta * ABS( gphif(:,:) - zphi0 ) * rad * ra )   ! f = f0 +beta* y ( y=0 at south) 
     496         ff_t(:,:) = ( zf0 + zbeta * ABS( gphit(:,:) - zphi0 ) * rad * ra )   ! clem: coriolis at T-point 
    423497         ! 
    424498         IF(lwp) THEN 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r6492 r7309  
    145145      ! Build the vertical coordinate system 
    146146      ! ------------------------------------ 
     147#if defined key_sas2D 
     148      WRITE(numout,*) ' domzgr: we use SAS2D (i.e. no ocean) with jpk=',jpk 
     149      mbathy(:,:) = 1   ;   bathy(:,:) = rn_hmin 
     150 
     151      gdept_0 (:,:,:) = rn_hmin 
     152      gdepw_0 (:,:,:) = rn_hmin   ;   gdep3w_0(:,:,:) = rn_hmin 
     153      gdept_1d(:)     = rn_hmin   ;   gdepw_1d(:)     = rn_hmin 
     154 
     155      e3t_0 (:,:,:) = rn_hmin 
     156      e3u_0 (:,:,:) = rn_hmin   ;   e3v_0 (:,:,:) = rn_hmin 
     157      e3f_0 (:,:,:) = rn_hmin   ;   e3w_0 (:,:,:) = rn_hmin 
     158      e3uw_0(:,:,:) = rn_hmin   ;   e3vw_0(:,:,:) = rn_hmin 
     159      e3t_1d(:)     = rn_hmin   ;   e3w_1d(:)     = rn_hmin 
     160 
     161      mikt(:,:) = 1   ;   mikv(:,:) = 1 
     162      miku(:,:) = 1   ;   mikf(:,:) = 1 
     163#else 
    147164                          CALL zgr_z            ! Reference z-coordinate system (always called) 
    148165                          CALL zgr_bat          ! Bathymetry fields (levels and meters) 
    … …  
    164181      END IF 
    165182      ! 
     183#endif 
     184       
    166185      IF( nprint == 1 .AND. lwp )   THEN 
    167186         WRITE(numout,*) ' MIN val mbathy  ', MINVAL(  mbathy(:,:) ), ' MAX ', MAXVAL( mbathy(:,:) ) 
    … …  
    476495            idta( 1    , :    ) = 0                 ;      zdta( 1    , :    ) =  0._wp 
    477496            idta(jpidta, :    ) = 0                 ;      zdta(jpidta, :    ) =  0._wp 
     497         ELSEIF( jperio == 7 ) THEN 
     498!           Nothing to do here 
    478499         ELSE 
    479500            ih = 0                                  ;      zh = 0._wp 
    … …  
    738759         IF( lk_mpp ) THEN 
    739760            IF( nbondi == -1 .OR. nbondi == 2 ) THEN 
    740                IF( jperio /= 1 )   mbathy(1,:) = 0 
     761               IF( jperio /= 1 .AND. jperio /= 7 )   mbathy(1,:) = 0 
    741762            ENDIF 
    742763            IF( nbondi == 1 .OR. nbondi == 2 ) THEN 
    743                IF( jperio /= 1 )   mbathy(nlci,:) = 0 
     764               IF( jperio /= 1 .AND. jperio /= 7 )   mbathy(nlci,:) = 0 
    744765            ENDIF 
    745766         ELSE 
    … …  
    756777         mbathy( 1 ,:) = mbathy(jpim1,:) 
    757778         mbathy(jpi,:) = mbathy(  2  ,:) 
     779         IF (jperio == 7) THEN 
     780            IF(lwp) WRITE(numout,*)' north south boundary conditions on mbathy: jperio = ', jperio 
     781            mbathy( : ,1) = mbathy(:, jpjm1) 
     782            mbathy(:, jpj)= mbathy(:,2) 
     783         ENDIF 
    758784      ELSEIF( nperio == 2 ) THEN 
    759785         IF(lwp) WRITE(numout,*) '   equatorial boundary conditions on mbathy: nperio = ', nperio 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90

    r6519 r7309  
    9696      CHARACTER(len=19) :: cldate  
    9797      CHARACTER(len=10) :: clname 
    98       INTEGER           ::   ji 
     98      INTEGER           :: ji, jkmin 
    9999      ! 
    100100      REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z_bnds 
    … …  
    169169 
    170170      ! Add vertical grid bounds 
     171      jkmin = MIN(2,jpk)  ! in case jpk=1 (i.e. sas2D) 
    171172      z_bnds(:      ,1) = gdepw_1d(:) 
    172       z_bnds(1:jpkm1,2) = gdepw_1d(2:jpk) 
     173      z_bnds(1:jpkm1,2) = gdepw_1d(jkmin:jpk) 
    173174      z_bnds(jpk:   ,2) = gdepw_1d(jpk) + e3t_1d(jpk) 
    174175      CALL iom_set_axis_attr( "deptht", bounds=z_bnds ) 
    175176      CALL iom_set_axis_attr( "depthu", bounds=z_bnds ) 
    176177      CALL iom_set_axis_attr( "depthv", bounds=z_bnds ) 
    177       z_bnds(:    ,2) = gdept_1d(:) 
    178       z_bnds(2:jpk,1) = gdept_1d(1:jpkm1) 
    179       z_bnds(1    ,1) = gdept_1d(1) - e3w_1d(1) 
     178      z_bnds(:        ,2) = gdept_1d(:) 
     179      z_bnds(jkmin:jpk,1) = gdept_1d(1:jpkm1) 
     180      z_bnds(1        ,1) = gdept_1d(1) - e3w_1d(1) 
    180181      CALL iom_set_axis_attr( "depthw", bounds=z_bnds ) 
    181182 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/IOM/iom_def.F90

    r6140 r7309  
    3939   INTEGER, PARAMETER, PUBLIC ::   jp_i1    = 204      !: write INTEGER(1) 
    4040 
    41    INTEGER, PARAMETER, PUBLIC ::   jpmax_files  = 100   !: maximum number of simultaneously opened file 
    42    INTEGER, PARAMETER, PUBLIC ::   jpmax_vars   = 600 !: maximum number of variables in one file 
     41   INTEGER, PARAMETER, PUBLIC ::   jpmax_files  = 100  !: maximum number of simultaneously opened file 
     42   INTEGER, PARAMETER, PUBLIC ::   jpmax_vars   = 1200 !: maximum number of variables in one file 
    4343   INTEGER, PARAMETER, PUBLIC ::   jpmax_dims   =  4   !: maximum number of dimensions for one variable 
    4444   INTEGER, PARAMETER, PUBLIC ::   jpmax_digits =  5   !: maximum number of digits for the cpu number in the file name 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90

    r6918 r7309  
    405405                                         ptab(nlci-jpreci+1:jpi   ,:,:) = zland    ! north 
    406406         ENDIF 
    407          !                                   ! North-South boundaries (always closed) 
    408          IF( .NOT. cd_type == 'F' )   ptab(:,     1       :jprecj,:) = zland       ! south except F-point 
    409                                       ptab(:,nlcj-jprecj+1:jpj   ,:) = zland       ! north 
     407                                          ! North-south cyclic 
     408         IF ( nbondj == 2 .AND. jperio == 7 )    THEN !* cyclic north south only with no mpp split in latitude 
     409            ptab(:,1 , :) = ptab(:, jpjm1,:) 
     410            ptab(:,jpj,:) = ptab(:,     2,:) 
     411         ELSE   !                                   ! North-South boundaries (closed) 
     412            IF( .NOT. cd_type == 'F' )   ptab(:,     1       :jprecj,:) = zland       ! south except F-point 
     413                                         ptab(:,nlcj-jprecj+1:jpj   ,:) = zland       ! north 
     414         ENDIF 
    410415         ! 
    411416      ENDIF 
    … …  
    608613                                                   pt2d_array(ii)%pt2d(nlci-jpreci+1 : jpi   ,:) = zland    ! north 
    609614            ENDIF 
    610             !                                   ! North-South boundaries (always closed) 
     615                                                ! Noth-South boundaries 
     616            IF ( nbondj == 2 .AND. jperio == 7 )    THEN !* cyclic north south 
     617               pt2d_array(ii)%pt2d(:, 1   ) =   pt2d_array(ii)%pt2d(:, jpjm1 ) 
     618               pt2d_array(ii)%pt2d(:, jpj ) =   pt2d_array(ii)%pt2d(:, 2 )           
     619            ELSE   !              
     620               !                                   ! North-South boundaries (closed) 
    611621               IF( .NOT. type_array(ii) == 'F' )   pt2d_array(ii)%pt2d(:,             1:jprecj ) = zland    ! south except F-point 
    612622                                                   pt2d_array(ii)%pt2d(:, nlcj-jprecj+1:jpj    ) = zland    ! north 
    613623            ! 
    614          ENDIF 
     624            ENDIF 
     625          ENDIF 
    615626      END DO 
    616627 
    … …  
    888899                                         pt2d(nlci-jpreci+1:jpi   ,:) = zland    ! north 
    889900         ENDIF 
    890          !                                   ! North-South boundaries (always closed) 
     901                                            ! North-South boudaries 
     902         IF ( nbondj == 2 .AND. jperio == 7 )    THEN !* cyclic north south 
     903            pt2d(:,  1 ) = pt2d(:,jpjm1) 
     904            pt2d(:, jpj) = pt2d(:,    2) 
     905         ELSE     
     906         !                                   ! North-South boundaries (closed) 
    891907            IF( .NOT. cd_type == 'F' )   pt2d(:,     1       :jprecj) = zland    !south except F-point 
    892908                                         pt2d(:,nlcj-jprecj+1:jpj   ) = zland    ! north 
    893          ! 
     909         ENDIF      
    894910      ENDIF 
    895911 
    … …  
    10711087                                       ptab2(nlci-jpreci+1:jpi   ,:,:) = 0.e0 
    10721088      ENDIF 
    1073  
    1074  
    1075       !                                      ! North-South boundaries 
     1089                                            ! North-South boundaries 
     1090      IF ( nbondj == 2 .AND. jperio == 7 )    THEN !* cyclic north south 
     1091         ptab1(:,     1       ,:) = ptab1(: ,  jpjm1 , :) 
     1092         ptab1(:,   jpj       ,:) = ptab1(: ,      2 , :) 
     1093         ptab2(:,     1       ,:) = ptab2(: ,  jpjm1 , :) 
     1094         ptab2(:,   jpj       ,:) = ptab2(: ,      2 , :) 
     1095      ELSE      
     1096      !                                      ! North-South boundaries closed 
    10761097      IF( .NOT. cd_type1 == 'F' )   ptab1(:,     1       :jprecj,:) = 0.e0    ! south except at F-point 
    10771098      IF( .NOT. cd_type2 == 'F' )   ptab2(:,     1       :jprecj,:) = 0.e0 
    10781099                                    ptab1(:,nlcj-jprecj+1:jpj   ,:) = 0.e0    ! north 
    10791100                                    ptab2(:,nlcj-jprecj+1:jpj   ,:) = 0.e0 
    1080  
     1101      ENDIF      
    10811102 
    10821103      ! 2. East and west directions exchange 
    … …  
    12671288      ! Order matters Here !!!! 
    12681289      ! 
    1269       !                                      !* North-South boundaries (always colsed) 
     1290                                           ! North-South cyclic 
     1291      IF ( nbondj == 2 .AND. jperio == 7 )    THEN !* cyclic north south 
     1292         pt2d(:, 1-jprj:  1     ) = pt2d ( :, jpjm1-jprj:jpjm1) 
     1293         pt2d(:, jpj   :jpj+jprj) = pt2d ( :, 2         :2+jprj) 
     1294      ELSE 
     1295         
     1296      !                                      !* North-South boundaries (closed) 
    12701297      IF( .NOT. cd_type == 'F' )   pt2d(:,  1-jprj   :  jprecj  ) = 0.e0    ! south except at F-point 
    12711298                                   pt2d(:,nlcj-jprecj+1:jpj+jprj) = 0.e0    ! north 
    1272  
     1299      ENDIF 
     1300                                 
    12731301      !                                      ! East-West boundaries 
    12741302      !                                           !* Cyclic east-west 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/LBC/mppini.F90

    r6412 r7309  
    7676          &              'the domain is lay out for distributed memory computing! ' ) 
    7777 
     78      IF( jperio == 7 ) CALL ctl_stop( ' jperio = 7 needs distributed memory computing ',   & 
     79          &              ' with 1 process. Add key_mpp_mpi in the list of active cpp keys ' ) 
    7880   END SUBROUTINE mpp_init 
    7981 
    … …  
    379381      ! w a r n i n g  narea (zone) /= nproc (processors)! 
    380382 
    381       IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) THEN 
     383      IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7 ) THEN 
    382384         IF( jpni == 1 )THEN 
    383385            nbondi = 2 
    … …  
    446448      ENDIF 
    447449 
     450      IF( jperio == 7 .AND. ( jpni /= 1 .OR. jpnj /= 1 ) ) & 
     451         &                  CALL ctl_stop( ' mpp_init: error jperio = 7 works only with jpni = jpnj = 1' ) 
    448452      IF( nperio == 1 .AND. jpni /= 1 ) CALL ctl_stop( ' mpp_init: error on cyclicity' ) 
    449453 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90

    r6140 r7309  
    44   !! Ocean forcing:  analytical momentum, heat and freshwater forcings 
    55   !!===================================================================== 
    6    !! History :  3.0   ! 2006-06  (G. Madec)  Original code 
    7    !!            3.2   ! 2009-07  (G. Madec)  Style only 
     6   !! History :  3.0   ! 2006-06  (G. Madec)    Original code 
     7   !!            3.2   ! 2009-07  (G. Madec)    Style only 
     8   !!            3.7   ! 2016-10  (C. Rousset)  Add analytic for LIM3 (ana_ice) 
    89   !!---------------------------------------------------------------------- 
    910 
    … …  
    1516   USE dom_oce         ! ocean space and time domain 
    1617   USE sbc_oce         ! Surface boundary condition: ocean fields 
     18   USE sbc_ice         ! Surface boundary condition: ice   fields 
    1719   USE phycst          ! physical constants 
    1820   USE in_out_manager  ! I/O manager 
    … …  
    2022   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2123   USE lib_fortran 
    22  
     24   USE wrk_nemo 
     25#if defined key_lim3 
     26   USE ice, ONLY       : pfrld, a_i_b 
     27   USE limthd_dh       ! for CALL lim_thd_snwblow 
     28#endif 
     29    
    2330   IMPLICIT NONE 
    2431   PRIVATE 
    2532 
    26    PUBLIC   sbc_ana    ! routine called in sbcmod module 
    27    PUBLIC   sbc_gyre   ! routine called in sbcmod module 
     33   PUBLIC   sbc_ana         ! routine called in sbcmod module 
     34   PUBLIC   sbc_gyre        ! routine called in sbcmod module 
     35#if defined key_lim3 
     36   PUBLIC   ana_ice_tau     ! routine called in sbc_ice_lim module 
     37   PUBLIC   ana_ice_flx     ! routine called in sbc_ice_lim module 
     38#endif 
    2839 
    2940   !                       !!* Namelist namsbc_ana * 
    30    INTEGER  ::   nn_tau000  ! nb of time-step during which the surface stress 
    31    !                        ! increase from 0 to its nominal value  
    32    REAL(wp) ::   rn_utau0   ! constant wind stress value in i-direction 
    33    REAL(wp) ::   rn_vtau0   ! constant wind stress value in j-direction 
    34    REAL(wp) ::   rn_qns0    ! non solar heat flux 
    35    REAL(wp) ::   rn_qsr0    !     solar heat flux 
    36    REAL(wp) ::   rn_emp0    ! net freshwater flux 
     41   ! --- oce variables --- ! 
     42   INTEGER  ::   nn_tau000 ! nb of time-step during which the surface stress 
     43   !                       ! increase from 0 to its nominal value  
     44   REAL(wp) ::   rn_utau0  ! constant wind stress value in i-direction 
     45   REAL(wp) ::   rn_vtau0  ! constant wind stress value in j-direction 
     46   REAL(wp) ::   rn_qns0   ! non solar heat flux 
     47   REAL(wp) ::   rn_qsr0   !     solar heat flux 
     48   REAL(wp) ::   rn_emp0   ! net freshwater flux 
     49   ! --- ice variables --- ! 
     50   REAL(wp) ::   rn_iutau0 ! constant wind stress value in i-direction over ice 
     51   REAL(wp) ::   rn_ivtau0 ! constant wind stress value in j-direction over ice 
     52   REAL(wp) ::   rn_iqns0  ! non solar heat flux over ice 
     53   REAL(wp) ::   rn_iqsr0  !     solar heat flux over ice 
     54   REAL(wp) ::   rn_sprec0 ! snow precip 
     55   REAL(wp) ::   rn_ievap0 ! sublimation 
    3756    
    3857   !! * Substitutions 
    … …  
    6887      REAL(wp) ::   zcoef, zty, zmod      !   -      - 
    6988      !! 
    70       NAMELIST/namsbc_ana/ nn_tau000, rn_utau0, rn_vtau0, rn_qns0, rn_qsr0, rn_emp0 
     89      NAMELIST/namsbc_ana/ nn_tau000, rn_utau0, rn_vtau0, rn_qns0, rn_qsr0, rn_emp0,  & 
     90         &                 rn_iutau0, rn_ivtau0, rn_iqsr0, rn_iqns0, rn_sprec0, rn_ievap0 
    7191      !!--------------------------------------------------------------------- 
    7292      ! 
    … …  
    85105         IF(lwp) WRITE(numout,*)' sbc_ana : Constant surface fluxes read in namsbc_ana namelist' 
    86106         IF(lwp) WRITE(numout,*)' ~~~~~~~ ' 
    87          IF(lwp) WRITE(numout,*)'              spin up of the stress  nn_tau000 = ', nn_tau000, ' time-steps' 
    88          IF(lwp) WRITE(numout,*)'              constant i-stress      rn_utau0  = ', rn_utau0 , ' N/m2' 
    89          IF(lwp) WRITE(numout,*)'              constant j-stress      rn_vtau0  = ', rn_vtau0 , ' N/m2' 
    90          IF(lwp) WRITE(numout,*)'              non solar heat flux    rn_qns0   = ', rn_qns0  , ' W/m2' 
    91          IF(lwp) WRITE(numout,*)'              solar heat flux        rn_qsr0   = ', rn_qsr0  , ' W/m2' 
    92          IF(lwp) WRITE(numout,*)'              net heat flux          rn_emp0   = ', rn_emp0  , ' Kg/m2/s' 
     107         IF(lwp) WRITE(numout,*)'              spin up of the stress         nn_tau000 = ', nn_tau000 , ' time-steps' 
     108         IF(lwp) WRITE(numout,*)'              constant i-stress             rn_utau0  = ', rn_utau0  , ' N/m2' 
     109         IF(lwp) WRITE(numout,*)'              constant j-stress             rn_vtau0  = ', rn_vtau0  , ' N/m2' 
     110         IF(lwp) WRITE(numout,*)'              non solar heat flux           rn_qns0   = ', rn_qns0   , ' W/m2' 
     111         IF(lwp) WRITE(numout,*)'              solar heat flux               rn_qsr0   = ', rn_qsr0   , ' W/m2' 
     112         IF(lwp) WRITE(numout,*)'              net freshwater flux           rn_emp0   = ', rn_emp0   , ' Kg/m2/s' 
     113         IF(lwp) WRITE(numout,*)'              constant ice-atm stress       rn_iutau0 = ', rn_iutau0 , ' N/m2' 
     114         IF(lwp) WRITE(numout,*)'              constant ice-atm stress       rn_ivtau0 = ', rn_ivtau0 , ' N/m2' 
     115         IF(lwp) WRITE(numout,*)'              solar heat flux over ice      rn_iqsr0  = ', rn_iqsr0  , ' W/m2' 
     116         IF(lwp) WRITE(numout,*)'              non solar heat flux over ice  rn_iqns0  = ', rn_iqns0  , ' W/m2' 
     117         IF(lwp) WRITE(numout,*)'              snow precip                   rn_sprec0 = ', rn_sprec0 , ' Kg/m2/s' 
     118         IF(lwp) WRITE(numout,*)'              sublimation                   rn_ievap0 = ', rn_ievap0 , ' Kg/m2/s' 
    93119         ! 
    94120         nn_tau000 = MAX( nn_tau000, 1 )     ! must be >= 1 
    … …  
    132158   END SUBROUTINE sbc_ana 
    133159 
    134  
     160#if defined key_lim3 
     161   SUBROUTINE ana_ice_tau 
     162      !!--------------------------------------------------------------------- 
     163      !!                     ***  ROUTINE ana_ice_tau  *** 
     164      !! 
     165      !! ** Purpose :   provide the surface boundary (momentum) condition over sea-ice 
     166      !!--------------------------------------------------------------------- 
     167      utau_ice(:,:) = rn_iutau0 
     168      vtau_ice(:,:) = rn_ivtau0 
     169      
     170   END SUBROUTINE ana_ice_tau 
     171    
     172   SUBROUTINE ana_ice_flx 
     173      !!--------------------------------------------------------------------- 
     174      !!                     ***  ROUTINE ana_ice_flx  *** 
     175      !! 
     176      !! ** Purpose :   provide the surface boundary (flux) condition over sea-ice 
     177      !!--------------------------------------------------------------------- 
     178      REAL(wp), DIMENSION(:,:), POINTER ::   zsnw       ! snw distribution after wind blowing 
     179      !!--------------------------------------------------------------------- 
     180      CALL wrk_alloc( jpi,jpj, zsnw )  
     181 
     182      ! ocean variables (renaming) 
     183      emp_oce (:,:)   = rn_emp0 
     184      qsr_oce (:,:)   = rn_qsr0 
     185      qns_oce (:,:)   = rn_qns0 
     186       
     187      ! ice variables 
     188      alb_ice (:,:,:) = 0.7_wp ! useless 
     189      qsr_ice (:,:,:) = rn_iqsr0 
     190      qns_ice (:,:,:) = rn_iqns0 
     191      sprecip (:,:)   = rn_sprec0 
     192      evap_ice(:,:,:) = rn_ievap0 
     193 
     194      ! ice variables deduced from above 
     195      zsnw(:,:) = 1._wp 
     196      !!CALL lim_thd_snwblow( pfrld, zsnw )  ! snow distribution over ice after wind blowing  
     197      emp_ice  (:,:)   = SUM( a_i_b(:,:,:) * evap_ice(:,:,:), dim=3 ) - sprecip(:,:) * zsnw(:,:) 
     198      emp_oce  (:,:)   = emp_oce(:,:) - sprecip(:,:) * (1._wp - zsnw(:,:) ) 
     199      qevap_ice(:,:,:) =   0._wp 
     200      qprec_ice(:,:)   =   rhosn * ( sst_m(:,:) * cpic - lfus ) * tmask(:,:,1) ! in J/m3 
     201      qemp_oce (:,:)   = - emp_oce(:,:) * sst_m(:,:) * rcp 
     202      qemp_ice (:,:)   =   sprecip(:,:) * zsnw * ( sst_m(:,:) * cpic - lfus ) * tmask(:,:,1) ! solid precip (only) 
     203 
     204      ! total fluxes 
     205      emp_tot (:,:) = emp_ice  + emp_oce  
     206      qns_tot (:,:) = pfrld(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:) 
     207      qsr_tot (:,:) = pfrld(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 ) 
     208 
     209      !-------------------------------------------------------------------- 
     210      ! FRACTIONs of net shortwave radiation which is not absorbed in the 
     211      ! thin surface layer and penetrates inside the ice cover 
     212      ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 ) 
     213      fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 
     214      fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 
     215 
     216      CALL wrk_dealloc( jpi,jpj, zsnw )  
     217       
     218   END SUBROUTINE ana_ice_flx 
     219#endif 
     220 
     221    
    135222   SUBROUTINE sbc_gyre( kt ) 
    136223      !!--------------------------------------------------------------------- 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r6813 r7309  
    1616   !!            3.4  !  2011-11  (C. Harris) Fill arrays required by CICE 
    1717   !!            3.7  !  2014-06  (L. Brodeau) simplification and optimization of CORE bulk 
     18   !!            4.0  !  2016-06  (C. Rousset) Add new param of drags with sea-ice (Lupkes at al 2012) 
    1819   !!---------------------------------------------------------------------- 
    1920 
    … …  
    3839   USE lib_fortran    ! to use key_nosignedzero 
    3940#if defined key_lim3 
    40    USE ice     , ONLY :   u_ice, v_ice, jpl, pfrld, a_i_b 
     41   USE ice, ONLY      : u_ice, v_ice, jpl, pfrld, a_i_b, at_i_b 
    4142   USE limthd_dh      ! for CALL lim_thd_snwblow 
    4243#elif defined key_lim2 
    43    USE ice_2   , ONLY :  u_ice, v_ice 
    44    USE par_ice_2      ! LIM-2 parameters 
     44   USE ice_2, ONLY    : u_ice, v_ice 
     45   USE par_ice_2 
    4546#endif 
    4647   ! 
    … …  
    6162   PUBLIC   blk_ice_core_flx     ! routine called in sbc_ice_lim module 
    6263#endif 
    63    PUBLIC   turb_core_2z         ! routine calles in sbcblk_mfs module 
     64   PUBLIC   turb_core_2z         ! routine called in sbcblk_mfs module 
    6465 
    6566   INTEGER , PARAMETER ::   jpfld   = 9           ! maximum number of files to read 
    … …  
    7778 
    7879   !                                             !!! CORE bulk parameters 
    79    REAL(wp), PARAMETER ::   rhoa =    1.22        ! air density 
    80    REAL(wp), PARAMETER ::   cpa  = 1000.5         ! specific heat of air 
    81    REAL(wp), PARAMETER ::   Lv   =    2.5e6       ! latent heat of vaporization 
    82    REAL(wp), PARAMETER ::   Ls   =    2.839e6     ! latent heat of sublimation 
    83    REAL(wp), PARAMETER ::   Stef =    5.67e-8     ! Stefan Boltzmann constant 
    84    REAL(wp), PARAMETER ::   Cice =    1.4e-3      ! iovi 1.63e-3     ! transfer coefficient over ice 
    85    REAL(wp), PARAMETER ::   albo =    0.066       ! ocean albedo assumed to be constant 
     80   REAL(wp), PARAMETER ::   rhoa   =    1.22        ! air density 
     81   REAL(wp), PARAMETER ::   cpa    = 1000.5         ! specific heat of air 
     82   REAL(wp), PARAMETER ::   Lv     =    2.5e6       ! latent heat of vaporization 
     83   REAL(wp), PARAMETER ::   Ls     =    2.839e6     ! latent heat of sublimation 
     84   REAL(wp), PARAMETER ::   Stef   =    5.67e-8     ! Stefan Boltzmann constant 
     85   REAL(wp), PARAMETER ::   Cd_ice =    1.4e-3      ! transfer coefficient over ice 
     86   REAL(wp), PARAMETER ::   albo   =    0.066       ! ocean albedo assumed to be constant 
    8687 
    8788   !                        !!* Namelist namsbc_core : CORE bulk parameters 
    … …  
    9293   REAL(wp) ::   rn_zqt      ! z(q,t) : height of humidity and temperature measurements 
    9394   REAL(wp) ::   rn_zu       ! z(u)   : height of wind measurements 
    94  
     95   LOGICAL  ::   ln_Cd_L12 = .FALSE. !  Modify the drag ice-atm and oce-atm depending on ice concentration (from Lupkes et al. JGR2012) 
     96 
     97   ! 
     98   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_oce   ! air-ocean drag (clem) 
     99    
    95100   !! * Substitutions 
    96101#  include "vectopt_loop_substitute.h90" 
    … …  
    102107CONTAINS 
    103108 
     109   INTEGER FUNCTION sbc_blk_core_alloc() 
     110      !!------------------------------------------------------------------- 
     111      !!             ***  ROUTINE sbc_blk_core_alloc (clem) *** 
     112      !!------------------------------------------------------------------- 
     113      ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_core_alloc ) 
     114      ! 
     115      IF( lk_mpp                  )   CALL mpp_sum( sbc_blk_core_alloc ) 
     116      IF( sbc_blk_core_alloc /= 0 )   CALL ctl_warn('sbc_blk_core_alloc: failed to allocate arrays') 
     117   END FUNCTION sbc_blk_core_alloc 
     118 
     119    
    104120   SUBROUTINE sbc_blk_core( kt ) 
    105121      !!--------------------------------------------------------------------- 
    … …  
    149165      TYPE(FLD_N) ::   sn_tdif                                 !   "                                 " 
    150166      NAMELIST/namsbc_core/ cn_dir , ln_taudif, rn_pfac, rn_efac, rn_vfac,  & 
    151          &                  sn_wndi, sn_wndj  , sn_humi, sn_qsr ,           & 
    152          &                  sn_qlw , sn_tair  , sn_prec, sn_snow,           & 
    153          &                  sn_tdif, rn_zqt   ,  rn_zu 
     167         &                  sn_wndi, sn_wndj, sn_humi  , sn_qsr ,           & 
     168         &                  sn_qlw , sn_tair, sn_prec  , sn_snow,           & 
     169         &                  sn_tdif, rn_zqt,  rn_zu    , ln_Cd_L12 
    154170      !!--------------------------------------------------------------------- 
    155171      ! 
    … …  
    157173      IF( kt == nit000 ) THEN                   !  First call kt=nit000  ! 
    158174         !                                      ! ====================== ! 
     175         ! 
     176         !                                      ! allocate sbc_blk_core array (clem) 
     177         IF( sbc_blk_core_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_blk_core : unable to allocate standard arrays' ) 
    159178         ! 
    160179         REWIND( numnam_ref )              ! Namelist namsbc_core in reference namelist : CORE bulk parameters 
    … …  
    321340         &               Cd, Ch, Ce, zt_zu, zq_zu ) 
    322341     
     342      Cd_oce(:,:) = Cd(:,:)  ! record value of pure ocean-atm. drag (clem) 
     343      
    323344      ! ... tau module, i and j component 
    324345      DO jj = 1, jpj 
    … …  
    439460      !!--------------------------------------------------------------------- 
    440461      INTEGER  ::   ji, jj    ! dummy loop indices 
    441       REAL(wp) ::   zcoef_wnorm, zcoef_wnorm2 
    442462      REAL(wp) ::   zwnorm_f, zwndi_f , zwndj_f               ! relative wind module and components at F-point 
    443463      REAL(wp) ::             zwndi_t , zwndj_t               ! relative wind components at T-point 
     464      REAL(wp), DIMENSION(:,:), POINTER ::   Cd               ! transfer coefficient for momentum      (tau) 
    444465      !!--------------------------------------------------------------------- 
    445466      ! 
    446467      IF( nn_timing == 1 )  CALL timing_start('blk_ice_core_tau') 
    447468      ! 
    448       ! local scalars ( place there for vector optimisation purposes) 
    449       zcoef_wnorm  = rhoa * Cice 
    450       zcoef_wnorm2 = rhoa * Cice * 0.5 
     469      CALL wrk_alloc( jpi,jpj, Cd ) 
     470 
     471      Cd(:,:) = Cd_ice 
     472       
     473      ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 
     474#if defined key_lim3 
     475      IF( ln_Cd_L12 ) THEN 
     476         CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 
     477      ENDIF 
     478#endif 
    451479 
    452480!!gm brutal.... 
    … …  
    469497               zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
    470498                  &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
    471                zwnorm_f = zcoef_wnorm * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
     499               zwnorm_f = rhoa * Cd(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
    472500               ! ... ice stress at I-point 
    473501               utau_ice(ji,jj) = zwnorm_f * zwndi_f 
    … …  
    478506               zwndj_t = sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.25 * (  v_ice(ji,jj+1) + v_ice(ji+1,jj+1)   & 
    479507                  &                                                    + v_ice(ji,jj  ) + v_ice(ji+1,jj  )  ) 
    480                wndm_ice(ji,jj)  = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 
     508               wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 
    481509            END DO 
    482510         END DO 
    … …  
    495523         DO jj = 2, jpjm1 
    496524            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    497                utau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )                          & 
     525               utau_ice(ji,jj) = rhoa * Cd(ji,jj) * 0.5_wp * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )                          & 
    498526                  &          * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 
    499                vtau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )                          & 
     527               vtau_ice(ji,jj) = rhoa * Cd(ji,jj) * 0.5_wp * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )                          & 
    500528                  &          * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 
    501529            END DO 
    … …  
    511539         CALL prt_ctl(tab2d_1=wndm_ice  , clinfo1=' blk_ice_core: wndm_ice : ') 
    512540      ENDIF 
     541 
     542      CALL wrk_dealloc( jpi,jpj, Cd ) 
    513543 
    514544      IF( nn_timing == 1 )  CALL timing_stop('blk_ice_core_tau') 
    … …  
    542572      REAL(wp), DIMENSION(:,:,:), POINTER ::   z_dqsb            ! sensible  heat sensitivity over ice 
    543573      REAL(wp), DIMENSION(:,:)  , POINTER ::   zevap, zsnw       ! evaporation and snw distribution after wind blowing (LIM3) 
     574      REAL(wp), DIMENSION(:,:)  , POINTER ::   Cd                ! transfer coefficient for momentum      (tau) 
    544575      !!--------------------------------------------------------------------- 
    545576      ! 
    … …  
    547578      ! 
    548579      CALL wrk_alloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb )  
     580      CALL wrk_alloc( jpi,jpj, Cd ) 
     581 
     582      Cd(:,:) = Cd_ice 
     583 
     584      ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 
     585#if defined key_lim3 
     586      IF( ln_Cd_L12 ) THEN 
     587         CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 
     588      ENDIF 
     589#endif 
    549590 
    550591      ! local scalars ( place there for vector optimisation purposes) 
    551592      zcoef_dqlw   = 4.0 * 0.95 * Stef 
    552       zcoef_dqla   = -Ls * Cice * 11637800. * (-5897.8) 
    553       zcoef_dqsb   = rhoa * cpa * Cice 
     593      zcoef_dqla   = -Ls * 11637800. * (-5897.8) 
     594      zcoef_dqsb   = rhoa * cpa 
    554595 
    555596      zztmp = 1. / ( 1. - albo ) 
    … …  
    577618               ! ... turbulent heat fluxes 
    578619               ! Sensible Heat 
    579                z_qsb(ji,jj,jl) = rhoa * cpa * Cice * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 
     620               z_qsb(ji,jj,jl) = rhoa * cpa * Cd(ji,jj) * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 
    580621               ! Latent Heat 
    581                qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls  * Cice * wndm_ice(ji,jj)   &                            
     622               qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls  * Cd(ji,jj) * wndm_ice(ji,jj)   &                            
    582623                  &                         * (  11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / rhoa - sf(jp_humi)%fnow(ji,jj,1)  ) ) 
    583624              ! Latent heat sensitivity for ice (Dqla/Dt) 
    584625               IF( qla_ice(ji,jj,jl) > 0._wp ) THEN 
    585                   dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * wndm_ice(ji,jj) / ( zst2 ) * EXP( -5897.8 / ptsu(ji,jj,jl) ) 
     626                  dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * Cd(ji,jj) * wndm_ice(ji,jj) / ( zst2 ) * EXP( -5897.8 / ptsu(ji,jj,jl) ) 
    586627               ELSE 
    587628                  dqla_ice(ji,jj,jl) = 0._wp 
    … …  
    589630 
    590631               ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 
    591                z_dqsb(ji,jj,jl) = zcoef_dqsb * wndm_ice(ji,jj) 
     632               z_dqsb(ji,jj,jl) = zcoef_dqsb * Cd(ji,jj) * wndm_ice(ji,jj) 
    592633 
    593634               ! ----------------------------! 
    … …  
    668709 
    669710      CALL wrk_dealloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 
     711      CALL wrk_dealloc( jpi,jpj, Cd ) 
    670712      ! 
    671713      IF( nn_timing == 1 )  CALL timing_stop('blk_ice_core_flx') 
    … …  
    905947   END FUNCTION psi_h 
    906948 
     949 
     950#if defined key_lim3 
     951   SUBROUTINE Cdn10_Lupkes2012( Cd ) 
     952      !!---------------------------------------------------------------------- 
     953      !!                      ***  ROUTINE  Cdn10_Lupkes2012  *** 
     954      !! 
     955      !! ** Purpose :    Recompute the ice-atm drag at 10m height to make 
     956      !!                 it dependent on edges at leads, melt ponds and flows. 
     957      !!                 After some approximations, this can be resumed to a dependency 
     958      !!                 on ice concentration. 
     959      !!                 
     960      !! ** Method :     The parameterization is taken from Lupkes et al. (2012) eq.(50) 
     961      !!                 with the highest level of approximation: level4, eq.(59) 
     962      !!                 The generic drag over a cell partly covered by ice can be re-written as follows: 
     963      !! 
     964      !!                 Cd = Cdw * (1-A) + Cdi * A + Ce * (1-A)**(nu+1/(10*beta)) * A**mu 
     965      !! 
     966      !!                    Ce = 2.23e-3       , as suggested by Lupkes (eq. 59) 
     967      !!                    nu = mu = beta = 1 , as suggested by Lupkes (eq. 59) 
     968      !!                    A is the concentration of ice minus melt ponds (if any) 
     969      !! 
     970      !!                 This new drag has a parabolic shape (as a function of A) starting at 
     971      !!                 Cdw(say 1.5e-3) for A=0, reaching 1.97e-3 for A~0.5  
     972      !!                 and going down to Cdi(say 1.4e-3) for A=1 
     973      !! 
     974      !!                 It is theoretically applicable to all ice conditions (not only MIZ) 
     975      !!                 => see Lupkes et al (2013) 
     976      !! 
     977      !! ** References : Lupkes et al. JGR 2012 (theory) 
     978      !!                 Lupkes et al. GRL 2013 (application to GCM) 
     979      !! 
     980      !!---------------------------------------------------------------------- 
     981      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Cd 
     982      REAL(wp), PARAMETER ::   zCe   = 2.23e-03_wp 
     983      REAL(wp), PARAMETER ::   znu   = 1._wp 
     984      REAL(wp), PARAMETER ::   zmu   = 1._wp 
     985      REAL(wp), PARAMETER ::   zbeta = 1._wp 
     986      REAL(wp)            ::   zcoef 
     987      !!---------------------------------------------------------------------- 
     988      zcoef = znu + 1._wp / ( 10._wp * zbeta ) 
     989 
     990      ! generic drag over a cell partly covered by ice 
     991      !!Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) +  &                        ! pure ocean drag 
     992      !!   &      Cd_ice      *           at_i_b(:,:)   +  &                        ! pure ice drag 
     993      !!   &      zCe         * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu   ! change due to sea-ice morphology 
     994 
     995      ! ice-atm drag 
     996      Cd(:,:) = Cd_ice +  &                                                          ! pure ice drag 
     997         &      zCe    * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp)   ! change due to sea-ice morphology 
     998       
     999   END SUBROUTINE Cdn10_Lupkes2012 
     1000#endif 
     1001    
    9071002   !!====================================================================== 
    9081003END MODULE sbcblk_core 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r6722 r7309  
    168168#  include "vectopt_loop_substitute.h90" 
    169169   !!---------------------------------------------------------------------- 
    170    !! NEMO/OPA 3.7 , NEMO Consortium (2015) 
     170   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    171171   !! $Id$ 
    172172   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    15991599      ENDIF 
    16001600 
    1601       !! clem: we should output qemp_oce and qemp_ice (at least) 
    1602       IF( iom_use('hflx_snow_cea') )   CALL iom_put( 'hflx_snow_cea', sprecip(:,:) * ( zcptn(:,:) - Lfus ) )   ! heat flux from snow (cell average) 
    1603       !! these diags are not outputed yet 
    1604 !!      IF( iom_use('hflx_rain_cea') )   CALL iom_put( 'hflx_rain_cea', ( tprecip(:,:) - sprecip(:,:) ) * zcptn(:,:) )   ! heat flux from rain (cell average) 
    1605 !!      IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put( 'hflx_snow_ao_cea', sprecip(:,:) * ( zcptn(:,:) - Lfus ) * (1._wp - zsnw(:,:)) ) ! heat flux from snow (cell average) 
    1606 !!      IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put( 'hflx_snow_ai_cea', sprecip(:,:) * ( zcptn(:,:) - Lfus ) * zsnw(:,:) ) ! heat flux from snow (cell average) 
     1601      ! some more outputs 
     1602      IF( iom_use('hflx_snow_cea') )    CALL iom_put('hflx_snow_cea',   sprecip(:,:) * ( zcptn(:,:) - Lfus ) )                       ! heat flux from snow (cell average) 
     1603      IF( iom_use('hflx_rain_cea') )    CALL iom_put('hflx_rain_cea', ( tprecip(:,:) - sprecip(:,:) ) * zcptn(:,:) )                 ! heat flux from rain (cell average) 
     1604      IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea',sprecip(:,:) * ( zcptn(:,:) - Lfus ) * (1._wp - zsnw(:,:)) ) ! heat flux from snow (cell average) 
     1605      IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea',sprecip(:,:) * ( zcptn(:,:) - Lfus ) * zsnw(:,:) )           ! heat flux from snow (cell average) 
    16071606 
    16081607#else 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r6416 r7309  
    2424   USE ice             ! LIM-3: ice variables 
    2525   USE thd_ice         ! LIM-3: thermodynamical variables 
    26    USE dom_ice         ! LIM-3: ice domain 
    2726   ! 
    2827   USE sbc_oce         ! Surface boundary condition: ocean fields 
    … …  
    3130   USE sbcblk_clio     ! Surface boundary condition: CLIO bulk 
    3231   USE sbccpl          ! Surface boundary condition: coupled interface 
     32   USE sbcana          ! Surface boundary condition: analytic formulation 
    3333   USE albedo          ! ocean & ice albedo 
    3434   ! 
    … …  
    4848   USE limvar          ! Ice variables switch 
    4949   USE limctl          !  
    50    USE limmsh          ! LIM mesh 
    5150   USE limistate       ! LIM initial state 
    5251   USE limthd_sal      ! LIM ice thermodynamics: salinity 
    … …  
    6564   USE bdyice_lim       ! unstructured open boundary data  (bdy_ice_lim routine) 
    6665#endif 
     66# if defined key_agrif 
     67   USE agrif_ice 
     68   USE agrif_lim3_update 
     69   USE agrif_lim3_interp 
     70# endif 
    6771 
    6872   IMPLICIT NONE 
    … …  
    102106      !!--------------------------------------------------------------------- 
    103107      INTEGER, INTENT(in) ::   kt      ! ocean time step 
    104       INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE, =5 COUPLED) 
    105       !! 
    106       INTEGER  ::    jl                 ! dummy loop index 
     108      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=1 ANALYTIC, =3 CLIO, =4 CORE, =5 COUPLED) 
     109      !! 
     110      INTEGER  ::   jl                 ! dummy loop index 
    107111      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky 
    108112      REAL(wp), POINTER, DIMENSION(:,:  )   ::   zutau_ice, zvtau_ice  
    … …  
    110114 
    111115      IF( nn_timing == 1 )  CALL timing_start('sbc_ice_lim') 
     116 
     117      ! clem: it is important to initialize agrif_lim3 variables here and not in sbc_lim_init 
     118# if defined key_agrif 
     119      IF( kt == nit000 ) THEN 
     120         IF( .NOT. Agrif_Root() )   CALL Agrif_InitValues_cont_lim3 
     121      ENDIF 
     122# endif 
    112123 
    113124      !-----------------------! 
    … …  
    115126      !-----------------------! 
    116127      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
     128 
     129# if defined key_agrif 
     130         IF( .NOT. Agrif_Root() )  lim_nbstep = MOD( lim_nbstep, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) + 1 
     131# endif 
    117132 
    118133         ! mean surface ocean current at ice velocity point (C-grid dynamics :  U- & V-points as the ocean) 
    … …  
    136151         !----------------------------------------------------------------- 
    137152         SELECT CASE( kblk ) 
    138          CASE( jp_clio    )   ;   CALL blk_ice_clio_tau                         ! CLIO bulk formulation             
    139          CASE( jp_core    )   ;   CALL blk_ice_core_tau                         ! CORE bulk formulation 
    140          CASE( jp_purecpl )   ;   CALL sbc_cpl_ice_tau( utau_ice , vtau_ice )   ! Coupled   formulation 
     153            CASE( jp_ana     )   ;    CALL ana_ice_tau                              ! analytic formulation             
     154            CASE( jp_clio    )   ;    CALL blk_ice_clio_tau                         ! CLIO bulk formulation             
     155            CASE( jp_core    )   ;    CALL blk_ice_core_tau                         ! CORE bulk formulation 
     156            CASE( jp_purecpl )   ;    CALL sbc_cpl_ice_tau( utau_ice , vtau_ice )   ! Coupled   formulation 
    141157         END SELECT 
    142158          
    143          IF( ln_mixcpl) THEN   ! Case of a mixed Bulk/Coupled formulation 
     159         IF( ln_mixcpl) THEN                                                       ! Case of a mixed Bulk/Coupled formulation 
    144160            CALL wrk_alloc( jpi,jpj    , zutau_ice, zvtau_ice) 
    145             CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice ) 
     161                                      CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice ) 
    146162            utau_ice(:,:) = utau_ice(:,:) * xcplmask(:,:,0) + zutau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) 
    147163            vtau_ice(:,:) = vtau_ice(:,:) * xcplmask(:,:,0) + zvtau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) 
    … …  
    154170         numit = numit + nn_fsbc                  ! Ice model time step 
    155171         !                                                    
    156          CALL sbc_lim_bef                         ! Store previous ice values 
    157          CALL sbc_lim_diag0                       ! set diag of mass, heat and salt fluxes to 0 
    158          CALL lim_rst_opn( kt )                   ! Open Ice restart file 
    159          ! 
    160          IF( .NOT. lk_c1d ) THEN 
     172                                      CALL sbc_lim_bef         ! Store previous ice values 
     173                                      CALL sbc_lim_diag0       ! set diag of mass, heat and salt fluxes to 0 
     174                                      CALL lim_rst_opn( kt )   ! Open Ice restart file 
     175         ! 
     176         ! --- zap this if no ice dynamics --- ! 
     177         IF( .NOT. lk_c1d .AND. ln_limdyn ) THEN 
    161178            ! 
    162             CALL lim_dyn( kt )                    ! Ice dynamics    ( rheology/dynamics )    
    163             ! 
    164             CALL lim_trp( kt )                    ! Ice transport   ( Advection/diffusion ) 
    165             ! 
    166             IF( nn_monocat /= 2 ) CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) 
    167             ! 
    168 #if defined key_bdy 
    169             CALL bdy_ice_lim( kt )                ! bdy ice thermo  
    170             IF( ln_icectl )       CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 
    171 #endif 
    172             ! 
    173             CALL lim_update1( kt )                ! Corrections 
     179            IF( nn_limdyn /= 0 ) THEN                          ! -- Ice dynamics 
     180                                      CALL lim_dyn( kt )       !     rheology   
     181            ELSE 
     182               u_ice(:,:) = rn_uice * umask(:,:,1)             !     or prescribed velocity 
     183               v_ice(:,:) = rn_vice * vmask(:,:,1) 
     184            ENDIF 
     185                                      CALL lim_trp( kt )       ! -- Ice transport (Advection/diffusion) 
     186            IF( nn_limdyn == 2 .AND. nn_monocat /= 2 )  &      ! -- Mechanical redistribution (ridging/rafting) 
     187               &                      CALL lim_itd_me          
     188            IF( nn_limdyn == 2 )      CALL lim_update1( kt )   ! -- Corrections 
    174189            ! 
    175190         ENDIF 
    176           
     191 
     192         ! --- 
     193#if defined key_agrif 
     194         IF( .NOT. Agrif_Root() )     CALL agrif_interp_lim3('T') 
     195#endif 
     196#if defined key_bdy 
     197         IF( ln_limthd )              CALL bdy_ice_lim( kt )   ! -- bdy ice thermo  
     198#endif 
     199 
    177200         ! previous lead fraction and ice volume for flux calculations 
    178          CALL sbc_lim_bef                         
    179          CALL lim_var_glo2eqv                     ! ht_i and ht_s for ice albedo calculation 
    180          CALL lim_var_agg(1)                      ! at_i for coupling (via pfrld)  
     201                                      CALL sbc_lim_bef                         
     202                                      CALL lim_var_glo2eqv     ! ht_i and ht_s for ice albedo calculation 
     203                                      CALL lim_var_agg(1)      ! at_i for coupling (via pfrld)  
     204         ! 
    181205         pfrld(:,:)   = 1._wp - at_i(:,:) 
    182206         phicif(:,:)  = vt_i(:,:) 
    … …  
    193217         !---------------------------------------------------------------------------------------- 
    194218         CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs ) 
    195          CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 
    196  
     219          
     220                                      CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 
    197221         SELECT CASE( kblk ) 
    198          CASE( jp_clio )                                       ! CLIO bulk formulation 
    199             ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo  
    200             ! (alb_ice) is computed within the bulk routine 
    201                                  CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, alb_ice ) 
    202             IF( ln_mixcpl      ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
    203             IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
    204          CASE( jp_core )                                       ! CORE bulk formulation 
    205             ! albedo depends on cloud fraction because of non-linear spectral effects 
    206             alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
    207                                  CALL blk_ice_core_flx( t_su, alb_ice ) 
    208             IF( ln_mixcpl      ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
    209             IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
    210          CASE ( jp_purecpl ) 
    211             ! albedo depends on cloud fraction because of non-linear spectral effects 
    212             alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
    213                                  CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
    214             IF( nn_limflx == 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
     222 
     223            CASE( jp_ana )                                        ! analytic formulation 
     224                                      CALL ana_ice_flx 
     225                
     226            CASE( jp_clio )                                       ! CLIO bulk formulation 
     227               ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo  
     228               ! (alb_ice) is computed within the bulk routine 
     229                                      CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, alb_ice ) 
     230               IF( ln_mixcpl      )   CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
     231               IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
     232                
     233            CASE( jp_core )                                       ! CORE bulk formulation 
     234               ! albedo depends on cloud fraction because of non-linear spectral effects 
     235               alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     236                                      CALL blk_ice_core_flx( t_su, alb_ice ) 
     237               IF( ln_mixcpl      )   CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
     238               IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
     239                
     240            CASE ( jp_purecpl )                                    ! Coupled formulation 
     241               ! albedo depends on cloud fraction because of non-linear spectral effects 
     242               alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     243                                      CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 
     244               IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 
     245 
    215246         END SELECT 
     247 
    216248         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs ) 
    217249 
    … …  
    219251         ! --- ice thermodynamics --- ! 
    220252         !----------------------------! 
    221          CALL lim_thd( kt )                         ! Ice thermodynamics       
    222          ! 
    223          CALL lim_update2( kt )                     ! Corrections 
    224          ! 
    225          CALL lim_sbc_flx( kt )                     ! Update surface ocean mass, heat and salt fluxes 
    226          ! 
    227          IF(ln_limdiaout) CALL lim_diahsb           ! Diagnostics and outputs  
    228          ! 
    229          CALL lim_wri( 1 )                          ! Ice outputs  
     253         ! --- zap this if no ice thermo --- ! 
     254         IF( ln_limthd )              CALL lim_thd( kt )        ! -- Ice thermodynamics       
     255         IF( ln_limthd )              CALL lim_update2( kt )    ! -- Corrections 
     256         ! --- 
     257# if defined key_agrif 
     258         IF( .NOT. Agrif_Root() )     CALL agrif_update_lim3( kt ) 
     259# endif 
     260                                      CALL lim_var_glo2eqv      ! necessary calls (at least for coupling) 
     261                                      CALL lim_var_agg( 2 )     ! necessary calls (at least for coupling) 
     262                                      ! 
     263# if defined key_agrif 
     264!!         IF( .NOT. Agrif_Root() )     CALL Agrif_ChildGrid_To_ParentGrid()  ! clem: should be called at the update frequency only (cf agrif_lim3_update) 
     265# endif 
     266                                      CALL lim_sbc_flx( kt )    ! -- Update surface ocean mass, heat and salt fluxes 
     267# if defined key_agrif 
     268!!         IF( .NOT. Agrif_Root() )     CALL Agrif_ParentGrid_To_ChildGrid()  ! clem: should be called at the update frequency only (cf agrif_lim3_update) 
     269# endif 
     270         IF( ln_limdiahsb )           CALL lim_diahsb( kt )     ! -- Diagnostics and outputs  
     271         ! 
     272                                      CALL lim_wri( 1 )         ! -- Ice outputs  
    230273         ! 
    231274         IF( kt == nit000 .AND. ln_rstart )   & 
    232             &             CALL iom_close( numrir )  ! close input ice restart file 
    233          ! 
    234          IF( lrst_ice )   CALL lim_rst_write( kt )  ! Ice restart file  
    235          ! 
    236          IF( ln_icectl )  CALL lim_ctl( kt )        ! alerts in case of model crash 
     275            &                         CALL iom_close( numrir )  ! close input ice restart file 
     276         ! 
     277         IF( lrst_ice )               CALL lim_rst_write( kt )  ! -- Ice restart file  
     278         ! 
     279         IF( ln_limctl )              CALL lim_ctl( kt )        ! alerts in case of model crash 
    237280         ! 
    238281      ENDIF   ! End sea-ice time step only 
    … …  
    242285      !-------------------------! 
    243286      ! Update surface ocean stresses (only in ice-dynamic case) otherwise the atm.-ocean stresses are used everywhere 
    244       IF( ln_limdyn )     CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents 
     287      !    using before instantaneous surf. currents 
     288      IF( ln_limdyn )                 CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) ) 
    245289!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!! 
    246290      ! 
    … …  
    259303      !!---------------------------------------------------------------------- 
    260304      IF(lwp) WRITE(numout,*) 
    261       IF(lwp) WRITE(numout,*) 'sbc_ice_lim : update ocean surface boudary condition'  
     305      IF(lwp) WRITE(numout,*) 'sbc_lim_init : update ocean surface boudary condition'  
    262306      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   via Louvain la Neuve Ice Model (LIM-3) time stepping' 
    263307      ! 
    … …  
    271315      !                                ! Allocate the ice arrays 
    272316      ierr =        ice_alloc        ()      ! ice variables 
    273       ierr = ierr + dom_ice_alloc    ()      ! domain 
    274317      ierr = ierr + sbc_ice_alloc    ()      ! surface forcing 
    275318      ierr = ierr + thd_ice_alloc    ()      ! thermodynamics 
    276       ierr = ierr + lim_itd_me_alloc ()      ! ice thickness distribution - mechanics 
     319      IF( ln_limdyn )   ierr = ierr + lim_itd_me_alloc ()      ! ice thickness distribution - mechanics 
    277320      ! 
    278321      IF( lk_mpp    )   CALL mpp_sum( ierr ) 
    279322      IF( ierr /= 0 )   CALL ctl_stop('STOP', 'sbc_lim_init : unable to allocate ice arrays') 
    280323      ! 
    281       !                                ! adequation jpk versus ice/snow layers/categories 
    282       IF( jpl > jpk .OR. (nlay_i+1) > jpk .OR. nlay_s > jpk )   & 
    283          &      CALL ctl_stop( 'STOP',                          & 
    284          &     'sbc_lim_init: the 3rd dimension of workspace arrays is too small.',   & 
    285          &     'use more ocean levels or less ice/snow layers/categories.' ) 
     324      CALL lim_dyn_init                ! set ice dynamics parameters 
    286325      ! 
    287326      CALL lim_itd_init                ! ice thickness distribution initialization 
    … …  
    293332      CALL lim_thd_sal_init            ! set ice salinity parameters 
    294333      ! 
    295       CALL lim_msh                     ! ice mesh initialization 
    296       ! 
    297       CALL lim_itd_me_init             ! ice thickness distribution initialization for mecanical deformation 
     334      IF( ln_limdyn )   CALL lim_itd_me_init             ! ice thickness distribution initialization for mecanical deformation 
    298335      !                                ! Initial sea-ice state 
    299336      IF( .NOT. ln_rstart ) THEN              ! start from rest: sea-ice deduced from sst 
    … …  
    305342         numit = nit000 - 1 
    306343      ENDIF 
    307       CALL lim_var_agg(1) 
     344      CALL lim_var_agg(2) 
    308345      CALL lim_var_glo2eqv 
    309346      ! 
    310347      CALL lim_sbc_init                 ! ice surface boundary condition    
     348      ! 
     349      IF( ln_limdiahsb) CALL lim_diahsb_init  ! initialization for diags 
    311350      ! 
    312351      fr_i(:,:)     = at_i(:,:)         ! initialisation of sea-ice fraction 
    … …  
    342381      !!------------------------------------------------------------------- 
    343382      INTEGER  ::   ios                 ! Local integer output status for namelist read 
    344       NAMELIST/namicerun/ jpl, nlay_i, nlay_s, cn_icerst_in, cn_icerst_indir, cn_icerst_out, cn_icerst_outdir,  & 
    345          &                ln_limdyn, rn_amax_n, rn_amax_s, ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt   
     383      NAMELIST/namicerun/ jpl, nlay_i, nlay_s, rn_amax_n, rn_amax_s, cn_icerst_in, cn_icerst_indir,   & 
     384         &                cn_icerst_out, cn_icerst_outdir, ln_limthd, ln_limdyn, nn_limdyn, rn_uice, rn_vice   
     385      NAMELIST/namicediag/ ln_limdiachk, ln_limdiahsb, ln_limctl, iiceprt, jiceprt   
    346386      !!------------------------------------------------------------------- 
    347387      !                     
    348388      REWIND( numnam_ice_ref )              ! Namelist namicerun in reference namelist : Parameters for ice 
    349389      READ  ( numnam_ice_ref, namicerun, IOSTAT = ios, ERR = 901) 
    350 901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namicerun in reference namelist', lwp ) 
    351       ! 
     390901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicerun in reference namelist', lwp ) 
     391 
    352392      REWIND( numnam_ice_cfg )              ! Namelist namicerun in configuration namelist : Parameters for ice 
    353393      READ  ( numnam_ice_cfg, namicerun, IOSTAT = ios, ERR = 902 ) 
    354 902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namicerun in configuration namelist', lwp ) 
    355       IF(lwm) WRITE( numoni, namicerun ) 
     394902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicerun in configuration namelist', lwp ) 
     395      IF(lwm) WRITE ( numoni, namicerun ) 
     396      ! 
     397      REWIND( numnam_ice_ref )              ! Namelist namicediag in reference namelist : Parameters for ice 
     398      READ  ( numnam_ice_ref, namicediag, IOSTAT = ios, ERR = 903) 
     399903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicediag in reference namelist', lwp ) 
     400 
     401      REWIND( numnam_ice_cfg )              ! Namelist namicediag in configuration namelist : Parameters for ice 
     402      READ  ( numnam_ice_cfg, namicediag, IOSTAT = ios, ERR = 904 ) 
     403904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicediag in configuration namelist', lwp ) 
     404      IF(lwm) WRITE ( numoni, namicediag ) 
    356405      ! 
    357406      IF(lwp) THEN                        ! control print 
    … …  
    362411         WRITE(numout,*) '   number of ice  layers                                   = ', nlay_i 
    363412         WRITE(numout,*) '   number of snow layers                                   = ', nlay_s 
    364          WRITE(numout,*) '   switch for ice dynamics (1) or not (0)      ln_limdyn   = ', ln_limdyn 
    365413         WRITE(numout,*) '   maximum ice concentration for NH                        = ', rn_amax_n  
    366414         WRITE(numout,*) '   maximum ice concentration for SH                        = ', rn_amax_s 
    367          WRITE(numout,*) '   Diagnose heat/salt budget or not          ln_limdiahsb  = ', ln_limdiahsb 
    368          WRITE(numout,*) '   Output   heat/salt budget or not          ln_limdiaout  = ', ln_limdiaout 
    369          WRITE(numout,*) '   control prints in ocean.out for (i,j)=(iiceprt,jiceprt) = ', ln_icectl 
    370          WRITE(numout,*) '   i-index for control prints (ln_icectl=true)             = ', iiceprt 
    371          WRITE(numout,*) '   j-index for control prints (ln_icectl=true)             = ', jiceprt 
     415         WRITE(numout,*) '   Ice thermodynamics (T) or not (F)            ln_limthd  = ', ln_limthd 
     416         WRITE(numout,*) '   Ice dynamics       (T) or not (F)            ln_limdyn  = ', ln_limdyn 
     417         WRITE(numout,*) '     (ln_limdyn=T) Ice dynamics switch          nn_limdyn  = ', nn_limdyn 
     418         WRITE(numout,*) '       2: total' 
     419         WRITE(numout,*) '       1: advection only (no diffusion, no ridging/rafting)' 
     420         WRITE(numout,*) '       0: advection only (as 1 + prescribed velocity, bypass rheology)' 
     421         WRITE(numout,*) '     (ln_limdyn=T) prescribed u-vel (case nn_limdyn=0)     = ', rn_uice 
     422         WRITE(numout,*) '     (ln_limdyn=T) prescribed v-vel (case nn_limdyn=0)     = ', rn_vice 
     423         WRITE(numout,*) 
     424         WRITE(numout,*) '...and ice diagnostics' 
     425         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~' 
     426         WRITE(numout,*) '   Diagnose online heat/mass/salt budget     ln_limdiachk  = ', ln_limdiachk 
     427         WRITE(numout,*) '   Output          heat/mass/salt budget     ln_limdiahsb  = ', ln_limdiahsb 
     428         WRITE(numout,*) '   control prints in ocean.out for (i,j)=(iiceprt,jiceprt) = ', ln_limctl 
     429         WRITE(numout,*) '   i-index for control prints (ln_limctl=true)             = ', iiceprt 
     430         WRITE(numout,*) '   j-index for control prints (ln_limctl=true)             = ', jiceprt 
    372431      ENDIF 
    373432      ! 
    374433      ! sea-ice timestep and inverse 
    375       rdt_ice   = nn_fsbc * rdt   
     434      rdt_ice   = REAL(nn_fsbc) * rdt   
    376435      r1_rdtice = 1._wp / rdt_ice  
    377436 
    … …  
    381440      ! 
    382441#if defined key_bdy 
    383       IF( lwp .AND. ln_limdiahsb )  CALL ctl_warn('online conservation check activated but it does not work with BDY') 
     442      IF( lwp .AND. ln_limdiachk )  CALL ctl_warn('online conservation check activated but it does not work with BDY') 
    384443#endif 
     444      ! 
     445      IF( lwp ) WRITE(numout,*) '   ice timestep rdt_ice  = ', rdt_ice 
    385446      ! 
    386447   END SUBROUTINE ice_run 
    … …  
    404465      ! 
    405466      REWIND( numnam_ice_ref )              ! Namelist namiceitd in reference namelist : Parameters for ice 
    406       READ  ( numnam_ice_ref, namiceitd, IOSTAT = ios, ERR = 903) 
    407 903   IF( ios /= 0 )  CALL ctl_nam ( ios , 'namiceitd in reference namelist', lwp ) 
    408       ! 
     467      READ  ( numnam_ice_ref, namiceitd, IOSTAT = ios, ERR = 905) 
     468905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in reference namelist', lwp ) 
     469 
    409470      REWIND( numnam_ice_cfg )              ! Namelist namiceitd in configuration namelist : Parameters for ice 
    410       READ  ( numnam_ice_cfg, namiceitd, IOSTAT = ios, ERR = 904 ) 
    411 904   IF( ios /= 0 )  CALL ctl_nam ( ios , 'namiceitd in configuration namelist', lwp ) 
    412       IF(lwm) WRITE( numoni, namiceitd ) 
     471      READ  ( numnam_ice_cfg, namiceitd, IOSTAT = ios, ERR = 906 ) 
     472906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in configuration namelist', lwp ) 
     473      IF(lwm) WRITE ( numoni, namiceitd ) 
    413474      ! 
    414475      IF(lwp) THEN                        ! control print 
    415476         WRITE(numout,*) 
    416          WRITE(numout,*) 'ice_itd : ice cat distribution' 
    417          WRITE(numout,*) ' ~~~~~~' 
    418          WRITE(numout,*) '   shape of ice categories distribution                     nn_catbnd = ', nn_catbnd 
    419          WRITE(numout,*) '   mean ice thickness in the domain (used if nn_catbnd=2)  rn_himean = ', rn_himean 
     477         WRITE(numout,*) 'lim_itd_init : Initialization of ice cat distribution ' 
     478         WRITE(numout,*) '~~~~~~~~~~~~' 
     479         WRITE(numout,*) '   shape of ice categories distribution                          nn_catbnd = ', nn_catbnd 
     480         WRITE(numout,*) '   mean ice thickness in the domain (only active if nn_catbnd=2) rn_himean = ', rn_himean 
    420481      ENDIF 
    421482      ! 
    … …  
    423484      !- Thickness categories boundaries  
    424485      !---------------------------------- 
    425       IF(lwp) WRITE(numout,*) 
    426       IF(lwp) WRITE(numout,*) 'lim_itd_init : Initialization of ice cat distribution ' 
    427       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 
    428       ! 
    429486      hi_max(:) = 0._wp 
    430487      ! 
    … …  
    463520 
    464521    
    465    SUBROUTINE ice_lim_flx( ptn_ice , palb_ice, pqns_ice ,    & 
    466       &                    pqsr_ice, pdqn_ice, pevap_ice, pdevap_ice, k_limflx ) 
     522   SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice, pdqn_ice, pevap_ice, pdevap_ice, k_limflx ) 
    467523      !!--------------------------------------------------------------------- 
    468524      !!                  ***  ROUTINE ice_lim_flx  *** 
    … …  
    557613      u_ice_b(:,:)     = u_ice(:,:) 
    558614      v_ice_b(:,:)     = v_ice(:,:) 
    559       !       
     615      at_i_b (:,:)     = SUM( a_i_b(:,:,:), dim=3 ) 
     616       
    560617   END SUBROUTINE sbc_lim_bef 
    561618 
    … …  
    569626      !!---------------------------------------------------------------------- 
    570627      sfx    (:,:) = 0._wp   ; 
    571       sfx_bri(:,:) = 0._wp   ;  
     628      sfx_bri(:,:) = 0._wp   ;   sfx_lam(:,:) = 0._wp 
    572629      sfx_sni(:,:) = 0._wp   ;   sfx_opw(:,:) = 0._wp 
    573630      sfx_bog(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp 
    … …  
    580637      wfx_bom(:,:) = 0._wp   ;   wfx_sum(:,:) = 0._wp 
    581638      wfx_res(:,:) = 0._wp   ;   wfx_sub(:,:) = 0._wp 
    582       wfx_spr(:,:) = 0._wp   ;    
    583       ! 
     639      wfx_spr(:,:) = 0._wp   ;   wfx_lam(:,:) = 0._wp   
     640       
    584641      hfx_thd(:,:) = 0._wp   ;    
    585642      hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp 
    … …  
    597654      diag_heat(:,:) = 0._wp ;   diag_smvi(:,:) = 0._wp ; 
    598655      diag_vice(:,:) = 0._wp ;   diag_vsnw(:,:) = 0._wp ; 
    599       ! 
     656 
     657      tau_icebfr(:,:) = 0._wp; ! landfast ice param only (clem: important to keep the init here) 
     658       
    600659   END SUBROUTINE sbc_lim_diag0 
    601660 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r6460 r7309  
    115115      ! 
    116116      !                          ! overwrite namelist parameter using CPP key information 
    117       IF( Agrif_Root() ) THEN                ! AGRIF zoom 
    118         IF( lk_lim2 )   nn_ice      = 2 
    119         IF( lk_lim3 )   nn_ice      = 3 
    120         IF( lk_cice )   nn_ice      = 4 
    121       ENDIF 
    122       IF( cp_cfg == 'gyre' ) THEN            ! GYRE configuration 
    123           ln_ana      = .TRUE.    
    124           nn_ice      =   0 
    125       ENDIF 
    126       ! 
     117#if defined key_agrif 
     118      IF( Agrif_Root() ) THEN                ! AGRIF zoom (cf r1242: possibility to run without ice in fine grid) 
     119         IF( lk_lim2 )   nn_ice      = 2 
     120         IF( lk_lim3 )   nn_ice      = 3 
     121         IF( lk_cice )   nn_ice      = 4 
     122      ENDIF 
     123#else 
     124      IF( lk_lim2 )   nn_ice      = 2 
     125      IF( lk_lim3 )   nn_ice      = 3 
     126      IF( lk_cice )   nn_ice      = 4      
     127#endif 
     128 
     129      IF( cp_cfg == 'gyre' )   ln_ana = .TRUE.          ! GYRE configuration 
     130              
    127131      IF(lwp) THEN               ! Control print 
    128132         WRITE(numout,*) '        Namelist namsbc (partly overwritten with CPP key setting)' 
    … …  
    200204 
    201205      !                                            ! restartability    
    202       IF( ( nn_ice == 2 .OR. nn_ice ==3 ) .AND. .NOT.( ln_blk_clio .OR. ln_blk_core .OR. ln_cpl ) )   & 
    203          &   CALL ctl_stop( 'LIM sea-ice model requires a bulk formulation or coupled configuration' ) 
    204206      IF( nn_ice == 4 .AND. .NOT.( ln_blk_core .OR. ln_cpl ) )   & 
    205207         &   CALL ctl_stop( 'CICE sea-ice model requires ln_blk_core or ln_cpl' ) 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfric.F90

    r7048 r7309  
    205205         DO jj = 2, jpjm1 
    206206            DO ji = fs_2, fs_jpim1 
    207                IF( fsdept(ji,jj,jk) < ekm_dep(ji,jj) ) THEN 
     207               IF( gdept_n(ji,jj,jk) < ekm_dep(ji,jj) ) THEN 
    208208                  avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk), rn_wvmix ) 
    209209                  avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk), rn_wvmix ) 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/lib_fortran.F90

    r6140 r7309  
    2525 
    2626   PUBLIC   glob_sum      ! used in many places (masked with tmask_i) 
    27    PUBLIC   glob_sum_full ! used in many places (masked with tmask_h, ie omly over the halos) 
     27   PUBLIC   glob_sum_full ! used in many places (masked with tmask_h, ie only over the halos) 
    2828   PUBLIC   DDPDD         ! also used in closea module 
    2929   PUBLIC   glob_min, glob_max 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r6152 r7309  
    104104 
    105105   !!---------------------------------------------------------------------- 
    106    !! NEMO/OPA 3.7 , NEMO Consortium (2015) 
     106   !! NEMO/OPA 3.7 , NEMO Consortium (2016) 
    107107   !! $Id$ 
    108108   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    141141# endif 
    142142# if defined key_lim2 
    143       CALL Agrif_Declare_Var_lim2  !  "      "   "   "      "  LIM 
     143      CALL Agrif_Declare_Var_lim2  !  "      "   "   "      "  LIM2 
     144# endif 
     145# if defined key_lim3 
     146      CALL Agrif_Declare_Var_lim3  !  "      "   "   "      "  LIM3 
    144147# endif 
    145148#endif 

  • branches/2016/dev_CNRS_AGRIF_2016/NEMOGCM/NEMO/OPA_SRC/step.F90

    r6464 r7309  
    295295      IF(.NOT.ln_linssh) CALL dom_vvl_sf_swp( kstp )  ! swap of vertical scale factors 
    296296      ! 
     297      IF( ln_diahsb        )   CALL dia_hsb( kstp )         ! - ML - global conservation diagnostics 
    297298 
    298299!!gm : This does not only concern the dynamics ==>>> add a new title 
    … …  
    316317      ENDIF 
    317318#endif 
    318       IF( ln_diahsb  )   CALL dia_hsb( kstp )         ! - ML - global conservation diagnostics 
    319319      IF( ln_diaobs  )   CALL dia_obs( kstp )         ! obs-minus-model (assimilation) diagnostics (call after dynamics update) 
    320320 
Note: See TracChangeset for help on using the changeset viewer.