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

Changeset 8637


Ignore:
Timestamp:
2017-10-18T19:14:32+02:00 (7 years ago)
Author:
gm
Message:

#1911 (ENHANCE-09): PART I.3 - phasing with updated branch dev_r8183_ICEMODEL revision 8626

Location:
branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM
Files:
3 added
4 deleted
63 edited

Legend:

Unmodified
Added
Removed

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/EXP00/1_namelist_cfg

    r8586 r8637  
    4343                           !  =1 use observed ice-cover      , 
    4444                           !  =2 or 3 automatically for LIM3 or CICE    ("key_lim3" or "key_cice") 
     45                           !          except in AGRIF zoom where it has to be specified 
    4546   ln_rnf      = .false.   !  runoffs                                   (T => fill namsbc_rnf) 
    4647   ln_ssr      = .false.   !  Sea Surface Restoring on T and/or S       (T => fill namsbc_ssr) 
    … …  
    126127!----------------------------------------------------------------------- 
    127128   !                       !  Operator type: 
     129   ln_traldf_NONE  =  .false.  !  No explicit diffusion 
    128130   ln_traldf_lap   =  .true.   !    laplacian operator 
    129131   ln_traldf_blp   =  .false.  !  bilaplacian operator 
    130132   !                       !  Direction of action: 
    131133   ln_traldf_lev   =  .false.  !  iso-level 
    132    ln_traldf_hor   =  .false.  !  horizontal (geopotential) 
     134   ln_traldf_hor   =  .false.  !  horizontal  (geopotential) 
    133135   ln_traldf_iso   =  .true.   !  iso-neutral (Standard operator) 
    134136   ln_traldf_triad =  .false.  !  iso-neutral (Triads   operator) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/EXP00/file_def_nemo.xml

    r8215 r8637  
    350350   </file_definition> 
    351351     
    352    
    353 </simulation> 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/EXP00/namelist_ice_cfg

    r8586 r8637  
    1313!!             11 - Ice growth in open water           (namthd_do) 
    1414!!             12 - Ice salinity                       (namthd_sal) 
    15 !!             13 - Ice melt ponds                     (nammp) 
     15!!             13 - Ice melt ponds                     (namthd_pnd) 
    1616!!             14 - Ice initialization                 (namini) 
    1717!!             15 - Ice/snow albedos                   (namalb) 
    … …  
    2020! 
    2121!------------------------------------------------------------------------------ 
    22 &nampar     !   Generic parameters 
     22&nampar         !   Generic parameters 
    2323!------------------------------------------------------------------------------ 
    2424/ 
    2525!------------------------------------------------------------------------------ 
    26 &namitd     !   Ice discretization 
     26&namitd         !   Ice discretization 
    2727!------------------------------------------------------------------------------ 
    2828/ 
    2929!------------------------------------------------------------------------------ 
    30 &namdyn     !   Ice dynamics 
     30&namdyn         !   Ice dynamics 
    3131!------------------------------------------------------------------------------ 
    3232/ 
    … …  
    4848/ 
    4949!------------------------------------------------------------------------------ 
    50 &namthd     !   Ice thermodynamics 
     50&namthd         !   Ice thermodynamics 
    5151!------------------------------------------------------------------------------ 
    5252/ 
    … …  
    5656/ 
    5757!------------------------------------------------------------------------------ 
    58 &namthd_da     !   Ice lateral melting 
     58&namthd_da      !   Ice lateral melting 
    5959!------------------------------------------------------------------------------ 
    6060/ 
    6161!------------------------------------------------------------------------------ 
    62 &namthd_do     !   Ice growth in open water 
     62&namthd_do      !   Ice growth in open water 
    6363!------------------------------------------------------------------------------ 
    6464/ 
    … …  
    6868/ 
    6969!------------------------------------------------------------------------------ 
    70 &nammp      !   Melt ponds 
     70&namthd_pnd     !   Melt ponds 
    7171!------------------------------------------------------------------------------ 
    7272/ 
    7373!------------------------------------------------------------------------------ 
    74 &namini     !   Ice initialization 
     74&namini         !   Ice initialization 
    7575!------------------------------------------------------------------------------ 
    7676/ 
    7777!------------------------------------------------------------------------------ 
    78 &namalb     !   albedo parameters 
     78&namalb         !   albedo parameters 
    7979!------------------------------------------------------------------------------ 
    8080/ 
    8181!------------------------------------------------------------------------------ 
    82 &namdia     !   Diagnostics 
     82&namdia         !   Diagnostics 
    8383!------------------------------------------------------------------------------ 
    8484/ 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_LIM3_PISCES/cpp_ORCA2_LIM3_PISCES.fcm

    r8215 r8637  
    1 bld::tool::fppkeys   key_lim3 key_top key_iomput key_mpp_mpi key_nosignedzero 
     1bld::tool::fppkeys   key_lim3 key_iomput key_mpp_mpi key_nosignedzero key_agrif 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_SAS_LIM3/EXP00/namelist_cfg

    r8568 r8637  
    7676&namtra_adv    !   advection scheme for tracer                          (default: NO selection) 
    7777!----------------------------------------------------------------------- 
     78   ln_traadv_NONE =  .true.    !  No advection 
     79/ 
     80!---------------------------------------------------------------------------------- 
     81&namtra_ldf    !   lateral diffusion scheme for tracers 
     82!---------------------------------------------------------------------------------- 
     83   ln_traldf_lap  =  .NONE.   ! No explicit diffusion 
     84/ 
     85!----------------------------------------------------------------------- 
     86&namdyn_adv    !   formulation of the momentum advection 
     87!----------------------------------------------------------------------- 
     88   ln_dynadv_NONE =  .true.    !  No advection 
    7889/ 
    7990!----------------------------------------------------------------------- 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/ORCA2_SAS_LIM3/EXP00/namelist_ice_cfg

    r8586 r8637  
    11!!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    22!! ESIM configuration namelist: Overwrites SHARED/namelist_ice_lim3_ref 
    3 !!              1 - Generic parameters                 (namice_run) 
    4 !!              2 - Ice thickness discretization       (namice_itd) 
    5 !!              3 - Ice dynamics                       (namice_dyn) 
    6 !!              4 - Ice ridging/rafting                (namice_rdgrft) 
    7 !!              5 - Ice rheology                       (namice_rhg) 
    8 !!              6 - Ice advection                      (namice_adv) 
    9 !!              7 - Ice thermodynamics                 (namice_thd) 
    10 !!              8 - Ice salinity                       (namice_sal) 
    11 !!              9 - Ice melt ponds                     (namice_mp) 
    12 !!             10 - Ice initialization                 (namice_ini) 
    13 !!             11 - Ice/snow albedos                   (namice_alb) 
    14 !!             12 - Ice diagnostics                    (namice_dia) 
     3!!              1 - Generic parameters                 (nampar) 
     4!!              2 - Ice thickness discretization       (namitd) 
     5!!              3 - Ice dynamics                       (namdyn) 
     6!!              4 - Ice ridging/rafting                (namdyn_rdgrft) 
     7!!              5 - Ice rheology                       (namdyn_rhg) 
     8!!              6 - Ice advection                      (namdyn_adv) 
     9!!              7 - Ice surface forcing                (namforcing) 
     10!!              8 - Ice thermodynamics                 (namthd) 
     11!!              9 - Ice heat diffusion                 (namthd_zdf) 
     12!!             10 - Ice lateral melting                (namthd_da) 
     13!!             11 - Ice growth in open water           (namthd_do) 
     14!!             12 - Ice salinity                       (namthd_sal) 
     15!!             13 - Ice melt ponds                     (namthd_pnd) 
     16!!             14 - Ice initialization                 (namini) 
     17!!             15 - Ice/snow albedos                   (namalb) 
     18!!             16 - Ice diagnostics                    (namdia) 
    1519!!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    1620! 
    1721!------------------------------------------------------------------------------ 
    18 &namice_run     !   Generic parameters 
     22&nampar         !   Generic parameters 
    1923!------------------------------------------------------------------------------ 
    2024/ 
    2125!------------------------------------------------------------------------------ 
    22 &namice_itd     !   Ice discretization 
     26&namitd         !   Ice discretization 
    2327!------------------------------------------------------------------------------ 
    2428/ 
    2529!------------------------------------------------------------------------------ 
    26 &namice_dyn     !   Ice dynamics 
     30&namdyn         !   Ice dynamics 
    2731!------------------------------------------------------------------------------ 
    2832/ 
    2933!------------------------------------------------------------------------------ 
    30 &namice_rdgrft  !   Ice ridging/rafting 
     34&namdyn_rdgrft  !   Ice ridging/rafting 
    3135!------------------------------------------------------------------------------ 
    3236/ 
    3337!------------------------------------------------------------------------------ 
    34 &namice_rhg     !   Ice rheology 
     38&namdyn_rhg     !   Ice rheology 
    3539!------------------------------------------------------------------------------ 
    3640/ 
    3741!------------------------------------------------------------------------------ 
    38 &namice_adv     !   Ice advection 
     42&namdyn_adv     !   Ice advection 
    3943!------------------------------------------------------------------------------ 
    4044/ 
    4145!------------------------------------------------------------------------------ 
    42 &namice_thd     !   Ice thermodynamics 
     46&namforcing     !   Ice surface forcing 
    4347!------------------------------------------------------------------------------ 
    4448/ 
    4549!------------------------------------------------------------------------------ 
    46 &namice_sal     !   Ice salinity 
     50&namthd         !   Ice thermodynamics 
    4751!------------------------------------------------------------------------------ 
    4852/ 
    4953!------------------------------------------------------------------------------ 
    50 &namicemp      !   Melt ponds 
     54&namthd_zdf     !   Ice heat diffusion 
    5155!------------------------------------------------------------------------------ 
    5256/ 
    5357!------------------------------------------------------------------------------ 
    54 &namice_ini     !   Ice initialization 
     58&namthd_da      !   Ice lateral melting 
    5559!------------------------------------------------------------------------------ 
    5660/ 
    5761!------------------------------------------------------------------------------ 
    58 &namice_alb     !   albedo parameters 
     62&namthd_do      !   Ice growth in open water 
    5963!------------------------------------------------------------------------------ 
    6064/ 
    6165!------------------------------------------------------------------------------ 
    62 &namice_dia     !   Diagnostics 
     66&namthd_sal     !   Ice salinity 
    6367!------------------------------------------------------------------------------ 
    6468/ 
     69!------------------------------------------------------------------------------ 
     70&namthd_pnd     !   Melt ponds 
     71!------------------------------------------------------------------------------ 
     72/ 
     73!------------------------------------------------------------------------------ 
     74&namini         !   Ice initialization 
     75!------------------------------------------------------------------------------ 
     76/ 
     77!------------------------------------------------------------------------------ 
     78&namalb         !   albedo parameters 
     79!------------------------------------------------------------------------------ 
     80/ 
     81!------------------------------------------------------------------------------ 
     82&namdia         !   Diagnostics 
     83!------------------------------------------------------------------------------ 
     84/ 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/SHARED/field_def_nemo-lim.xml

    r8586 r8637  
    44    <!--  
    55============================================================================================================ 
    6 =                                  definition of all existing variables                                    = 
     6=                                  definition of all existing ESIM variables                                    = 
    77=                                            DO NOT CHANGE                                                 = 
    88============================================================================================================ 
    … …  
    132132         <field id="icevmp"       long_name="melt pond volume"                                             standard_name="sea_ice_meltpond_volume"                            unit="m"            />  
    133133 
    134          <field id="iceconc_cat"  long_name="Sea-ice area fractions in thickness categories"               unit=""   grid_ref="grid_T_3D_ncatice" /> 
     134         <field id="iceconc_cat"  long_name="Sea-ice concentration in thickness categories"                unit=""   grid_ref="grid_T_3D_ncatice" /> 
    135135         <field id="icethic_cat"  long_name="Sea-ice thickness in thickness categories"                    unit="m"  grid_ref="grid_T_3D_ncatice" /> 
    136136         <field id="snowthic_cat" long_name="Snow thickness in thickness categories"                       unit="m"  grid_ref="grid_T_3D_ncatice" /> 
    … …  
    140140         <field id="icetemp_cat"  long_name="Ice temperature for categories"                               unit="degC"   grid_ref="grid_T_3D_ncatice" /> 
    141141         <field id="snwtemp_cat"  long_name="Snow temperature for categories"                              unit="degC"   grid_ref="grid_T_3D_ncatice" /> 
    142          <field id="iceamp_cat"   long_name="Ice melt pond fraction for categories"                        unit="%"      grid_ref="grid_T_3D_ncatice" />  
     142         <field id="iceamp_cat"   long_name="Ice melt pond concentration for categories"                   unit="%"      grid_ref="grid_T_3D_ncatice" />  
    143143         <field id="icevmp_cat"   long_name="Ice melt pond volume for categories"                          unit="m"      grid_ref="grid_T_3D_ncatice" />  
     144         <field id="icehmp_cat"   long_name="Ice melt pond thickness for categories"                       unit="m"      grid_ref="grid_T_3D_ncatice" />  
     145         <field id="iceafp_cat"   long_name="Ice melt pond fraction for categories"                        unit="m"      grid_ref="grid_T_3D_ncatice" />  
    144146 
    145147         <field id="micet"        long_name="Mean ice temperature"                                         unit="degC"     /> 
    … …  
    186188         <field id="vfxspr"       long_name="snw precipitation on ice"                                     unit="kg/m2/s"   /> 
    187189         <field id="vfxthin"      long_name="thermo ice prod. for thin ice(20cm) + open water"             unit="kg/m2/s"   /> 
     190         <field id="vfxpnd"       long_name="melt pond water flux to ocean"                                unit="kg/m2/s"   /> 
    188191 
    189192         <field id="afxtot"       long_name="area tendency (total)"                                        unit="s-1"   /> 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/SHARED/field_def_nemo-opa.xml

    r8586 r8637  
    298298         <field id="uoces"        long_name="ocean transport along i-axis times salinity (CRS)"                                                  unit="1e-3*m/s"   grid_ref="grid_U_3D" /> 
    299299 
    300          <!-- u-eddy coefficients (ldftra) --> 
     300         <!-- u-eddy diffusivity coefficients (available if ln_traldf_NONE=F) --> 
    301301         <field id="ahtu_2d"      long_name=" surface u-eddy diffusivity coefficient"   unit="m2/s or m4/s" /> 
    302302         <field id="ahtu_3d"      long_name=" 3D u-EIV coefficient"                     unit="m2/s or m4/s"      grid_ref="grid_U_3D"/> 
     303         <!-- v-eiv diffusivity coefficients (available if ln_ldfeiv=F) --> 
    303304         <field id="aeiu_2d"      long_name=" surface u-EIV coefficient"                unit="m2/s" /> 
    304305         <field id="aeiu_3d"      long_name=" 3D u-EIV coefficient"                     unit="m2/s"              grid_ref="grid_U_3D"/> 
    … …  
    307308         <field id="psiu_mle"     long_name="MLE streamfunction along i-axis"   unit="m3/s"   grid_ref="grid_U_3D" /> 
    308309 
    309          <!-- uoce_eiv: available EIV --> 
     310         <!-- uoce_eiv: available EIV (ln_ldfeiv=T and ln_ldfeiv_dia=T) --> 
    310311         <field id="uoce_eiv"     long_name="EIV ocean current along i-axis"   standard_name="bolus_sea_water_x_velocity"   unit="m/s"   grid_ref="grid_U_3D" /> 
    311312 
    312          <!-- uoce_eiv: available with key_trabbl --> 
     313         <!-- uoce_bbl: available with ln_trabbl=T and nn_bbl_adv=1 --> 
    313314         <field id="uoce_bbl"     long_name="BBL ocean current along i-axis"    unit="m/s"  /> 
     315         <!-- ahu_bbl : available with ln_trabbl=T and nn_bbl_ldf=1 --> 
    314316         <field id="ahu_bbl"      long_name="BBL diffusive flux along i-axis"   unit="m3/s" /> 
    315317 
    … …  
    346348         <field id="voces"        long_name="ocean transport along j-axis times salinity (CRS)"                                                  unit="1e-3*m/s"   grid_ref="grid_V_3D" /> 
    347349 
    348          <!-- v-eddy coefficients (ldftra, ldfdyn) --> 
     350         <!-- v-eddy diffusivity coefficients (available if ln_traldf_NONE=F) --> 
    349351         <field id="ahtv_2d"      long_name=" surface v-eddy diffusivity coefficient"     unit="m2/s or (m4/s)^1/2" /> 
    350352         <field id="ahtv_3d"      long_name=" 3D v-eddy diffusivity coefficient"          unit="m2/s or (m4/s)^1/2"           grid_ref="grid_V_3D"/> 
     353         <!-- v-eiv diffusivity coefficients (available if ln_ldfeiv=F) --> 
    351354         <field id="aeiv_2d"      long_name=" surface v-EIV coefficient"                  unit="m2/s" /> 
    352355         <field id="aeiv_3d"      long_name=" 3D v-EIV coefficient"                       unit="m2/s"                         grid_ref="grid_V_3D" /> 
    353356 
    354          <!-- variables available with MLE --> 
     357         <!-- variables available with MLE (ln_mle=T) --> 
    355358         <field id="psiv_mle"     long_name="MLE streamfunction along j-axis"   unit="m3/s"   grid_ref="grid_V_3D" /> 
    356359 
    357          <!-- voce_eiv: available with EIV --> 
     360         <!-- voce_eiv: available EIV (ln_ldfeiv=T and ln_ldfeiv_dia=T) --> 
    358361         <field id="voce_eiv"     long_name="EIV ocean current along j-axis"   standard_name="bolus_sea_water_y_velocity"   unit="m/s"   grid_ref="grid_V_3D" /> 
    359362 
    360          <!-- voce_eiv: available with key_trabbl --> 
     363         <!-- voce_bbl: available with ln_trabbl=T and nn_bbl_adv=1 --> 
    361364         <field id="voce_bbl"     long_name="BBL ocean current along j-axis"    unit="m/s"  /> 
     365         <!-- ahu_bbl : available with ln_trabbl=T and nn_bbl_ldf=1 --> 
    362366         <field id="ahv_bbl"      long_name="BBL diffusive flux along j-axis"   unit="m3/s" /> 
    363367 
    … …  
    386390         <field id="wocetr_eff"   long_name="effective ocean vertical transport"                                               unit="m3/s" /> 
    387391 
    388          <!-- woce_eiv: available with EIV --> 
     392         <!-- woce_eiv: available with EIV  (ln_ldfeiv=T and ln_ldfeiv_dia=T)  --> 
    389393         <field id="woce_eiv"     long_name="EIV ocean vertical velocity"   standard_name="bolus_upward_sea_water_velocity"   unit="m/s" /> 
    390394 
    … …  
    393397         <field id="avm"          long_name="vertical eddy viscosity"     standard_name="ocean_vertical_momentum_diffusivity"   unit="m2/s" /> 
    394398 
    395          <!-- avs: available with key_zdfddm --> 
     399         <!-- avs: /= avt with ln_zdfddm=T --> 
    396400         <field id="avs"          long_name="salt vertical eddy diffusivity"   standard_name="ocean_vertical_salt_diffusivity"   unit="m2/s" /> 
    397401         <field id="logavs"       long_name="logarithm of salt vertical eddy diffusivity"   standard_name="ocean_vertical_heat_diffusivity"       unit="m2/s" /> 
    … …  
    401405         <field id="avm_evd"      long_name="convective enhancement of vertical viscosity"     standard_name="ocean_vertical_momentum_diffusivity_due_to_convection"   unit="m2/s" /> 
    402406 
    403          <!-- avt_tide: available with key_zdftmx --> 
    404          <field id="av_tide"      long_name="tidal vertical diffusivity"   standard_name="ocean_vertical_tracer_diffusivity_due_to_tides"   unit="m2/s" /> 
    405  
    406          <!-- variables available with key_zdftmx_new --> 
    407          <field id="av_ratio"     long_name="S over T diffusivity ratio"            standard_name="salinity_over_temperature_diffusivity_ratio"                     unit="1"    /> 
    408          <field id="av_wave"      long_name="wave-induced vertical diffusivity"     standard_name="ocean_vertical_tracer_diffusivity_due_to_internal_waves"         unit="m2/s" /> 
    409          <field id="bflx_tmx"     long_name="wave-induced buoyancy flux"            standard_name="buoyancy_flux_due_to_internal_waves"                             unit="W/kg" /> 
    410          <field id="pcmap_tmx"    long_name="power consumed by wave-driven mixing"  standard_name="vertically_integrated_power_consumption_by_wave_driven_mixing"   unit="W/m2"      grid_ref="grid_W_2D" /> 
    411          <field id="emix_tmx"     long_name="power density available for mixing"    standard_name="power_available_for_mixing_from_breaking_internal_waves"         unit="W/kg" /> 
     407         <!-- avt_tide: available with ln_zdfiwm=T --> 
     408         <field id="av_ratio"     long_name="S over T diffusivity ratio"                     standard_name="salinity_over_temperature_diffusivity_ratio"                     unit="1"    /> 
     409         <field id="av_wave"      long_name="internal wave-induced vertical diffusivity"     standard_name="ocean_vertical_tracer_diffusivity_due_to_internal_waves"         unit="m2/s" /> 
     410         <field id="bflx_tmx"     long_name="internal wave-induced buoyancy flux"            standard_name="buoyancy_flux_due_to_internal_waves"                             unit="W/kg" /> 
     411         <field id="pcmap_tmx"    long_name="power consumed by wave-driven mixing"           standard_name="vertically_integrated_power_consumption_by_wave_driven_mixing"   unit="W/m2"      grid_ref="grid_W_2D" /> 
     412         <field id="emix_tmx"     long_name="power density available for mixing"             standard_name="power_available_for_mixing_from_breaking_internal_waves"         unit="W/kg" /> 
    412413 
    413414         <!-- variables available with diaar5 -->    
    … …  
    614615       <field id="ttrd_ad_e3t"       unit="degC/s * m" >  ttrd_ad  * e3t </field> 
    615616       <field id="strd_ad_e3t"       unit="1e-3/s * m" >  strd_ad  * e3t </field> 
    616        <field id="ttrd_totad_e3t"    unit="degC/s * m" >  ttrd_totad  * e3t </field> 
    617        <field id="strd_totad_e3t"    unit="1e-3/s * m" >  strd_totad  * e3t </field> 
     617       <field id="ttrd_totad_e3t"    unit="degC/s * m" >  ttrd_totad * e3t </field> 
     618       <field id="strd_totad_e3t"    unit="1e-3/s * m" >  strd_totad * e3t </field> 
    618619       <field id="ttrd_ldf_e3t"      unit="degC/s * m" >  ttrd_ldf * e3t </field> 
    619620       <field id="strd_ldf_e3t"      unit="1e-3/s * m" >  strd_ldf * e3t </field> 
    … …  
    624625 
    625626       <!-- ln_traldf_iso=T only (iso-neutral diffusion) --> 
    626        <field id="ttrd_iso_e3t"      unit="degC/s * m"  >  ttrd_iso * e3t </field> 
    627        <field id="strd_iso_e3t"      unit="1e-3/s * m"  >  strd_iso * e3t </field> 
     627       <field id="ttrd_iso_e3t"      unit="degC/s * m"  >  ttrd_iso  * e3t </field> 
     628       <field id="strd_iso_e3t"      unit="1e-3/s * m"  >  strd_iso  * e3t </field> 
    628629       <field id="ttrd_zdfp_e3t"     unit="degC/s * m"  >  ttrd_zdfp * e3t </field> 
    629630       <field id="strd_zdfp_e3t"     unit="1e-3/s * m"  >  strd_zdfp * e3t </field> 
    … …  
    669670       <field id="KE"            long_name="kinetic energy: u(n)*u(n+1)/2"                    unit="W/s^2"                        />    
    670671 
     672      <!-- variables available when explicit lateral mixing is used (ln_dynldf_NONE=F) --> 
     673      <field id="dispkexyfo"    long_name="KE-trend: lateral  mixing induced dissipation"   standard_name="ocean_kinetic_energy_dissipation_per_unit_area_due_to_xy_friction"                   unit="W/m^2" grid_ref="grid_T_2D" />    
     674      <field id="dispkevfo"     long_name="KE-trend: vertical mixing induced dissipation"   standard_name="ocean_kinetic_energy_dissipation_per_unit_area_due_to_vertical_friction"             unit="W/m^2" grid_ref="grid_T_2D" />    
     675      <!-- variables available with ln_traadv_eiv=T and ln_diaeiv=T --> 
     676      <field id="eketrd_eiv"    long_name="EKE-trend due to parameterized eddy advection"   standard_name="tendency_of_ocean_eddy_kinetic_energy_content_due_to_parameterized_eddy_advection"   unit="W/m^2" grid_ref="grid_T_2D" />    
     677 
    671678       <!-- variables available with ln_PE_trd --> 
    672679       <field id="petrd_xad"     long_name="pe-trend: i-advection"                unit="W/m^3"                        /> 
    673680       <field id="petrd_yad"     long_name="pe-trend: j-advection"                unit="W/m^3"                        /> 
    674681       <field id="petrd_zad"     long_name="pe-trend: k-advection"                unit="W/m^3"                        /> 
    675        <field id="petrd_sad"     long_name="pe-trend: surface adv. (linssh true)"      unit="W/m^3"   grid_ref="grid_T_2D" /> 
     682       <field id="petrd_sad"     long_name="pe-trend: surface adv. (linssh true)" unit="W/m^3"   grid_ref="grid_T_2D" /> 
    676683       <field id="petrd_ldf"     long_name="pe-trend: lateral  diffusion"         unit="W/m^3"                        /> 
    677684       <field id="petrd_zdf"     long_name="pe-trend: vertical diffusion"         unit="W/m^3"                        /> 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/SHARED/namelist_ice_lim3_ref

    r8586 r8637  
    1313!!             11 - Ice growth in open water           (namthd_do) 
    1414!!             12 - Ice salinity                       (namthd_sal) 
    15 !!             13 - Ice melt ponds                     (nammp) 
     15!!             13 - Ice melt ponds                     (namthd_pnd) 
    1616!!             14 - Ice initialization                 (namini) 
    1717!!             15 - Ice/snow albedos                   (namalb) 
    … …  
    2020! 
    2121!------------------------------------------------------------------------------ 
    22 &nampar     !   Generic parameters 
     22&nampar         !   Generic parameters 
    2323!------------------------------------------------------------------------------ 
    2424   jpl              =   5             !  number of ice  categories 
    … …  
    3939/ 
    4040!------------------------------------------------------------------------------ 
    41 &namitd     !   Ice discretization 
     41&namitd         !   Ice discretization 
    4242!------------------------------------------------------------------------------ 
    4343   rn_himean        =   2.0           !  expected domain-average ice thickness (m) 
    … …  
    4545/ 
    4646!------------------------------------------------------------------------------ 
    47 &namdyn     !   Ice dynamics 
     47&namdyn         !   Ice dynamics 
    4848!------------------------------------------------------------------------------ 
    4949   ln_dynFULL       = .true.          !  dyn.: full ice dynamics               (rheology + advection + ridging/rafting + correction) 
    … …  
    120120/ 
    121121!------------------------------------------------------------------------------ 
    122 &namthd     !   Ice thermodynamics 
     122&namthd         !   Ice thermodynamics 
    123123!------------------------------------------------------------------------------ 
    124124   ln_icedH         = .true.          !  activate ice thickness change from growing/melting (T) or not (F) 
    125125   ln_icedA         = .true.          !  activate lateral melting param. (T) or not (F) 
    126126   ln_icedO         = .true.          !  activate ice growth in open-water (T) or not (F) 
    127    ln_icedS         = .true.          !  activate gravity drainage and flushing (T) or not (F) 
     127   ln_icedS         = .true.          !  activate brine drainage (T) or not (F) 
    128128/ 
    129129!------------------------------------------------------------------------------ 
    130130&namthd_zdf     !   Ice heat diffusion 
    131131!------------------------------------------------------------------------------ 
    132    ln_zdf_Beer      = .true.          !  Heat diffusion follows a Beer law 
     132   ln_zdf_BL99      = .true.          !  Heat diffusion follows Bitz and Lipscomb 1999 
    133133   ln_cndi_U64      = .false.         !  sea ice thermal conductivity: k = k0 + beta.S/T            (Untersteiner, 1964) 
    134134   ln_cndi_P07      = .true.          !  sea ice thermal conductivity: k = k0 + beta1.S/T - beta2.T (Pringle et al., 2007) 
    … …  
    150150/ 
    151151!------------------------------------------------------------------------------ 
    152 &namthd_do     !   Ice growth in open water 
     152&namthd_do      !   Ice growth in open water 
    153153!------------------------------------------------------------------------------ 
    154154   rn_hinew         =   0.1           !  thickness for new ice formation in open water (m), must be larger than rn_hnewice 
    … …  
    174174/ 
    175175!------------------------------------------------------------------------------ 
    176 &nammp      !   Melt ponds 
    177 !------------------------------------------------------------------------------ 
    178    ln_pnd           = .false.         !  active melt ponds 
    179    ln_pnd_rad       = .false.         !  active melt ponds radiative coupling 
    180    ln_pnd_fw        = .false.         !  active melt ponds freshwater coupling 
    181    nn_pnd_scheme    =   0             !  type of melt pond scheme  : =0 prescribed ( Tsu=0 ), =1 empirical, =2 topographic 
    182    rn_apnd          =   0.2           !  prescribed pond fraction, at Tsu=0  : (0<rn_apnd<1, nn_pnd_scheme = 0) 
    183    rn_hpnd          =   0.05          !  prescribed pond depth, at Tsu=0     : (0<rn_apnd<1, nn_pnd_scheme = 0) 
    184 / 
    185 !------------------------------------------------------------------------------ 
    186 &namini     !   Ice initialization 
     176&namthd_pnd     !   Melt ponds 
     177!------------------------------------------------------------------------------ 
     178   ln_pnd_H12       = .false.         !  activate evolutive melt ponds (from Holland et al 2012) 
     179      ln_pnd_fwb    = .false.         !     melt ponds store freshwater or not 
     180   ln_pnd_CST       = .false.         !  activate constant melt ponds 
     181      rn_apnd       =   0.2           !     prescribed pond fraction, at Tsu=0 
     182      rn_hpnd       =   0.05          !     prescribed pond depth, at Tsu=0 
     183   ln_pnd_alb       = .false.         !  melt ponds affect albedo or not 
     184/ 
     185!------------------------------------------------------------------------------ 
     186&namini         !   Ice initialization 
    187187!------------------------------------------------------------------------------ 
    188188   ln_iceini        = .true.          !  activate ice initialization (T) or not (F) 
    … …  
    209209/ 
    210210!------------------------------------------------------------------------------ 
    211 &namalb     !   albedo parameters 
    212 !------------------------------------------------------------------------------ 
    213    nn_ice_alb       =    1            !  parameterization of ice/snow albedo 
    214       !                               !  =0 : Shine & Henderson-Sellers (JGR 1985), giving clear-sky albedo 
    215       !                               !  =1 : "home made" based on Brandt et al. (JClim 2005) and 
    216       !                               !        Grenfell & Perovich (JGR 2004), giving cloud-sky albedo 
    217       !                               !  =2 : as 1 with melt ponds 
    218       !                               ! 
    219       !                               !  recommanded values depending on the param chosen: 
    220       !                               !             nn_ice_alb   !   =0   ! = 1 or 2 !  obs range (cloud-sky) ! 
    221       rn_alb_sdry      =   0.85       !  dry snow albedo         :  0.80  !   0.85   !  0.85 -- 0.87  
    222       rn_alb_smlt      =   0.75       !  melting snow albedo     :  0.65  !   0.75   !  0.72 -- 0.82  
    223       rn_alb_idry      =   0.60       !  dry ice albedo          :  0.72  !   0.60   !  0.54 -- 0.65  
    224       rn_alb_imlt      =   0.50       !  bare puddled ice albedo :  0.53  !   0.50   !  0.49 -- 0.58  
    225       rn_alb_dpnd      =   0.27       !  ponded ice albedo       :  0.25  !   0.27   !  0.10 -- 0.30  
    226 / 
    227 !------------------------------------------------------------------------------ 
    228 &namdia     !   Diagnostics 
     211&namalb         !   albedo parameters 
     212!------------------------------------------------------------------------------ 
     213   !                                  !                          !  obs range (cloud-sky) 
     214   rn_alb_sdry      =   0.85          !  dry snow albedo         :  0.85 -- 0.87 
     215   rn_alb_smlt      =   0.75          !  melting snow albedo     :  0.72 -- 0.82 
     216   rn_alb_idry      =   0.60          !  dry ice albedo          :  0.54 -- 0.65 
     217   rn_alb_imlt      =   0.50          !  bare puddled ice albedo :  0.49 -- 0.58 
     218   rn_alb_dpnd      =   0.27          !  ponded ice albedo       :  0.10 -- 0.30  
     219/ 
     220!------------------------------------------------------------------------------ 
     221&namdia         !   Diagnostics 
    229222!------------------------------------------------------------------------------ 
    230223   ln_icediachk     = .false.         !  check online the heat, mass & salt budgets (T) or not (F) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/SHARED/namelist_ref

    r8586 r8637  
    195195                           !  =1 use observed ice-cover      , 
    196196                           !  =2 or 3 automatically for LIM3 or CICE    ("key_lim3" or "key_cice") 
     197                           !          except in AGRIF zoom where it has to be specified 
    197198   ln_ice_embd = .false.   !  =T embedded sea-ice (pressure + mass and salt exchanges) 
    198199                           !  =F levitating ice (no pressure, mass and salt exchanges) 
    … …  
    257258   rn_vfac     = 0.        !  multiplicative factor for ocean & ice velocity used to 
    258259   !                       !  calculate the wind stress (0.=absolute or 1.=relative winds) 
    259    ln_Cd_L12   = .false.   !  air-ice and ocean-ice drags = F(ice concentration) (Lupkes et al. JGR 2012) 
     260   ln_Cd_L12   = .false.   !  air-ice drags = F(ice concentration) (Lupkes et al. 2012) 
     261   ln_Cd_L15   = .false.   !  air-ice drags = F(ice concentration) (Lupkes et al. 2015) 
    260262/ 
    261263!----------------------------------------------------------------------- 
    … …  
    312314   sn_frq      = 'sas_grid_T',     120           , 'frq_m'   ,   .true.    , .true. , 'yearly'  ,    ''    ,    ''    ,    '' 
    313315 
    314    l_sasread   = .true.    !  =T Read the above fields in a file, =F initialize to 0. in sbcssm.F90 if .FALSE. 
     316   l_sasread   = .true.    !  =T Read the above fields in a file, =F initialize to 0. in sbcssm.F90 
    315317   ln_3d_uve   = .true.    !  specify whether we are supplying a 3D u,v and e3 field 
    316318   ln_read_frq = .false.   !  specify whether we must read frq or not 
    … …  
    768770! 
    769771!----------------------------------------------------------------------- 
    770 &nam_vvl    !   vertical coordinate options                             (default: z-star) 
     772&nam_vvl       !   vertical coordinate options                             (default: z-star) 
    771773!----------------------------------------------------------------------- 
    772774   ln_vvl_zstar  = .true.           !  z-star vertical coordinate 
    … …  
    775777   ln_vvl_ztilde_as_zstar = .false. !  ztilde vertical coordinate emulating zstar 
    776778   ln_vvl_zstar_at_eqtor  = .false. !  ztilde near the equator 
    777    rn_ahe3       =  0.0e0           !  thickness diffusion coefficient 
    778    rn_rst_e3t    = 30.0e0           !  ztilde to zstar restoration timescale [days] 
    779    rn_lf_cutoff  =  5.0e0           !  cutoff frequency for low-pass filter  [days] 
    780    rn_zdef_max   =  0.9e0           !  maximum fractional e3t deformation 
     779   rn_ahe3       =  0.0             !  thickness diffusion coefficient 
     780   rn_rst_e3t    = 30.0             !  ztilde to zstar restoration timescale [days] 
     781   rn_lf_cutoff  =  5.0             !  cutoff frequency for low-pass filter  [days] 
     782   rn_zdef_max   =  0.9             !  maximum fractional e3t deformation 
    781783   ln_vvl_dbg    = .true.           !  debug prints    (T/F) 
    782784/ 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/TEST_CASES/ISOMIP/EXP00/namelist_cfg

    r8568 r8637  
    6666                           !  =1 use observed ice-cover      ,  
    6767                           !  =2 ice-model used                 
    68    nn_ice_embd = 0         !  =0 levitating ice (no mass exchange, concentration/dilution effect) 
    69                            !  =1 levitating ice with mass and salt exchange but no presure effect 
    70                            !  =2 embedded sea-ice (full salt and mass exchanges and pressure) 
     68   ln_ice_embd = .false.   !  =F levitating ice with mass and salt exchange but no presure effect 
     69                           !  =T embedded sea-ice (full salt and mass exchanges and pressure) 
     70                     ! Misc. options of sbc : 
    7171   ln_traqsr   = .false.   !  Light penetration (T) or not (F) 
    7272   ln_rnf      = .false.   !  runoffs                                   (T => fill namsbc_rnf) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/TEST_CASES/SAS_BIPER/EXP00/1_namelist_cfg

    r8568 r8637  
    151151   ln_dynadv_cen2= .false. !  flux form - 2nd order centered scheme 
    152152   ln_dynadv_ubs = .true.  !  flux form - 3rd order UBS      scheme 
    153 !----------------------------------------------------------------------- 
    154153/ 
    155154!----------------------------------------------------------------------- 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/TEST_CASES/SAS_BIPER/EXP00/namelist_cfg

    r8568 r8637  
    2626&namdom        !   space and time domain (bathymetry, mesh, timestep) 
    2727!----------------------------------------------------------------------- 
    28    ln_linssh   = .false.   !  =T  linear free surface  ==>>  model level are fixed in time 
     28   ln_linssh   = .true.    !  =T  linear free surface  ==>>  model level are fixed in time 
    2929   ! 
    3030   nn_msh      =    1      !  create (>0) a mesh file or not (=0) 
    … …  
    8484&namdrg            !   top/bottom drag coefficient                      (default: NO selection) 
    8585!----------------------------------------------------------------------- 
    86    ln_NONE    = .false.    !  free-slip       : Cd = 0                   
     86   ln_NONE    = .true.     !  free-slip       : Cd = 0                   
    8787/ 
    8888!----------------------------------------------------------------------- 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/TEST_CASES/SAS_BIPER/EXP00/namelist_ice_cfg

    r7820 r8637  
    11!!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    2 !! LIM3 namelist:   
    3 !!              1 - Generic parameters                 (namicerun) 
    4 !!              2 - Diagnostics                        (namicediag) 
    5 !!              3 - Ice initialization                 (namiceini) 
    6 !!              4 - Ice discretization                 (namiceitd) 
    7 !!              5 - Ice dynamics and transport         (namicedyn) 
    8 !!              6 - Ice diffusion                      (namicehdf) 
    9 !!              7 - Ice thermodynamics                 (namicethd) 
    10 !!              8 - Ice salinity                       (namicesal) 
    11 !!              9 - Ice mechanical redistribution      (namiceitdme) 
     2!! ESIM namelist:   
     3!!              1 - Generic parameters                 (nampar) 
     4!!              2 - Ice thickness discretization       (namitd) 
     5!!              3 - Ice dynamics                       (namdyn) 
     6!!              4 - Ice ridging/rafting                (namdyn_rdgrft) 
     7!!              5 - Ice rheology                       (namdyn_rhg) 
     8!!              6 - Ice advection                      (namdyn_adv) 
     9!!              7 - Ice surface forcing                (namforcing) 
     10!!              8 - Ice thermodynamics                 (namthd) 
     11!!              9 - Ice heat diffusion                 (namthd_zdf) 
     12!!             10 - Ice lateral melting                (namthd_da) 
     13!!             11 - Ice growth in open water           (namthd_do) 
     14!!             12 - Ice salinity                       (namthd_sal) 
     15!!             13 - Ice melt ponds                     (namthd_pnd) 
     16!!             14 - Ice initialization                 (namini) 
     17!!             15 - Ice/snow albedos                   (namalb) 
     18!!             16 - Ice diagnostics                    (namdia) 
    1219!!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    1320! 
    1421!------------------------------------------------------------------------------ 
    15 &namicerun     !   Generic parameters 
     22&nampar         !   Generic parameters 
    1623!------------------------------------------------------------------------------ 
    17    jpl              =    1          !  number of ice  categories 
    18    nlay_i           =    1          !  number of ice  layers 
    19    ln_limthd        =  .false.       !  ice thermo   (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    20    ln_limdyn        =  .true.       !  ice dynamics (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    21    nn_limdyn        =   0           !     (ln_limdyn=T) switch for ice dynamics 
    22                                     !      2: total 
    23                                     !      1: advection only (no diffusion, no ridging/rafting) 
    24                                     !      0: advection only (as 1 but with prescribed velocity, bypass rheology) 
    25    rn_uice          =   0.5    !     (nn_limdyn=0) ice u-velocity 
    26    rn_vice          =   0.0    !     (nn_limdyn=0) ice v-velocity 
     24   jpl              =   1             !  number of ice  categories 
     25   nlay_i           =   1             !  number of ice  layers 
     26   ln_icedyn        = .true.          !  ice dynamics (T) or not (F) 
     27   ln_icethd        = .false.         !  ice thermo   (T) or not (F) 
    2728/ 
    2829!------------------------------------------------------------------------------ 
    29 &namicediag    !   Diagnostics 
     30&namitd         !   Ice discretization 
     31!------------------------------------------------------------------------------ 
     32   rn_himin         =   0.1           !  minimum ice thickness (m) used in remapping 
     33/ 
     34!------------------------------------------------------------------------------ 
     35&namdyn         !   Ice dynamics 
     36!------------------------------------------------------------------------------ 
     37   ln_dynFULL       = .false.         !  dyn.: full ice dynamics               (rheology + advection + ridging/rafting + correction) 
     38   ln_dynRHGADV     = .false.         !  dyn.: no ridge/raft & no corrections  (rheology + advection) 
     39   ln_dynADV        = .true.          !  dyn.: only advection w prescribed vel.(rn_uvice + advection) 
     40      rn_uice       =   0.5           !        prescribed ice u-velocity 
     41      rn_vice       =   0.            !        prescribed ice v-velocity 
     42/ 
     43!------------------------------------------------------------------------------ 
     44&namdyn_rdgrft  !   Ice ridging/rafting 
    3045!------------------------------------------------------------------------------ 
    3146/ 
    3247!------------------------------------------------------------------------------ 
    33 &namiceini     !   Ice initialization 
    34 !------------------------------------------------------------------------------ 
    35                   ! -- limistate -- ! 
    36    ln_limini      = .true.         !  activate ice initialization (T) or not (F) 
    37    ln_limini_file = .true.         !  netcdf file provided for initialization (T) or not (F) 
    38    cn_dir="./" 
    39    sn_hti = 'initice'    , -12 ,'hti'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    40    sn_hts = 'initice'    , -12 ,'hts'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    41    sn_ati = 'initice'    , -12 ,'ati'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    42    sn_tsu = 'initice'    , -12 ,'tsu'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    43    sn_tmi = 'initice'    , -12 ,'tmi'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    44    sn_smi = 'initice'    , -12 ,'smi'   ,  .false.  , .true., 'yearly'  , '' , '', '' 
    45 / 
    46 !------------------------------------------------------------------------------ 
    47 &namiceitd     !   Ice discretization 
     48&namdyn_rhg     !   Ice rheology 
    4849!------------------------------------------------------------------------------ 
    4950/ 
    5051!------------------------------------------------------------------------------ 
    51 &namicedyn     !   Ice dynamics and transport 
    52 !------------------------------------------------------------------------------ 
    53                   ! -- limtrp & limadv -- ! 
    54    nn_limadv      =    0            !  choose the advection scheme (-1=Prather ; 0=Ultimate-Macho) 
    55    nn_limadv_ord  =    5            !  choose the order of the advection scheme (if nn_limadv=0) 
    56 / 
    57 !------------------------------------------------------------------------------ 
    58 &namicehdf     !   Ice horizontal diffusion 
     52&namdyn_adv     !   Ice advection 
    5953!------------------------------------------------------------------------------ 
    6054/ 
    6155!------------------------------------------------------------------------------ 
    62 &namicethd     !   Ice thermodynamics 
     56&namforcing     !   Ice surface forcing 
    6357!------------------------------------------------------------------------------ 
    64                   ! -- limthd_dh -- ! 
    65    ln_limdH       = .true.          !  activate ice thickness change from growing/melting (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    66                   ! -- limthd_da -- ! 
    67    ln_limdA       = .true.          !  activate lateral melting param. (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    68                  ! -- limthd_lac -- ! 
    69    ln_limdO       = .true.          !  activate ice growth in open-water (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    70    rn_hnewice     = 0.02            !  thickness for new ice formation in open water (m) 
    71                   ! -- limitd_th -- ! 
    72    rn_himin       = 0.01            !  minimum ice thickness (m) used in remapping, must be smaller than rn_hnewice 
    7358/ 
    7459!------------------------------------------------------------------------------ 
    75 &namicesal     !   Ice salinity 
     60&namthd         !   Ice thermodynamics 
    7661!------------------------------------------------------------------------------ 
    77                  ! -- limthd_sal -- ! 
    78    ln_limdS       = .true.          !  activate gravity drainage and flushing (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    7962/ 
    8063!------------------------------------------------------------------------------ 
    81 &namiceitdme   !   Ice mechanical redistribution (ridging and rafting) 
     64&namthd_zdf     !   Ice heat diffusion 
    8265!------------------------------------------------------------------------------ 
    83                   ! -- limitd_me -- ! 
    84    ln_ridging     =   .true.        !  ridging activated (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    85    ln_rafting     =   .true.        !  rafting activated (T) or not (F) => DO NOT TOUCH UNLESS U KNOW WHAT U DO 
    8666/ 
     67!------------------------------------------------------------------------------ 
     68&namthd_da      !   Ice lateral melting 
     69!------------------------------------------------------------------------------ 
     70/ 
     71!------------------------------------------------------------------------------ 
     72&namthd_do      !   Ice growth in open water 
     73!------------------------------------------------------------------------------ 
     74/ 
     75!------------------------------------------------------------------------------ 
     76&namthd_sal     !   Ice salinity 
     77!------------------------------------------------------------------------------ 
     78/ 
     79!------------------------------------------------------------------------------ 
     80&namthd_pnd     !   Melt ponds 
     81!------------------------------------------------------------------------------ 
     82/ 
     83!------------------------------------------------------------------------------ 
     84&namini         !   Ice initialization 
     85!------------------------------------------------------------------------------ 
     86/ 
     87!------------------------------------------------------------------------------ 
     88&namalb         !   albedo parameters 
     89!------------------------------------------------------------------------------ 
     90/ 
     91!------------------------------------------------------------------------------ 
     92&namdia         !   Diagnostics 
     93!------------------------------------------------------------------------------ 
     94/ 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/CONFIG/cfg.txt

    r8568 r8637  
    77ORCA2_OFF_TRC OPA_SRC OFF_SRC TOP_SRC 
    88GYRE_PISCES_XIOS OPA_SRC TOP_SRC 
     9ORCA2_LIM3 OPA_SRC LIM_SRC_3 NST_SRC 
    910ORCA2_LIM3_PISCES OPA_SRC LIM_SRC_3 TOP_SRC NST_SRC 
    1011GYRE_PISCES_RK3 RK3_SRC TOP_SRC 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/ice.F90

    r8586 r8637  
    204204   REAL(wp), PUBLIC ::   rn_simin         !: minimum ice salinity [PSU] 
    205205 
    206    ! MV MP 2016 
    207    !                                     !!** melt pond namelist (nammp) 
    208    LOGICAL , PUBLIC ::   ln_pnd           !: activate ponds or not 
    209    LOGICAL , PUBLIC ::   ln_pnd_rad       !: ponds radiatively active or not 
    210    LOGICAL , PUBLIC ::   ln_pnd_fw        !: ponds active wrt meltwater or not 
    211    INTEGER , PUBLIC ::   nn_pnd_scheme    !: type of melt pond scheme:   =0 prescribed, =1 empirical, =2 topographic 
    212    REAL(wp), PUBLIC ::   rn_apnd          !: prescribed pond fraction (0<rn_apnd<1), only if nn_pnd_scheme = 0 
    213    REAL(wp), PUBLIC ::   rn_hpnd          !: prescribed pond depth    (0<rn_hpnd<1), only if nn_pnd_scheme = 0 
    214    ! END MV MP 2016 
     206   !                                     !!** namelist namthd_pnd 
     207   LOGICAL , PUBLIC ::   ln_pnd_H12       !: Melt ponds scheme from Holland et al 2012 
     208   LOGICAL , PUBLIC ::   ln_pnd_fwb       !: melt ponds store freshwater 
     209   LOGICAL , PUBLIC ::   ln_pnd_CST       !: Melt ponds scheme with constant fraction and depth 
     210   REAL(wp), PUBLIC ::   rn_apnd          !: prescribed pond fraction (0<rn_apnd<1) 
     211   REAL(wp), PUBLIC ::   rn_hpnd          !: prescribed pond depth    (0<rn_hpnd<1) 
     212   LOGICAL , PUBLIC ::   ln_pnd_alb       !: melt ponds affect albedo 
     213 
    215214   !                                     !!** ice-diagnostics namelist (namdia) ** 
    216215   LOGICAL , PUBLIC ::   ln_icediachk     !: flag for ice diag (T) or not (F) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/ice1D.F90

    r8586 r8637  
    8989 
    9090   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   sprecip_1d    !: <==> the 2D  sprecip 
    91    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   at_i_1d        !: <==> the 2D  at_i 
     91   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   at_i_1d       !: <==> the 2D  at_i 
    9292   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   ato_i_1d      !: <==> the 2D  ato_i 
    9393   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   fhtur_1d      !: <==> the 2D  fhtur 
    … …  
    103103   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_i_1d        !: <==> the 2D  a_i 
    104104   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_ib_1d       !: <==> the 2D  a_i_b 
    105    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_i_1d       !: 
    106    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_ib_1d      !: 
    107    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_s_1d       !: 
     105   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_i_1d        !: 
     106   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_ib_1d       !: 
     107   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_s_1d        !: 
    108108   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   fc_su         !: Surface Conduction flux  
    109109   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   fc_bo_i       !: Bottom  Conduction flux  
    … …  
    112112   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   dh_i_sub      !: Ice surface sublimation [m] 
    113113   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   dh_i_bott     !: Ice bottom accretion/ablation  [m] 
     114   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   dh_s_mlt      !: Snow melt [m] 
    114115   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   dh_snowice    !: Snow ice formation             [m of ice] 
    115    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   s_i_1d       !: Ice bulk salinity [ppt] 
     116   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   s_i_1d        !: Ice bulk salinity [ppt] 
    116117   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   s_i_new       !: Salinity of new ice at the bottom 
     118   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   v_i_1d        !: 
     119   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   v_s_1d        !: 
     120   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   sv_i_1d       !: 
     121   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   oa_i_1d       !: 
    117122   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_ip_1d       !: 
    118123   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   v_ip_1d       !: 
    119    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   v_i_1d        !: 
    120    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   v_s_1d        !: 
    121    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   sv_i_1d      !: 
    122    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   oa_i_1d       !: 
    123  
    124    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_s_1d   !: corresponding to the 2D var  t_s 
    125    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_i_1d   !: corresponding to the 2D var  t_i 
    126    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sz_i_1d  !: profiled ice salinity 
    127    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e_i_1d   !:    Ice  enthalpy per unit volume 
    128    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e_s_1d   !:    Snow enthalpy per unit volume 
    129  
    130    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   eh_i_old !: ice heat content (q*h, J.m-2) 
    131    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   h_i_old  !: ice thickness layer (m) 
     124   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_ip_1d       !: 
     125   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_ip_frac_1d  !: 
     126 
     127   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_s_1d      !: corresponding to the 2D var  t_s 
     128   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_i_1d      !: corresponding to the 2D var  t_i 
     129   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sz_i_1d     !: profiled ice salinity 
     130   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e_i_1d      !:    Ice  enthalpy per unit volume 
     131   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e_s_1d      !:    Snow enthalpy per unit volume 
     132 
     133   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   eh_i_old    !: ice heat content (q*h, J.m-2) 
     134   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   h_i_old     !: ice thickness layer (m) 
    132135 
    133136   ! Conduction flux diagnostics (SIMIP) 
    … …  
    197200      ! 
    198201      ii = ii + 1 
    199       ALLOCATE( t_su_1d  (jpij) , t_si_1d   (jpij) , a_i_1d  (jpij) , a_ib_1d(jpij) ,                  & 
    200          &      h_i_1d  (jpij) , h_ib_1d  (jpij) , h_s_1d (jpij) , fc_su  (jpij) , fc_bo_i(jpij) ,  &     
    201          &      dh_s_tot (jpij) , dh_i_surf (jpij) , dh_i_sub(jpij) ,                                  &     
    202          &      dh_i_bott(jpij) , dh_snowice(jpij) , s_i_1d (jpij) , s_i_new(jpij) , & 
     202      ALLOCATE( t_su_1d  (jpij) , t_si_1d   (jpij) , a_i_1d  (jpij) , a_ib_1d(jpij) ,                & 
     203         &      h_i_1d  (jpij) , h_ib_1d  (jpij) , h_s_1d (jpij) , fc_su  (jpij) , fc_bo_i(jpij) ,   &     
     204         &      dh_s_tot (jpij) , dh_i_surf (jpij) , dh_i_sub(jpij) ,                                &     
     205         &      dh_i_bott(jpij) , dh_s_mlt(jpij), dh_snowice(jpij) , s_i_1d (jpij) , s_i_new(jpij) , & 
    203206         &      a_ip_1d  (jpij) , v_ip_1d   (jpij) , v_i_1d  (jpij) , v_s_1d (jpij) , & 
     207         &      h_ip_1d  (jpij) , a_ip_frac_1d(jpij) , & 
    204208         &      sv_i_1d (jpij) , oa_i_1d   (jpij) , STAT=ierr(ii) ) 
    205209      ! 
    206210      ii = ii + 1 
    207211      ALLOCATE( t_s_1d  (jpij,nlay_s)     , t_i_1d (jpij,nlay_i)     , sz_i_1d(jpij,nlay_i) ,  &             
    208          &      e_i_1d  (jpij,nlay_i)     , e_s_1d (jpij,nlay_s)     ,                        & 
     212         &      e_i_1d  (jpij,nlay_i)     , e_s_1d (jpij,nlay_s)     ,                         & 
    209213         &      eh_i_old(jpij,0:nlay_i+1) , h_i_old(jpij,0:nlay_i+1) , STAT=ierr(ii) ) 
    210214      ! 
    … …  
    216220      ! 
    217221      ii = ii + 1 
    218       ALLOCATE( a_i_2d(jpij,jpl) , a_ib_2d(jpij,jpl) , h_i_2d(jpij,jpl) , h_ib_2d(jpij,jpl) , & 
    219          &      v_i_2d(jpij,jpl) ,v_s_2d(jpij,jpl) ,oa_i_2d(jpij,jpl) ,sv_i_2d(jpij,jpl) ,  & 
    220          &      a_ip_2d(jpij,jpl) ,v_ip_2d(jpij,jpl) ,t_su_2d(jpij,jpl) ,  & 
    221          &      STAT=ierr(ii) ) 
     222      ALLOCATE( a_i_2d (jpij,jpl) , a_ib_2d(jpij,jpl) , h_i_2d (jpij,jpl) , h_ib_2d(jpij,jpl) ,   & 
     223         &      v_i_2d (jpij,jpl) , v_s_2d (jpij,jpl) , oa_i_2d(jpij,jpl) , sv_i_2d(jpij,jpl) ,   & 
     224         &      a_ip_2d(jpij,jpl) , v_ip_2d(jpij,jpl) , t_su_2d(jpij,jpl) ,                   STAT=ierr(ii) ) 
    222225 
    223226      ice1D_alloc = MAXVAL( ierr(:) ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icealb.F90

    r8586 r8637  
    3535   REAL(wp) , PARAMETER ::   pp1_c2 = 1. / ppc2 
    3636   ! 
    37    ! ** albedo namelist (namalb) 
    38    INTEGER  ::   nn_ice_alb       ! type of albedo scheme: 0: Shine & Henderson-Sellers (JGR 1985) 
    39    !                                      !                         1: "home made" based on Brandt et al. (JClim 2005) 
    40    !                                      !                                             and Grenfell & Perovich (JGR 2004) 
    41    !                                      !                         2: Same as 1 but with melt ponds 
     37   !                             !!* albedo namelist (namalb) 
    4238   REAL(wp) ::   rn_alb_sdry      ! dry snow albedo 
    4339   REAL(wp) ::   rn_alb_smlt      ! melting snow albedo 
    … …  
    5349CONTAINS 
    5450 
    55    SUBROUTINE ice_alb( pt_ice, ph_ice, ph_snw, pafrac_pnd, ph_pnd, ld_pnd, pa_ice_cs, pa_ice_os ) 
     51   SUBROUTINE ice_alb( pt_su, ph_ice, ph_snw, ld_pnd_alb, pafrac_pnd, ph_pnd, palb_cs, palb_os ) 
    5652      !!---------------------------------------------------------------------- 
    5753      !!               ***  ROUTINE ice_alb  *** 
    … …  
    5955      !! ** Purpose :   Computation of the albedo of the snow/ice system  
    6056      !!        
    61       !! ** Method  :   Two schemes are available (from namelist parameter nn_ice_alb) 
    62       !!                  0: the scheme is that of Shine & Henderson-Sellers (JGR 1985) for clear-skies 
    63       !!                  1: the scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005) 
    64       !!                                                                           and Grenfell & Perovich (JGR 2004) 
    65       !!                  2: fractional surface-based formulation of scheme 1 (NEW) 
    66       !!                Description of scheme 1: 
     57      !! ** Method  :   The scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005) 
     58      !!                                                                      and Grenfell & Perovich (JGR 2004) 
    6759      !!                  1) Albedo dependency on ice thickness follows the findings from Brandt et al (2005) 
    6860      !!                     which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999 
    … …  
    7668      !!                  4) The needed 4 parameters are: dry and melting snow, freezing ice and bare puddled ice 
    7769      !! 
    78       !! ** Note    :   The parameterization from Shine & Henderson-Sellers presents several misconstructions: 
     70      !!                     compilation of values from literature (reference overcast sky values) 
     71      !!                        rn_alb_sdry = 0.85      ! dry snow 
     72      !!                        rn_alb_smlt = 0.75      ! melting snow 
     73      !!                        rn_alb_idry = 0.60      ! bare frozen ice 
     74      !!                        rn_alb_imlt = 0.50      ! bare puddled ice albedo 
     75      !!                        rn_alb_dpnd = 0.36      ! ponded-ice overcast albedo (Lecomte et al, 2015) 
     76      !!                                                ! early melt pnds 0.27, late melt ponds 0.14 Grenfell & Perovich 
     77      !!                     Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved 
     78      !!                        rn_alb_sdry = 0.85      ! dry snow 
     79      !!                        rn_alb_smlt = 0.72      ! melting snow 
     80      !!                        rn_alb_idry = 0.65      ! bare frozen ice 
     81      !!                     Brandt et al 2005 (East Antarctica) 
     82      !!                        rn_alb_sdry = 0.87      ! dry snow 
     83      !!                        rn_alb_smlt = 0.82      ! melting snow 
     84      !!                        rn_alb_idry = 0.54      ! bare frozen ice 
     85      !! 
     86      !! ** Note    :   The old parameterization from Shine & Henderson-Sellers (not here anymore) presented several misconstructions 
    7987      !!                  1) ice albedo when ice thick. tends to 0 is different than ocean albedo 
    8088      !!                  2) for small ice thick. covered with some snow (<3cm?), albedo is larger  
    … …  
    8795      !!                Grenfell & Perovich 2004, JGR, vol 109  
    8896      !!---------------------------------------------------------------------- 
    89       REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pt_ice       !  ice surface temperature (Kelvin) 
     97      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pt_su        !  ice surface temperature (Kelvin) 
    9098      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_ice       !  sea-ice thickness 
    9199      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_snw       !  snow depth 
     100      LOGICAL , INTENT(in   )                   ::   ld_pnd_alb   !  effect of melt ponds on albedo 
    92101      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pafrac_pnd   !  melt pond relative fraction (per unit ice area) 
    93102      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_pnd       !  melt pond depth 
    94       LOGICAL , INTENT(in   )                   ::   ld_pnd       !  melt ponds radiatively active or not 
    95       REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_cs    !  albedo of ice under clear    sky 
    96       REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_os    !  albedo of ice under overcast sky 
     103      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   palb_cs      !  albedo of ice under clear    sky 
     104      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   palb_os      !  albedo of ice under overcast sky 
    97105      ! 
    98106      INTEGER  ::   ji, jj, jl                ! dummy loop indices 
    99       REAL(wp) ::   zswitch, z1_c1, z1_c2     ! local scalar 
    100       REAL(wp) ::   z1_href_pnd               !   -      -                  
    101       REAL(wp) ::   zalb_sm, zalb_sf, zalb_st ! albedo of snow melting, freezing, total 
    102       ! 
    103       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb, zalb_it   ! intermediate variable & albedo of ice (snow free) 
    104 !! MV MP 
    105       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_pnd        ! ponded sea ice albedo 
    106       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_ice        ! bare sea ice albedo 
    107       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb_snw        ! snow-covered sea ice albedo 
    108       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_snw      ! relative snow fraction 
    109       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_ice      ! relative ice fraction 
    110       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zafrac_pnd      ! relative ice fraction (effective) 
     107      REAL(wp) ::   z1_c1, z1_c2,z1_c3, z1_c4 ! local scalar 
     108      REAL(wp) ::   z1_href_pnd               ! inverse of the characteristic length scale (Lecomte et al. 2015) 
     109      REAL(wp) ::   zalb_pnd, zafrac_pnd      ! ponded sea ice albedo & relative pound fraction 
     110      REAL(wp) ::   zalb_ice, zafrac_ice      ! bare sea ice albedo & relative ice fraction 
     111      REAL(wp) ::   zalb_snw, zafrac_snw      ! snow-covered sea ice albedo & relative snow fraction 
    111112      !!--------------------------------------------------------------------- 
    112113      ! 
    113114      IF( nn_timing == 1 )   CALL timing_start('icealb') 
    114115      ! 
    115       !----------------------------------------------------- 
    116       !  Snow-free albedo (no ice thickness dependence yet) 
    117       !----------------------------------------------------- 
    118       ! 
    119       ! Part common to nn_ice_alb = 0, 1, 2 
    120       ! 
    121       IF ( .NOT. ld_pnd ) THEN   !--- No melt ponds OR radiatively inactive melt ponds 
    122          ! Bare ice albedo is prescribed, with implicit assumption on pond fraction and depth 
    123          WHERE     ( ph_snw == 0._wp .AND. pt_ice >= rt0 )   ;   zalb(:,:,:) = rn_alb_imlt 
    124          ELSEWHERE                                           ;   zalb(:,:,:) = rn_alb_idry 
    125          END WHERE 
    126       ENDIF 
    127  
    128       SELECT CASE ( nn_ice_alb )       ! select a parameterization 
    129       ! 
    130       !              !------------------------------------------ 
    131       CASE( 0 )      !  Shine and Henderson-Sellers (1985)     !   (based on clear sky values) 
    132          !           !------------------------------------------ 
    133          ! 
    134          !                       ! Thickness-dependent bare ice albedo 
    135          WHERE    ( 1.5  < ph_ice                     )   ;   zalb_it = zalb 
    136          ELSEWHERE( 1.0  < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_it = 0.472  + 2.0 * ( zalb - 0.472 ) * ( ph_ice - 1.0 ) 
    137          ELSEWHERE( 0.05 < ph_ice .AND. ph_ice <= 1.0 )   ;   zalb_it = 0.2467 + 0.7049 * ph_ice              & 
    138             &                                                                  - 0.8608 * ph_ice * ph_ice     & 
    139             &                                                                  + 0.3812 * ph_ice * ph_ice * ph_ice 
    140          ELSEWHERE                                        ;   zalb_it = 0.1    + 3.6    * ph_ice 
    141          END WHERE 
    142          ! 
    143          IF ( ld_pnd ) THEN      ! Depth-dependent ponded ice albedo 
    144             z1_href_pnd = 1. / 0.05    ! inverse of the characteristic length scale (Lecomte et al. 2015) 
    145             zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_it ) * EXP( - ph_pnd * z1_href_pnd )  
    146             !                          ! Snow-free ice albedo is a function of pond fraction 
    147             WHERE ( ph_snw == 0._wp .AND. pt_ice >= rt0 )    
    148                zalb_it = zalb_it * ( 1. - pafrac_pnd  ) + zalb_pnd * pafrac_pnd  
    149             END WHERE 
    150          ENDIF  
    151          ! 
    152          DO jl = 1, jpl 
    153             DO jj = 1, jpj 
    154                DO ji = 1, jpi 
    155                   !                    ! Freezing snow 
    156                   ! no effect of underlying ice layer IF snow thickness > c1. Albedo does not depend on snow thick if > ppc2 
    157                   zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( ph_snw(ji,jj,jl) - ppc1 ) ) ) 
    158                   zalb_sf = ( 1._wp - zswitch ) * (  zalb_it(ji,jj,jl)  & 
    159                      &                                   + ph_snw(ji,jj,jl) * ( rn_alb_sdry - zalb_it(ji,jj,jl) ) * pp1_c1  )   & 
    160                      &    +           zswitch   * rn_alb_sdry   
    161                      ! 
    162                   !                    ! Melting snow 
    163                   ! no effect of underlying ice layer. Albedo does not depend on snow thick IF > ppc2 
    164                   zswitch = MAX( 0._wp , SIGN( 1._wp , ph_snw(ji,jj,jl) - ppc2 ) ) 
    165                   zalb_sm = ( 1._wp - zswitch ) * ( rn_alb_imlt + ph_snw(ji,jj,jl) * ( rn_alb_smlt - rn_alb_imlt ) * pp1_c2 )   & 
    166                      &      +         zswitch   *   rn_alb_smlt  
    167                      ! 
    168                   !                    ! Snow albedo 
    169                   zswitch =  MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )    
    170                   zalb_st =  zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf 
    171                   ! 
    172                   !                    ! Surface albedo 
    173                   zswitch             = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) ) 
    174                   pa_ice_cs(ji,jj,jl) = zswitch * zalb_st + ( 1._wp - zswitch ) * zalb_it(ji,jj,jl) 
    175                   ! 
    176                END DO 
     116      z1_href_pnd = 0.05 
     117      z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) )  
     118      z1_c2 = 1. / 0.05 
     119      z1_c3 = 1. / 0.02 
     120      z1_c4 = 1. / 0.03 
     121      ! 
     122      DO jl = 1, jpl 
     123         DO jj = 1, jpj 
     124            DO ji = 1, jpi 
     125               !                       !--- Specific snow, ice and pond fractions (for now, we prevent melt ponds and snow at the same time) 
     126               IF( ph_snw(ji,jj,jl) == 0._wp ) THEN 
     127                  zafrac_snw = 0._wp 
     128                  IF( ld_pnd_alb ) THEN 
     129                     zafrac_pnd = pafrac_pnd(ji,jj,jl) 
     130                  ELSE 
     131                     zafrac_pnd = 0._wp 
     132                  ENDIF 
     133                  zafrac_ice = 1._wp - zafrac_pnd 
     134               ELSE 
     135                  zafrac_snw = 1._wp      ! Snow fully "shades" melt ponds and ice 
     136                  zafrac_pnd = 0._wp 
     137                  zafrac_ice = 0._wp 
     138               ENDIF 
     139               ! 
     140               !                       !--- Bare ice albedo (for hi > 150cm) 
     141               IF( ld_pnd_alb ) THEN 
     142                  zalb_ice = rn_alb_idry 
     143               ELSE 
     144                  IF( ph_snw(ji,jj,jl) == 0._wp .AND. pt_su(ji,jj,jl) >= rt0 ) THEN  ;   zalb_ice = rn_alb_imlt 
     145                  ELSE                                                               ;   zalb_ice = rn_alb_idry   ;   ENDIF 
     146               ENDIF 
     147               !                       !--- Bare ice albedo (for hi < 150cm) 
     148               IF( 0.05 < ph_ice(ji,jj,jl) .AND. ph_ice(ji,jj,jl) <= 1.5 ) THEN      ! 5cm < hi < 150cm 
     149                  zalb_ice = zalb_ice    + ( 0.18 - zalb_ice   ) * z1_c1 * ( LOG(1.5) - LOG(ph_ice(ji,jj,jl)) ) 
     150               ELSEIF( ph_ice(ji,jj,jl) <= 0.05 ) THEN                               ! 0cm < hi < 5cm 
     151                  zalb_ice = rn_alb_oce  + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice(ji,jj,jl) 
     152               ENDIF 
     153               ! 
     154               !                       !--- Snow-covered ice albedo (freezing, melting cases) 
     155               IF( pt_su(ji,jj,jl) < rt0_snow ) THEN 
     156                  zalb_snw = rn_alb_sdry - ( rn_alb_sdry - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c3 ) 
     157               ELSE 
     158                  zalb_snw = rn_alb_smlt - ( rn_alb_smlt - zalb_ice ) * EXP( - ph_snw(ji,jj,jl) * z1_c4 ) 
     159               ENDIF 
     160               !                       !--- Ponded ice albedo 
     161               IF( ld_pnd_alb ) THEN 
     162                  zalb_pnd = rn_alb_dpnd - ( rn_alb_dpnd - zalb_ice ) * EXP( - ph_pnd(ji,jj,jl) * z1_href_pnd )  
     163               ELSE 
     164                  zalb_pnd = rn_alb_dpnd 
     165               ENDIF 
     166               !                       !--- Surface albedo is weighted mean of snow, ponds and bare ice contributions 
     167               palb_os(ji,jj,jl) = zafrac_snw * zalb_snw + zafrac_pnd * zalb_pnd + zafrac_ice * zalb_ice 
     168               ! 
    177169            END DO 
    178170         END DO 
    179          ! 
    180          pa_ice_os(:,:,:) = pa_ice_cs(:,:,:) + ppcloud       ! Oberhuber correction for overcast sky 
    181          ! 
    182          ! 
    183          !           !------------------------------------------ 
    184       CASE( 1 )      !  New parameterization (2016)            !    (based on overcast sky values) 
    185          !           !------------------------------------------ 
    186          ! 
    187          !                    compilation of values from literature (reference overcast sky values) 
    188          !                          rn_alb_sdry = 0.85      ! dry snow 
    189          !                          rn_alb_smlt = 0.75      ! melting snow 
    190          !                          rn_alb_idry = 0.60      ! bare frozen ice 
    191          !                          rn_alb_dpnd = 0.36      ! ponded-ice overcast albedo (Lecomte et al, 2015) 
    192          !                                                  ! early melt pnds 0.27, late melt ponds 0.14 Grenfell & Perovich 
    193          !                    Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved 
    194          !                          rn_alb_sdry = 0.85      ! dry snow 
    195          !                          rn_alb_smlt = 0.72      ! melting snow 
    196          !                          rn_alb_idry = 0.65      ! bare frozen ice 
    197          !                    Brandt et al 2005 (East Antarctica) 
    198          !                          rn_alb_sdry = 0.87      ! dry snow 
    199          !                          rn_alb_smlt = 0.82      ! melting snow 
    200          !                          rn_alb_idry = 0.54      ! bare frozen ice 
    201          !  
    202          !              !--- Computation of snow-free ice albedo  
    203          z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) )  
    204          z1_c2 = 1. / 0.05 
    205  
    206          !              !--- Accounting for the ice-thickness dependency 
    207          WHERE     ( 1.5  < ph_ice                     )   ;   zalb_it = zalb 
    208          ELSE WHERE( 0.05 < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_it = zalb + ( 0.18 - zalb ) * z1_c1 * ( LOG(1.5) - LOG(ph_ice) ) 
    209          ELSE WHERE                                        ;   zalb_it = rn_alb_oce + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice 
    210          END WHERE 
    211          ! 
    212          IF ( ld_pnd ) THEN      ! Depth-dependent ponded ice albedo 
    213             z1_href_pnd = 0.05        ! inverse of the characteristic length scale (Lecomte et al. 2015) 
    214             zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_it ) * EXP( - ph_pnd * z1_href_pnd )  
    215             ! 
    216             !                    ! Snow-free ice albedo is weighted mean of ponded ice and bare ice contributions 
    217             WHERE ( ph_snw == 0._wp .AND. pt_ice >= rt0 ) 
    218                zalb_it = zalb_it * ( 1. - pafrac_pnd  ) + zalb_pnd * pafrac_pnd 
    219             END WHERE 
    220          ENDIF  
    221          ! 
    222          !              !--- Overcast sky surface albedo (accounting for snow, ice melt ponds) 
    223          z1_c1 = 1. / 0.02 
    224          z1_c2 = 1. / 0.03 
    225          DO jl = 1, jpl 
    226             DO jj = 1, jpj 
    227                DO ji = 1, jpi 
    228                   ! Snow depth dependence of snow albedo 
    229                   zalb_sf = rn_alb_sdry - ( rn_alb_sdry - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c1 ); 
    230                   zalb_sm = rn_alb_smlt - ( rn_alb_smlt - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c2 ); 
    231  
    232                   ! Snow albedo (MV I guess we could use rt0 instead of rt0_snow) 
    233                   zswitch = MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )    
    234                   zalb_st = zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf 
    235  
    236                   ! Surface albedo    
    237                   zswitch             = MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) ) 
    238                   pa_ice_os(ji,jj,jl) = ( 1._wp - zswitch ) * zalb_st + zswitch *  zalb_it(ji,jj,jl) 
    239  
    240               END DO 
    241             END DO 
    242          END DO 
    243          ! 
    244          !              !--- Clear sky surface albedo 
    245          pa_ice_cs = pa_ice_os - ( - 0.1010 * pa_ice_os * pa_ice_os + 0.1933 * pa_ice_os - 0.0148 );  
    246          ! 
    247          ! 
    248          !           !------------------------------------------ 
    249       CASE( 2 )      !  Fractional surface-based parameterization (2017) 
    250          !           !------------------------------------------ 
    251          !              MV: I propose this formulation that is more elegant, and more easy to expand towards 
    252          !              varying pond and snow fraction. 
    253          !              Formulation 1 has issues to handle ponds and snow together that can't easily be fixed.  
    254          !              This one handles it better, I believe. 
    255          ! 
    256          !                 !==  Snow, bare ice and ponded ice fractions  ==! 
    257          ! 
    258          !                       ! Specific fractions (zafrac) refer to relative area covered by snow within each ice category 
    259          ! 
    260          !                       !--- Effective pond fraction (for now, we prevent melt ponds and snow at the same time) 
    261          zafrac_pnd = 0._wp 
    262          IF ( ld_pnd ) THEN 
    263             WHERE( ph_snw == 0._wp )   zafrac_pnd = pafrac_pnd   ! Snow fully "shades" melt ponds 
    264          ENDIF          
    265          ! 
    266          !                       !--- Specific snow fraction (for now, prescribed) 
    267          WHERE( ph_snw > 0._wp )   ;   zafrac_snw = 1. 
    268          ELSEWHERE                 ;   zafrac_snw = 0. 
    269          END WHERE 
    270          ! 
    271          !                       !--- Specific ice fraction  
    272          zafrac_ice = 1. - zafrac_snw - zafrac_pnd 
    273          ! 
    274          !                 !==  Snow-covered, pond-covered, and bare ice albedos  ==! 
    275          ! 
    276          !                       !--- Bare ice albedo 
    277          z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) )  
    278          z1_c2 = 1. / 0.05 
    279          WHERE    ( 1.5  < ph_ice                     )   ;   zalb_ice = zalb 
    280          ELSEWHERE( 0.05 < ph_ice .AND. ph_ice <= 1.5 )   ;   zalb_ice = zalb       + ( 0.18 - zalb       ) * z1_c1 *  & 
    281             &                                                                       ( LOG(1.5) - LOG(ph_ice) ) 
    282          ELSEWHERE                                        ;   zalb_ice = rn_alb_oce + ( 0.18 - rn_alb_oce ) * z1_c2 * ph_ice 
    283          END WHERE 
    284          ! 
    285          z1_c1 = 1. / 0.02       !--- Snow-covered ice albedo (freezing, melting cases) 
    286          z1_c2 = 1. / 0.03 
    287          ! 
    288          WHERE( pt_ice < rt0_snow )   ;   zalb_snw = rn_alb_sdry - ( rn_alb_sdry - zalb_ice ) * EXP( - ph_snw * z1_c1 ) 
    289          ELSEWHERE                    ;   zalb_snw = rn_alb_smlt - ( rn_alb_smlt - zalb_ice ) * EXP( - ph_snw * z1_c2 ) 
    290          END WHERE 
    291          ! 
    292          IF ( ld_pnd ) THEN      !--- Depth-dependent ponded ice albedo 
    293             z1_href_pnd = 0.05        ! inverse of the characteristic length scale (Lecomte et al. 2015) 
    294             zalb_pnd  = rn_alb_dpnd - ( rn_alb_dpnd - zalb_ice ) * EXP( - ph_pnd * z1_href_pnd )  
    295          ELSE 
    296             zalb_pnd  = rn_alb_dpnd 
    297          ENDIF 
    298          !                       !--- Surface albedo is weighted mean of snow, ponds and bare ice contributions 
    299          pa_ice_os = zafrac_snw * zalb_snw  +  zafrac_pnd * zalb_pnd  +  zafrac_ice * zalb_ice 
    300          pa_ice_cs = pa_ice_os - ( - 0.1010 * pa_ice_os * pa_ice_os + 0.1933 * pa_ice_os - 0.0148 ) 
    301  
    302       END SELECT 
     171      END DO 
     172      ! 
     173      palb_cs(:,:,:) = palb_os(:,:,:) - ( - 0.1010 * palb_os(:,:,:) * palb_os(:,:,:) + 0.1933 * palb_os(:,:,:) - 0.0148 ) 
    303174      ! 
    304175      IF( nn_timing == 1 )   CALL timing_stop('icealb') 
    … …  
    315186      !! ** Method  :   Read the namelist namalb 
    316187      !!---------------------------------------------------------------------- 
    317       INTEGER  ::   ios                 ! Local integer output status for namelist read 
    318       !! 
    319       NAMELIST/namalb/ nn_ice_alb, rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd 
     188      INTEGER ::   ios   ! Local integer output status for namelist read 
     189      !! 
     190      NAMELIST/namalb/ rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd 
    320191      !!---------------------------------------------------------------------- 
    321192      ! 
    … …  
    334205         WRITE(numout,*) '~~~~~~~~~~~~' 
    335206         WRITE(numout,*) '   Namelist namalb:' 
    336          WRITE(numout,*) '      choose the albedo parameterization   nn_ice_alb  = ', nn_ice_alb 
    337207         WRITE(numout,*) '      albedo of dry snow                   rn_alb_sdry = ', rn_alb_sdry 
    338208         WRITE(numout,*) '      albedo of melting snow               rn_alb_smlt = ', rn_alb_smlt 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90

    r8586 r8637  
    7676      IF( icount == 0 ) THEN 
    7777         !                          ! water flux 
    78          pdiag_fv = glob_sum( -( wfx_bog(:,:) + wfx_bom(:,:)     + wfx_sum(:,:)     + wfx_sni(:,:)     +                     & 
    79             &                    wfx_opw(:,:) + wfx_res(:,:)     + wfx_dyn(:,:)     + wfx_lam(:,:)     + wfx_ice_sub(:,:) +  & 
    80             &                    wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) + wfx_spr(:,:)    & 
     78         pdiag_fv = glob_sum( -( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) +                  & 
     79            &                    wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_lam(:,:) + wfx_pnd(:,:)  +  & 
     80            &                    wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) +  & 
     81            &                    wfx_ice_sub(:,:) + wfx_spr(:,:)  & 
    8182            &                  ) * e1e2t(:,:) ) * zconv 
    8283         ! 
    … …  
    9697 
    9798         pdiag_t = glob_sum( (  SUM( SUM( e_i(:,:,1:nlay_i,:), dim=4 ), dim=3 )     & 
    98             &                 + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 )   ) * e1e2t ) * zconv 
     99            &                 + SUM( SUM( e_s(:,:,1:nlay_s,:), dim=4 ), dim=3 ) ) * e1e2t ) * zconv 
    99100 
    100101      ELSEIF( icount == 1 ) THEN 
    101102 
    102103         ! water flux 
    103          zfv  = glob_sum( -( wfx_bog(:,:) + wfx_bom(:,:)     + wfx_sum(:,:)     + wfx_sni(:,:)     +                     & 
    104             &                wfx_opw(:,:) + wfx_res(:,:)     + wfx_dyn(:,:)     + wfx_lam(:,:)     + wfx_ice_sub(:,:) +  & 
    105             &                wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) + wfx_spr(:,:)        & 
     104         zfv  = glob_sum( -( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) +                  & 
     105            &                wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_lam(:,:) + wfx_pnd(:,:)  +  & 
     106            &                wfx_snw_sni(:,:) + wfx_snw_sum(:,:) + wfx_snw_dyn(:,:) + wfx_snw_sub(:,:) +  & 
     107            &                wfx_ice_sub(:,:) + wfx_spr(:,:)  & 
    106108            &              ) * e1e2t(:,:) ) * zconv - pdiag_fv 
    107109 
    … …  
    117119  
    118120         ! outputs 
    119          zv = ( ( glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t  ) * zconv  & 
     121         zv = ( ( glob_sum( SUM( v_i * rhoic + v_s * rhosn, dim=3 ) * e1e2t ) * zconv  & 
    120122            &     - pdiag_v ) * r1_rdtice - zfv ) * rday 
    121123 
    122          zs = ( ( glob_sum( SUM( sv_i * rhoic            , dim=3 ) * e1e2t ) * zconv  & 
     124         zs = ( ( glob_sum( SUM( sv_i * rhoic             , dim=3 ) * e1e2t ) * zconv  & 
    123125            &     - pdiag_s ) * r1_rdtice + zfs ) * rday 
    124126 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icedyn.F90

    r8586 r8637  
    172172            END DO 
    173173         END DO 
    174       END DO 
    175              
    176       IF ( nn_pnd_scheme > 0 ) THEN               !-- correct pond fraction to avoid a_ip > a_i 
    177          WHERE( a_ip(:,:,:) > a_i(:,:,:) )   a_ip(:,:,:) = a_i(:,:,:) 
    178       ENDIF 
     174      END DO  
     175      !                                           !-- correct pond fraction to avoid a_ip > a_i 
     176      WHERE( a_ip(:,:,:) > a_i(:,:,:) )   a_ip(:,:,:) = a_i(:,:,:) 
    179177      ! 
    180178   END SUBROUTINE Hbig 
    … …  
    203201      ! 
    204202   END SUBROUTINE Hpiling 
     203 
    205204 
    206205   SUBROUTINE ice_dyn_init 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icedyn_adv_pra.F90

    r8586 r8637  
    120120            z0ei(:,:,jk,jl) = pe_i(:,:,jk,jl) * e1e2t(:,:) ! Ice  heat content 
    121121         END DO 
    122          IF ( nn_pnd_scheme > 0 ) THEN 
     122         IF ( ln_pnd_H12 ) THEN 
    123123            z0ap(:,:,jl)  = pa_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond fraction 
    124124            z0vp(:,:,jl)  = pv_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond volume 
    … …  
    167167                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    168168               END DO 
    169                IF ( nn_pnd_scheme > 0 ) THEN 
     169               IF ( ln_pnd_H12 ) THEN 
    170170                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &    !--- melt pond fraction -- 
    171171                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  ) 
    … …  
    220220                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    221221               END DO 
    222                IF ( nn_pnd_scheme > 0 ) THEN 
     222               IF ( ln_pnd_H12 ) THEN 
    223223                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &   !--- melt pond fraction --- 
    224224                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  ) 
    … …  
    249249         END DO 
    250250         ! MV MP 2016 
    251          IF ( nn_pnd_scheme > 0 ) THEN 
     251         IF ( ln_pnd_H12 ) THEN 
    252252            pa_ip  (:,:,jl)   = z0ap (:,:,jl) * r1_e1e2t(:,:) 
    253253            pv_ip  (:,:,jl)   = z0vp (:,:,jl) * r1_e1e2t(:,:) 
    … …  
    755755                  sxyage(:,:,jl)= z2d(:,:) 
    756756               END DO 
    757                IF ( nn_pnd_scheme > 0 ) THEN 
     757               IF ( ln_pnd_H12 ) THEN 
    758758                  DO jl = 1, jpl  
    759759                     WRITE(zchar,'(I2.2)') jl 
    … …  
    833833               syyc0 (:,:,:) = 0._wp   ;   syye (:,:,:,:) = 0._wp   ;   syysal (:,:,:) = 0._wp   ;   syyage (:,:,:) = 0._wp 
    834834               sxyc0 (:,:,:) = 0._wp   ;   sxye (:,:,:,:) = 0._wp   ;   sxysal (:,:,:) = 0._wp   ;   sxyage (:,:,:) = 0._wp 
    835                IF ( nn_pnd_scheme > 0 ) THEN 
     835               IF ( ln_pnd_H12 ) THEN 
    836836                  sxap  (:,:,:) = 0._wp    ; sxvp  (:,:,:) = 0._wp  
    837837                  syap  (:,:,:) = 0._wp    ; syvp  (:,:,:) = 0._wp  
    … …  
    854854            syyc0 (:,:,:) = 0._wp   ;   syye (:,:,:,:) = 0._wp   ;   syysal (:,:,:) = 0._wp   ;   syyage (:,:,:) = 0._wp 
    855855            sxyc0 (:,:,:) = 0._wp   ;   sxye (:,:,:,:) = 0._wp   ;   sxysal (:,:,:) = 0._wp   ;   sxyage (:,:,:) = 0._wp 
    856             IF ( nn_pnd_scheme > 0 ) THEN 
     856            IF ( ln_pnd_H12 ) THEN 
    857857               sxap  (:,:,:) = 0._wp    ; sxvp  (:,:,:) = 0._wp  
    858858               syap  (:,:,:) = 0._wp    ; syvp  (:,:,:) = 0._wp  
    … …  
    989989            END DO 
    990990         END DO 
    991          IF ( nn_pnd_scheme > 0 ) THEN 
     991         IF ( ln_pnd_H12 ) THEN 
    992992            DO jl = 1, jpl  
    993993               WRITE(zchar,'(I2.2)') jl 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icedyn_adv_umx.F90

    r8586 r8637  
    124124            CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pv_s(:,:,jl) )         ! Snow volume 
    125125            CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pe_s(:,:,1,jl) )       ! Snow heat content 
    126             IF ( nn_pnd_scheme > 0 ) THEN 
     126            IF ( ln_pnd_H12 ) THEN 
    127127               CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pa_ip(:,:,jl) )     ! Melt pond fraction 
    128128               CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pv_ip(:,:,jl) )     ! Melt pond volume 
    … …  
    191191      SELECT CASE( k_order ) 
    192192      CASE ( 20 )                          ! centered second order 
    193          DO jj = 2, jpjm1 
    194             DO ji = fs_2, fs_jpim1   ! vector opt. 
     193         DO jj = 1, jpjm1 
     194            DO ji = 1, fs_jpim1   ! vector opt. 
    195195               zfu_ho(ji,jj) = 0.5 * puc(ji,jj) * ( ptc(ji,jj) + ptc(ji+1,jj) ) 
    196196               zfv_ho(ji,jj) = 0.5 * pvc(ji,jj) * ( ptc(ji,jj) + ptc(ji,jj+1) ) 
    … …  
    201201         CALL macho( k_order, kt, pdt, ptc, puc, pvc, pubox, pvbox, zt_u, zt_v ) 
    202202         ! 
    203          DO jj = 2, jpjm1 
    204             DO ji = fs_2, fs_jpim1   ! vector opt. 
     203         DO jj = 1, jpjm1 
     204            DO ji = 1, fs_jpim1   ! vector opt. 
    205205               zfu_ho(ji,jj) = puc(ji,jj) * zt_u(ji,jj) 
    206206               zfv_ho(ji,jj) = pvc(ji,jj) * zt_v(ji,jj) 
    … …  
    212212      ! antidiffusive flux : high order minus low order 
    213213      ! -------------------------------------------------- 
    214       DO jj = 2, jpjm1 
    215          DO ji = fs_2, fs_jpim1   ! vector opt. 
     214      DO jj = 1, jpjm1 
     215         DO ji = 1, fs_jpim1   ! vector opt. 
    216216            zfu_ho(ji,jj) = zfu_ho(ji,jj) - zfu_ups(ji,jj) 
    217217            zfv_ho(ji,jj) = zfv_ho(ji,jj) - zfv_ups(ji,jj) 
    218218         END DO 
    219219      END DO 
    220       CALL lbc_lnk_multi( zfu_ho, 'U', -1., zfv_ho, 'V', -1. )         ! Lateral bondary conditions 
    221220       
    222221      ! monotonicity algorithm 
    … …  
    360359      CASE( 1 )                                                   !==  1st order central TIM  ==! (Eq. 21) 
    361360         !         
    362          DO jj = 1, jpj 
     361         DO jj = 2, jpjm1 
    363362            DO ji = 1, fs_jpim1   ! vector opt. 
    364363               pt_u(ji,jj) = 0.5_wp * umask(ji,jj,1) * (                               pt(ji+1,jj) + pt(ji,jj)   & 
    … …  
    369368      CASE( 2 )                                                   !==  2nd order central TIM  ==! (Eq. 23) 
    370369         ! 
    371          DO jj = 1, jpj 
     370         DO jj = 2, jpjm1 
    372371            DO ji = 1, fs_jpim1   ! vector opt. 
    373372               zcu  = puc(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    … …  
    376375            END DO 
    377376         END DO 
    378          CALL lbc_lnk( pt_u(:,:) , 'U',  1. ) 
    379377         !   
    380378      CASE( 3 )                                                   !==  3rd order central TIM  ==! (Eq. 24) 
    381379         ! 
    382          DO jj = 1, jpj 
     380         DO jj = 2, jpjm1 
    383381            DO ji = 1, fs_jpim1   ! vector opt. 
    384382               zcu  = puc(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    … …  
    394392      CASE( 4 )                                                   !==  4th order central TIM  ==! (Eq. 27) 
    395393         ! 
    396          DO jj = 1, jpj 
     394         DO jj = 2, jpjm1 
    397395            DO ji = 1, fs_jpim1   ! vector opt. 
    398396               zcu  = puc(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    … …  
    408406      CASE( 5 )                                                   !==  5th order central TIM  ==! (Eq. 29) 
    409407         ! 
    410          DO jj = 1, jpj 
     408         DO jj = 2, jpjm1 
    411409            DO ji = 1, fs_jpim1   ! vector opt. 
    412410               zcu  = puc(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 
    … …  
    485483      CASE( 1 )                                                !==  1st order central TIM  ==! (Eq. 21) 
    486484         DO jj = 1, jpjm1 
    487             DO ji = 1, jpi 
     485            DO ji = fs_2, fs_jpim1 
    488486               pt_v(ji,jj) = 0.5_wp * vmask(ji,jj,1) * (                              ( pt(ji,jj+1) + pt(ji,jj) )  & 
    489487                  &                                     - SIGN( 1._wp, pvc(ji,jj) ) * ( pt(ji,jj+1) - pt(ji,jj) ) ) 
    … …  
    493491      CASE( 2 )                                                !==  2nd order central TIM  ==! (Eq. 23) 
    494492         DO jj = 1, jpjm1 
    495             DO ji = 1, jpi 
     493            DO ji = fs_2, fs_jpim1 
    496494               zcv  = pvc(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    497495               pt_v(ji,jj) = 0.5_wp * vmask(ji,jj,1) * (        ( pt(ji,jj+1) + pt(ji,jj) )  & 
    … …  
    503501      CASE( 3 )                                                !==  3rd order central TIM  ==! (Eq. 24) 
    504502         DO jj = 1, jpjm1 
    505             DO ji = 1, jpi 
     503            DO ji = fs_2, fs_jpim1 
    506504               zcv  = pvc(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    507505               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
    … …  
    516514      CASE( 4 )                                                !==  4th order central TIM  ==! (Eq. 27) 
    517515         DO jj = 1, jpjm1 
    518             DO ji = 1, jpi 
     516            DO ji = fs_2, fs_jpim1 
    519517               zcv  = pvc(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    520518               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
    … …  
    529527      CASE( 5 )                                                !==  5th order central TIM  ==! (Eq. 29) 
    530528         DO jj = 1, jpjm1 
    531             DO ji = 1, jpi 
     529            DO ji = fs_2, fs_jpim1 
    532530               zcv  = pvc(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 
    533531               zdy2 = e2v(ji,jj) * e2v(ji,jj) 
    … …  
    629627      ! ------------------------------------- 
    630628      DO jj = 2, jpjm1 
    631          DO ji = fs_2, fs_jpim1   ! vector opt. 
     629         DO ji = 1, fs_jpim1   ! vector opt. 
    632630            zau = MIN( 1._wp , zbetdo(ji,jj) , zbetup(ji+1,jj) ) 
    633631            zbu = MIN( 1._wp , zbetup(ji,jj) , zbetdo(ji+1,jj) ) 
    634632            zcu = 0.5  + SIGN( 0.5 , paa(ji,jj) ) 
    635633            ! 
     634            paa(ji,jj) = paa(ji,jj) * ( zcu * zau + ( 1._wp - zcu) * zbu ) 
     635         END DO 
     636      END DO 
     637      ! 
     638      DO jj = 1, jpjm1 
     639         DO ji = fs_2, fs_jpim1   ! vector opt. 
    636640            zav = MIN( 1._wp , zbetdo(ji,jj) , zbetup(ji,jj+1) ) 
    637641            zbv = MIN( 1._wp , zbetup(ji,jj) , zbetdo(ji,jj+1) ) 
    638642            zcv = 0.5  + SIGN( 0.5 , pbb(ji,jj) ) 
    639643            ! 
    640             paa(ji,jj) = paa(ji,jj) * ( zcu * zau + ( 1._wp - zcu) * zbu ) 
    641644            pbb(ji,jj) = pbb(ji,jj) * ( zcv * zav + ( 1._wp - zcv) * zbv ) 
    642             ! 
    643          END DO 
    644       END DO 
    645       CALL lbc_lnk_multi( paa, 'U', -1., pbb, 'V', -1. )   ! lateral boundary condition (changed sign) 
     645         END DO 
     646      END DO 
    646647      ! 
    647648!!      IF( nn_timing == 1 )  CALL timing_stop('nonosc_2d') 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rdgrft.F90

    r8586 r8637  
    253253         CALL tab_3d_2d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
    254254         CALL tab_3d_2d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i(:,:,:) ) 
    255          IF ( nn_pnd_scheme > 0 ) THEN 
    256             CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
    257             CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    258          ENDIF 
     255         CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
     256         CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    259257         DO jl = 1, jpl 
    260258            DO jk = 1, nlay_s 
    … …  
    270268         CALL tab_2d_1d( npti, nptidx(1:npti), hfx_dyn_1d    (1:npti), hfx_dyn    (:,:) ) 
    271269         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_dyn_1d(1:npti), wfx_snw_dyn(:,:) ) 
    272          IF ( nn_pnd_scheme > 0 ) THEN 
    273             CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    274          ENDIF 
     270         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    275271          
    276272         !-----------------------------------------------------------------------------! 
    … …  
    320316         CALL tab_2d_3d( npti, nptidx(1:npti), sv_i_2d(1:npti,1:jpl), sv_i(:,:,:) ) 
    321317         CALL tab_2d_3d( npti, nptidx(1:npti), oa_i_2d(1:npti,1:jpl), oa_i(:,:,:) ) 
    322          IF ( nn_pnd_scheme > 0 ) THEN 
    323             CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
    324             CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    325          ENDIF 
     318         CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 
     319         CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 
    326320         DO jl = 1, jpl 
    327321            DO jk = 1, nlay_s 
    … …  
    337331         CALL tab_1d_2d( npti, nptidx(1:npti), hfx_dyn_1d    (1:npti), hfx_dyn    (:,:) ) 
    338332         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_dyn_1d(1:npti), wfx_snw_dyn(:,:) ) 
    339          IF ( nn_pnd_scheme > 0 ) THEN 
    340             CALL tab_1d_2d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    341          ENDIF 
     333         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_pnd_1d(1:npti), wfx_pnd(:,:) ) 
    342334 
    343335      ENDIF ! npti > 0 
    … …  
    650642               esrft(ji)  = ze_s_2d(ji,1,jl1) * afrft 
    651643 
    652                !MV MP 2016 
    653                IF ( nn_pnd_scheme > 0 ) THEN 
     644               IF ( ln_pnd_H12 ) THEN 
    654645                  aprdg1     = a_ip_2d(ji,jl1) * afrdg 
    655646                  aprdg2(ji) = a_ip_2d(ji,jl1) * afrdg * hi_hrdg(ji,jl1) 
    … …  
    659650                  vprft (ji) = v_ip_2d(ji,jl1) * afrft 
    660651               ENDIF 
    661                ! END MV MP 2016 
    662652 
    663653               !----------------------------------------------------------------- 
    … …  
    679669               ! Put the melt pond water into the ocean 
    680670               !------------------------------------------             
    681                IF ( ( nn_pnd_scheme > 0 ) .AND. ln_pnd_fw ) THEN 
     671               IF ( ln_pnd_fwb ) THEN 
    682672                  wfx_pnd_1d(ji) = wfx_pnd_1d(ji) + ( rhofw * vprdg(ji) * ( 1._wp - rn_fpndrdg )   &        ! fresh water source for ocean 
    683673                     &                              + rhofw * vprft(ji) * ( 1._wp - rn_fpndrft ) ) * r1_rdtice 
    … …  
    701691               sv_i_2d(ji,jl1) = sv_i_2d(ji,jl1) - sirdg1    - sirft(ji) 
    702692               oa_i_2d(ji,jl1) = oa_i_2d(ji,jl1) - oirdg1    - oirft1 
    703                ! MV MP 2016 
    704                IF ( nn_pnd_scheme > 0 ) THEN 
     693               IF ( ln_pnd_H12 ) THEN 
    705694                  a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - aprdg1    - aprft1 
    706695                  v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - vprdg(ji) - vprft(ji) 
    707696               ENDIF 
    708                ! END MV MP 2016 
    709697               ze_s_2d(ji,1,jl1) = ze_s_2d(ji,1,jl1) - esrdg (ji) - esrft (ji) 
    710698            ENDIF 
    … …  
    765753                  v_s_2d (ji,jl2) = v_s_2d (ji,jl2) + ( vsrdg (ji) * rn_fsnwrdg * fvol(ji)  +  & 
    766754                     &                                    vsrft (ji) * rn_fsnwrft * zswitch(ji) ) 
    767                   ! MV MP 2016 
    768                   IF ( nn_pnd_scheme > 0 ) THEN 
     755                  IF ( ln_pnd_H12 ) THEN 
    769756                     v_ip_2d (ji,jl2) = v_ip_2d(ji,jl2) + (   vprdg (ji) * rn_fpndrdg * fvol   (ji)   & 
    770757                        &                                   + vprft (ji) * rn_fpndrft * zswitch(ji)   ) 
    … …  
    772759                        &                                   + aprft2(ji) * rn_fpndrft * zswitch(ji)   ) 
    773760                  ENDIF 
    774                   ! END MV MP 2016 
    775761                  ze_s_2d(ji,1,jl2) = ze_s_2d(ji,1,jl2) + ( esrdg (ji) * rn_fsnwrdg * fvol(ji)  +  & 
    776762                     &                                      esrft (ji) * rn_fsnwrft * zswitch(ji) ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rhg_evp.F90

    r8586 r8637  
    355355         CALL lbc_lnk( zds, 'F', 1. ) 
    356356 
    357          DO jj = 2, jpjm1 
    358             DO ji = 2, jpim1 ! no vector loop 
     357         DO jj = 2, jpj    ! loop to jpi,jpj to avoid making a communication for zs1,zs2,zs12 
     358            DO ji = 2, jpi ! no vector loop 
    359359 
    360360               ! shear**2 at T points (doc eq. A16) 
    … …  
    400400            END DO 
    401401         END DO 
    402          CALL lbc_lnk_multi( zs1, 'T', 1., zs2, 'T', 1., zs12, 'F', 1. ) 
    403   
    404402 
    405403         ! --- Ice internal stresses (Appendix C of Hunke and Dukowicz, 2002) --- ! 
    … …  
    651649         END DO 
    652650      END DO            
    653       CALL lbc_lnk( zds, 'F', 1. ) 
    654651       
    655652      DO jj = 2, jpjm1 
    … …  
    685682       
    686683      ! --- Store the stress tensor for the next time step --- ! 
     684      CALL lbc_lnk_multi( zs1, 'T', 1., zs2, 'T', 1., zs12, 'F', 1. ) 
    687685      pstress1_i (:,:) = zs1 (:,:) 
    688686      pstress2_i (:,:) = zs2 (:,:) 
    … …  
    700698 
    701699      ! --- divergence, shear and strength --- ! 
    702       IF( iom_use('idive')  )   CALL iom_put( "idive"  , pdivu_i(:,:)   * zswi(:,:) )   ! divergence 
    703       IF( iom_use('ishear') )   CALL iom_put( "ishear" , pshear_i(:,:)  * zswi(:,:) )   ! shear 
     700      IF( iom_use('idive')  )   CALL iom_put( "idive"  , pdivu_i(:,:)  * zswi(:,:) )   ! divergence 
     701      IF( iom_use('ishear') )   CALL iom_put( "ishear" , pshear_i(:,:) * zswi(:,:) )   ! shear 
    704702      IF( iom_use('icestr') )   CALL iom_put( "icestr" , strength(:,:) * zswi(:,:) )   ! Ice strength 
    705703 
    … …  
    724722!!                                                                                                               ! (scheme converges if this value is ~1, see Bouillon et al 2009 (eq. 11)) 
    725723               zsig1(ji,jj) = 0.5_wp * zdum2 * ( pstress1_i(ji,jj) )          ! compressive stress, see Bouillon et al. 2015 
    726                zsig2(ji,jj) = 0.5_wp * zdum2 * ( zshear )                    ! shear stress 
     724               zsig2(ji,jj) = 0.5_wp * zdum2 * ( zshear )                     ! shear stress 
    727725               zsig3(ji,jj) = zdum2**2 * ( ( pstress1_i(ji,jj) + strength(ji,jj) )**2 + ( rn_ecc * zshear )**2 ) 
    728726            END DO 
    … …  
    757755                
    758756               ! Stress terms of the momentum equation (N/m2) 
    759                zdiag_dssh_dx(ji,jj) = zspgU(ji,jj) * rswitch    ! sea surface slope stress term 
     757               zdiag_dssh_dx(ji,jj) = zspgU(ji,jj) * rswitch     ! sea surface slope stress term 
    760758               zdiag_dssh_dy(ji,jj) = zspgV(ji,jj) * rswitch 
    761759                
    762                zdiag_corstrx(ji,jj) = zCorx(ji,jj) * rswitch    ! Coriolis stress term 
     760               zdiag_corstrx(ji,jj) = zCorx(ji,jj) * rswitch     ! Coriolis stress term 
    763761               zdiag_corstry(ji,jj) = zCory(ji,jj) * rswitch 
    764762                
    765                zdiag_intstrx(ji,jj) = zfU(ji,jj)   * rswitch    ! internal stress term 
     763               zdiag_intstrx(ji,jj) = zfU(ji,jj)   * rswitch     ! internal stress term 
    766764               zdiag_intstry(ji,jj) = zfV(ji,jj)   * rswitch 
    767765                

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/iceforcing.F90

    r8586 r8637  
    130130 
    131131      ! --- cloud-sky and overcast-sky ice albedos --- ! 
    132       CALL ice_alb( t_su, h_i, h_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os ) 
     132      CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_frac, h_ip, zalb_cs, zalb_os ) 
    133133 
    134134      ! albedo depends on cloud fraction because of non-linear spectral effects 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/iceistate.F90

    r8586 r8637  
    9494      REAL(wp) ::   ztmelts, zdh 
    9595      INTEGER  ::   i_hemis, i_fill, jl0 
    96       REAL(wp) ::   zarg, zV, zconv, zdv 
     96      REAL(wp) ::   zarg, zV, zconv, zdv, zfac 
    9797      INTEGER , DIMENSION(4)           ::   itest 
    9898      REAL(wp), DIMENSION(jpi,jpj)     ::   z2d 
    … …  
    314314 
    315315         ! for constant salinity in time 
    316          IF( nn_icesal == 1 .OR. nn_icesal == 3 )  THEN 
     316         IF( nn_icesal /= 2 )  THEN 
    317317            CALL ice_var_salprof 
    318318            sv_i = s_i * v_i 
    … …  
    358358         tn_ice (:,:,:) = t_su (:,:,:) 
    359359 
    360          ! MV MP 2016 
    361360         ! Melt pond volume and fraction 
    362          IF ( ln_pnd ) THEN 
    363             DO jl = 1, jpl 
    364                a_ip_frac(:,:,jl) = 0.2 *  zswitch(:,:) 
    365                h_ip     (:,:,jl) = 0.05 * zswitch(:,:) 
    366                a_ip(:,:,jl)      = a_ip_frac(:,:,jl) * a_i (:,:,jl)  
    367                v_ip(:,:,jl)      = h_ip     (:,:,jl) * a_ip(:,:,jl) 
    368             END DO 
    369          ELSE 
    370             a_ip(:,:,:)      = 0._wp 
    371             v_ip(:,:,:)      = 0._wp 
    372             a_ip_frac(:,:,:) = 0._wp 
    373             h_ip     (:,:,:) = 0._wp 
    374          ENDIF 
    375          ! END MV MP 2016 
     361         IF ( ln_pnd_CST .OR. ln_pnd_H12 ) THEN   ;   zfac = 1._wp 
     362         ELSE                                     ;   zfac = 0._wp 
     363         ENDIF  
     364         DO jl = 1, jpl 
     365            a_ip_frac(:,:,jl) = rn_apnd * zswitch(:,:) * zfac 
     366            h_ip     (:,:,jl) = rn_hpnd * zswitch(:,:) * zfac 
     367         END DO 
     368         a_ip(:,:,:) = a_ip_frac(:,:,:) * a_i (:,:,:)  
     369         v_ip(:,:,:) = h_ip     (:,:,:) * a_ip(:,:,:) 
    376370 
    377371      ELSE ! if ln_iceini=false 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90

    r8586 r8637  
    293293            IF ( a_i_1d(ji) > epsi10 .AND. h_i_1d(ji) < rn_himin ) THEN 
    294294               a_i_1d (ji) = a_i_1d(ji) * h_i_1d(ji) / rn_himin  
    295                ! MV MP 2016 
    296                IF ( nn_pnd_scheme > 0 ) THEN 
    297                   a_ip_1d(ji) = a_ip_1d(ji) * h_i_1d(ji) / rn_himin 
    298                ENDIF 
    299                ! END MV MP 2016 
     295               IF ( ln_pnd_H12 )    a_ip_1d(ji) = a_ip_1d(ji) * h_i_1d(ji) / rn_himin 
    300296               h_i_1d(ji) = rn_himin 
    301297            ENDIF 
    … …  
    457453               zaTsfn(ji,jl2)  = zaTsfn(ji,jl2) + ztrans 
    458454               !   
    459                ! MV MP 2016  
    460                IF ( nn_pnd_scheme > 0 ) THEN 
     455               IF ( ln_pnd_H12 ) THEN 
    461456                  !                                                  ! Pond fraction 
    462457                  ztrans          = a_ip_2d(ji,jl1) * pdaice(ji,jl) !!clem: should be * zworka(ji) but it does not work 
    … …  
    468463                  v_ip_2d(ji,jl2) = v_ip_2d(ji,jl2) + ztrans 
    469464               ENDIF 
    470                ! END MV MP 2016 
    471465               ! 
    472466            ENDIF   ! jl1 >0 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icerst.F90

    r8586 r8637  
    147147      END DO 
    148148 
    149       ! MV MP 2016 
    150       IF ( nn_pnd_scheme > 0 ) THEN 
    151          DO jl = 1, jpl  
    152             WRITE(zchar,'(I2.2)') jl 
    153             znam = 'a_ip'//'_htc'//zchar 
    154             z2d(:,:) = a_ip(:,:,jl) 
    155             CALL iom_rstput( iter, nitrst, numriw, znam , z2d ) ! a_ip 
    156             znam = 'v_ip'//'_htc'//zchar 
    157             z2d(:,:) = v_ip(:,:,jl) 
    158             CALL iom_rstput( iter, nitrst, numriw, znam , z2d ) ! v_ip 
    159          END DO 
    160       ENDIF 
    161       ! END MV MP 2016 
     149      DO jl = 1, jpl  
     150         WRITE(zchar,'(I2.2)') jl 
     151         znam = 'a_ip'//'_htc'//zchar 
     152         z2d(:,:) = a_ip(:,:,jl) 
     153         CALL iom_rstput( iter, nitrst, numriw, znam , z2d ) ! a_ip 
     154         znam = 'v_ip'//'_htc'//zchar 
     155         z2d(:,:) = v_ip(:,:,jl) 
     156         CALL iom_rstput( iter, nitrst, numriw, znam , z2d ) ! v_ip 
     157      END DO 
    162158 
    163159      DO jl = 1, jpl  
    … …  
    193189      !! ** purpose  :   read restart file 
    194190      !!---------------------------------------------------------------------- 
    195       INTEGER  :: jk, jl 
    196       REAL(wp) :: zfice, ziter 
     191      INTEGER  ::   jk, jl 
     192      INTEGER  ::   id1            ! local integer 
     193      REAL(wp) ::   zfice, ziter 
    197194      REAL(wp), DIMENSION(jpi,jpj) ::   z2d 
    198195      CHARACTER(len=25) ::   znam 
    … …  
    248245      END DO 
    249246 
    250       ! MV MP 2016 
    251       IF ( nn_pnd_scheme > 0 ) THEN 
     247      id1 = iom_varid( numrir, 'a_ip_htc01' , ldstop = .FALSE. ) 
     248      IF( id1 > 0 ) THEN                       ! fields exist (melt ponds) 
    252249         DO jl = 1, jpl  
    253250            WRITE(zchar,'(I2.2)') jl 
    … …  
    259256            v_ip(:,:,jl) = z2d(:,:) 
    260257         END DO 
    261       ENDIF 
    262       ! END MV MP 2016 
     258      ELSE                                     ! start from rest 
     259         IF(lwp) WRITE(numout,*) '   ==>>   previous run without melt ponds output then set it' 
     260         a_ip(:,:,:) = 0._wp 
     261         v_ip(:,:,:) = 0._wp 
     262      ENDIF 
    263263 
    264264      DO jl = 1, jpl  

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90

    r8586 r8637  
    3535   USE icedyn         ! sea-ice: dynamics 
    3636   USE icethd         ! sea-ice: thermodynamics 
    37    USE limmp          ! sea-ice: melt ponds 
    3837   USE icecor         ! sea-ice: corrections 
    3938   USE iceupdate      ! sea-ice: sea surface boundary condition update 
    … …  
    170169         IF( ln_icethd )                CALL ice_thd( kt )            ! -- Ice thermodynamics       
    171170         ! 
    172          IF ( ln_pnd )                  CALL lim_mp( kt )             ! -- Melt ponds 
    173171         ! 
    174172         IF( ln_icethd )                CALL ice_cor( kt , 2 )        ! -- Corrections 
    … …  
    238236      CALL ice_itd_init                ! ice thickness distribution initialization 
    239237      ! 
    240       CALL lim_mp_init                 ! set melt ponds parameters (clem: important to be located here) 
    241       ! 
     238      IF( ln_icethd ) THEN 
     239         CALL ice_thd_init             ! set ice thermodynics parameters (clem: important to call it first for melt ponds) 
     240      ENDIF    
    242241      !                                ! Initial sea-ice state 
    243242      IF( .NOT. ln_rstart ) THEN              ! start from rest: sea-ice deduced from sst 
    … …  
    255254         CALL ice_dyn_init             ! set ice dynamics parameters 
    256255      ENDIF 
    257       ! 
    258       IF( ln_icethd ) THEN 
    259          CALL ice_thd_init             ! set ice thermodynics parameters 
    260       ENDIF    
    261256      ! 
    262257      CALL ice_update_init             ! ice surface boundary condition 
    … …  
    398393      wfx_snw_sub(:,:) = 0._wp ; wfx_ice_sub(:,:) = 0._wp 
    399394      wfx_snw_sni(:,:) = 0._wp  
    400       ! MV MP 2016 
    401395      wfx_pnd(:,:) = 0._wp 
    402       ! END MV MP 2016 
    403396 
    404397      hfx_thd(:,:) = 0._wp   ; 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icethd.F90

    r8586 r8637  
    3434   USE icethd_ent     ! sea-ice: enthalpy redistribution 
    3535   USE icethd_do      ! sea-ice: growth in open water 
     36   USE icethd_pnd     ! sea-ice: melt ponds 
    3637   USE iceitd         ! sea-ice: remapping thickness distribution 
    3738   USE icetab         ! sea-ice: 1D <==> 2D transformation 
    … …  
    8687      !!             - call ice_thd_rem  for remapping thickness distribution 
    8788      !!             - call ice_thd_do   for ice growth in leads 
    88       !!--------------------------------------------------------------------- 
     89      !!------------------------------------------------------------------- 
    8990      INTEGER, INTENT(in) :: kt    ! number of iteration 
    9091      ! 
    … …  
    230231            s_i_new   (1:npti) = 0._wp ; dh_s_tot (1:npti) = 0._wp  ! --- some init --- !  (important to have them here)  
    231232            dh_i_surf (1:npti) = 0._wp ; dh_i_bott(1:npti) = 0._wp 
    232             dh_snowice(1:npti) = 0._wp ; dh_i_sub (1:npti) = 0._wp 
     233            dh_snowice(1:npti) = 0._wp ; dh_i_sub (1:npti) = 0._wp ; dh_s_mlt(1:npti) = 0._wp 
    233234            ! 
    234235            IF( ln_icedH ) THEN                                     ! --- growing/melting --- ! 
    235236                              CALL ice_thd_zdf                             ! Ice/Snow Temperature profile 
    236237                              CALL ice_thd_dh                              ! Ice/Snow thickness    
     238                              CALL ice_thd_pnd                             ! Melt ponds formation 
    237239                              CALL ice_thd_ent( e_i_1d(1:npti,:) )         ! Ice enthalpy remapping 
    238240            ENDIF 
    … …  
    362364         CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti), at_i             ) 
    363365         CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,kl)     ) 
    364          CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d(1:npti), h_i(:,:,kl)     ) 
    365          CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d(1:npti), h_s(:,:,kl)     ) 
     366         CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i (:,:,kl)     ) 
     367         CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d (1:npti), h_s (:,:,kl)     ) 
    366368         CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti), t_su(:,:,kl)     ) 
    367          CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d(1:npti), s_i(:,:,kl)     ) 
     369         CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d (1:npti), s_i (:,:,kl)     ) 
    368370         DO jk = 1, nlay_s 
    369             CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl)   ) 
    370             CALL tab_2d_1d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl)   ) 
     371            CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl)    ) 
     372            CALL tab_2d_1d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl)    ) 
    371373         END DO 
    372374         DO jk = 1, nlay_i 
    373             CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d(1:npti,jk), t_i(:,:,jk,kl)   ) 
    374             CALL tab_2d_1d( npti, nptidx(1:npti), e_i_1d(1:npti,jk), e_i(:,:,jk,kl)   ) 
    375             CALL tab_2d_1d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl)   ) 
    376          END DO 
     375            CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d (1:npti,jk), t_i (:,:,jk,kl)  ) 
     376            CALL tab_2d_1d( npti, nptidx(1:npti), e_i_1d (1:npti,jk), e_i (:,:,jk,kl)  ) 
     377            CALL tab_2d_1d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl)  ) 
     378         END DO 
     379         CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d     (1:npti), a_ip     (:,:,kl) ) 
     380         CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d     (1:npti), h_ip     (:,:,kl) ) 
     381         CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_frac_1d(1:npti), a_ip_frac(:,:,kl) ) 
    377382         ! 
    378383         CALL tab_2d_1d( npti, nptidx(1:npti), qprec_ice_1d(1:npti), qprec_ice        ) 
    … …  
    403408         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr          ) 
    404409         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam          ) 
     410         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd          ) 
    405411         ! 
    406412         CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog          ) 
    … …  
    454460         ! 
    455461         ! Change thickness to volume (replaces routine ice_var_eqv2glo) 
    456          v_i_1d(1:npti)  = h_i_1d(1:npti) * a_i_1d(1:npti) 
    457          v_s_1d(1:npti)  = h_s_1d(1:npti) * a_i_1d(1:npti) 
    458          sv_i_1d(1:npti) = s_i_1d(1:npti) * v_i_1d(1:npti) 
     462         v_i_1d (1:npti) = h_i_1d (1:npti) * a_i_1d (1:npti) 
     463         v_s_1d (1:npti) = h_s_1d (1:npti) * a_i_1d (1:npti) 
     464         sv_i_1d(1:npti) = s_i_1d (1:npti) * v_i_1d (1:npti) 
     465         v_ip_1d(1:npti) = h_ip_1d(1:npti) * a_ip_1d(1:npti) 
    459466          
    460467         CALL tab_1d_2d( npti, nptidx(1:npti), at_i_1d(1:npti), at_i             ) 
    461468         CALL tab_1d_2d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,kl)     ) 
    462          CALL tab_1d_2d( npti, nptidx(1:npti), h_i_1d(1:npti), h_i(:,:,kl)     ) 
    463          CALL tab_1d_2d( npti, nptidx(1:npti), h_s_1d(1:npti), h_s(:,:,kl)     ) 
     469         CALL tab_1d_2d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i (:,:,kl)     ) 
     470         CALL tab_1d_2d( npti, nptidx(1:npti), h_s_1d (1:npti), h_s (:,:,kl)     ) 
    464471         CALL tab_1d_2d( npti, nptidx(1:npti), t_su_1d(1:npti), t_su(:,:,kl)     ) 
    465          CALL tab_1d_2d( npti, nptidx(1:npti), s_i_1d(1:npti), s_i(:,:,kl)     ) 
     472         CALL tab_1d_2d( npti, nptidx(1:npti), s_i_1d (1:npti), s_i (:,:,kl)     ) 
    466473         DO jk = 1, nlay_s 
    467             CALL tab_1d_2d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl) ) 
    468             CALL tab_1d_2d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl) ) 
     474            CALL tab_1d_2d( npti, nptidx(1:npti), t_s_1d(1:npti,jk), t_s(:,:,jk,kl)    ) 
     475            CALL tab_1d_2d( npti, nptidx(1:npti), e_s_1d(1:npti,jk), e_s(:,:,jk,kl)    ) 
    469476         END DO 
    470477         DO jk = 1, nlay_i 
    471             CALL tab_1d_2d( npti, nptidx(1:npti), t_i_1d(1:npti,jk), t_i(:,:,jk,kl) ) 
    472             CALL tab_1d_2d( npti, nptidx(1:npti), e_i_1d(1:npti,jk), e_i(:,:,jk,kl) ) 
    473             CALL tab_1d_2d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl) ) 
    474          END DO 
     478            CALL tab_1d_2d( npti, nptidx(1:npti), t_i_1d (1:npti,jk), t_i (:,:,jk,kl)  ) 
     479            CALL tab_1d_2d( npti, nptidx(1:npti), e_i_1d (1:npti,jk), e_i (:,:,jk,kl)  ) 
     480            CALL tab_1d_2d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl)  ) 
     481         END DO 
     482         CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_1d     (1:npti), a_ip     (:,:,kl) ) 
     483         CALL tab_1d_2d( npti, nptidx(1:npti), h_ip_1d     (1:npti), h_ip     (:,:,kl) ) 
     484         CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_frac_1d(1:npti), a_ip_frac(:,:,kl) ) 
    475485         ! 
    476486         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sni_1d(1:npti), wfx_snw_sni ) 
    … …  
    488498         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr        ) 
    489499         CALL tab_1d_2d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam        ) 
     500         CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd        ) 
    490501         ! 
    491502         CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog        ) 
    … …  
    523534         CALL tab_1d_2d( npti, nptidx(1:npti), v_s_1d (1:npti), v_s (:,:,kl) ) 
    524535         CALL tab_1d_2d( npti, nptidx(1:npti), sv_i_1d(1:npti), sv_i(:,:,kl) ) 
     536         CALL tab_1d_2d( npti, nptidx(1:npti), v_ip_1d(1:npti), v_ip(:,:,kl) ) 
    525537         ! 
    526538      END SELECT 
    … …  
    530542 
    531543   SUBROUTINE ice_thd_init 
    532       !!----------------------------------------------------------------------- 
     544      !!------------------------------------------------------------------- 
    533545      !!                   ***  ROUTINE ice_thd_init ***  
    534546      !!                  
    … …  
    570582      IF( ln_icedO )   CALL ice_thd_do_init    ! set ice growth in open water parameters 
    571583                       CALL ice_thd_sal_init   ! set ice salinity parameters 
    572       ! 
    573       IF( ln_icedS .AND. nn_icesal == 1 ) THEN 
    574          ln_icedS = .FALSE. 
    575          CALL ctl_warn('ln_icedS is set to false since constant ice salinity is chosen (nn_icesal=1)') 
    576       ENDIF 
     584                       CALL ice_thd_pnd_init   ! set melt ponds parameters 
    577585      ! 
    578586   END SUBROUTINE ice_thd_init 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90

    r8586 r8637  
    148148            ! Contribution to mass flux 
    149149            wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhosn * h_s_1d(ji) * a_i_1d(ji) * r1_rdtice 
     150            dh_s_mlt(ji)       = dh_s_mlt(ji) - h_s_1d(ji) 
    150151            ! updates 
    151             h_s_1d(ji)   = 0._wp 
     152            h_s_1d(ji)    = 0._wp 
    152153            e_s_1d (ji,1) = 0._wp 
    153154            t_s_1d (ji,1) = rt0 
    … …  
    211212         ! snow melting only = water into the ocean (then without snow precip), >0 
    212213         wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice 
     214         dh_s_mlt(ji)       = dh_s_mlt(ji) + zdeltah(ji,1) 
    213215         ! updates available heat + precipitations after melting 
    214216         zq_su     (ji) = MAX( 0._wp , zq_su (ji) + zdeltah(ji,1) * zqprec(ji) )       
    … …  
    233235            hfx_snw_1d(ji)   = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_1d(ji) * e_s_1d(ji,jk) * r1_rdtice  
    234236            ! snow melting only = water into the ocean (then without snow precip) 
    235             wfx_snw_sum_1d(ji)   = wfx_snw_sum_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice 
     237            wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice 
     238            dh_s_mlt(ji)       = dh_s_mlt(ji) + zdeltah(ji,jk) 
    236239            ! updates available heat + thickness 
    237             zq_su (ji)  = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) ) 
     240            zq_su (ji) = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) ) 
    238241            h_s_1d(ji) = MAX( 0._wp , h_s_1d(ji) + zdeltah(ji,jk) ) 
    239242         END DO 
    … …  
    571574         zdeltah  (ji,1) = MIN( 0._wp , MAX( zdeltah(ji,1) , - h_s_1d(ji) ) ) ! bound melting 
    572575         dh_s_tot (ji)   = dh_s_tot(ji)  + zdeltah(ji,1) 
    573          h_s_1d   (ji)  = h_s_1d(ji)   + zdeltah(ji,1) 
     576         h_s_1d   (ji)   = h_s_1d(ji)   + zdeltah(ji,1) 
    574577         
    575578         zq_rema(ji)     = zq_rema(ji) + zdeltah(ji,1) * e_s_1d(ji,1)                ! update available heat (J.m-2) 
    … …  
    577580         hfx_snw_1d(ji)  = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_1d(ji) * e_s_1d(ji,1) * r1_rdtice ! W.m-2 (>0) 
    578581         ! Contribution to mass flux 
    579          wfx_snw_sum_1d(ji)  =  wfx_snw_sum_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice 
     582         wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhosn * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice 
     583         dh_s_mlt(ji)       = dh_s_mlt(ji) + zdeltah(ji,1) 
    580584         !     
    581585         ! Remaining heat flux (W.m-2) is sent to the ocean heat budget 
    … …  
    611615 
    612616         ! Case constant salinity in time: virtual salt flux to keep salinity constant 
    613          IF( nn_icesal == 1 .OR. nn_icesal == 3 )  THEN 
     617         IF( nn_icesal /= 2 )  THEN 
    614618            sfx_bri_1d(ji) = sfx_bri_1d(ji) - sss_1d (ji) * a_i_1d(ji) * zfmdt                  * r1_rdtice  & ! put back sss_m     into the ocean 
    615619               &                            - s_i_1d(ji)  * a_i_1d(ji) * dh_snowice(ji) * rhoic * r1_rdtice    ! and get  rn_icesal from the ocean  

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icethd_zdf.F90

    r8586 r8637  
    3131 
    3232   !!** namelist (namthd_zdf) ** 
    33    LOGICAL  ::   ln_zdf_Beer      ! Heat diffusion follows a Beer Law 
     33   LOGICAL  ::   ln_zdf_BL99      ! Heat diffusion follows Bitz and Lipscomb (1999) 
    3434   LOGICAL  ::   ln_cndi_U64      ! thermal conductivity: Untersteiner (1964) 
    3535   LOGICAL  ::   ln_cndi_P07      ! thermal conductivity: Pringle et al (2007) 
    … …  
    677677      INTEGER  ::   ios   ! Local integer output status for namelist read 
    678678      !! 
    679       NAMELIST/namthd_zdf/ ln_zdf_Beer, ln_cndi_U64, ln_cndi_P07, rn_cnd_s, rn_kappa_i, ln_dqns_i 
     679      NAMELIST/namthd_zdf/ ln_zdf_BL99, ln_cndi_U64, ln_cndi_P07, rn_cnd_s, rn_kappa_i, ln_dqns_i 
    680680      !!------------------------------------------------------------------- 
    681681      ! 
    … …  
    694694         WRITE(numout,*) '~~~~~~~~~~~~~~~~' 
    695695         WRITE(numout,*) '   Namelist namthd_zdf:' 
    696          WRITE(numout,*) '      Diffusion follows a Beer Law                            ln_zdf_Beer  = ', ln_zdf_Beer 
     696         WRITE(numout,*) '      Diffusion follows a Bitz and Lipscomb (1999)            ln_zdf_BL99  = ', ln_zdf_BL99 
    697697         WRITE(numout,*) '      thermal conductivity in the ice (Untersteiner 1964)     ln_cndi_U64  = ', ln_cndi_U64 
    698698         WRITE(numout,*) '      thermal conductivity in the ice (Pringle et al 2007)    ln_cndi_P07  = ', ln_cndi_P07 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/iceupdate.F90

    r8586 r8637  
    164164            ! mass flux from ice/ocean 
    165165            wfx_ice(ji,jj) = wfx_bog(ji,jj) + wfx_bom(ji,jj) + wfx_sum(ji,jj) + wfx_sni(ji,jj)   & 
    166                &           + wfx_opw(ji,jj) + wfx_dyn(ji,jj) + wfx_res(ji,jj) + wfx_lam(ji,jj)  
    167  
    168             IF ( ln_pnd_fw )   wfx_ice(ji,jj) = wfx_ice(ji,jj) + wfx_pnd(ji,jj) 
     166               &           + wfx_opw(ji,jj) + wfx_dyn(ji,jj) + wfx_res(ji,jj) + wfx_lam(ji,jj) + wfx_pnd(ji,jj) 
    169167 
    170168            ! add the snow melt water to snow mass flux to the ocean 
    … …  
    199197      ! Snow/ice albedo (only if sent to coupler, useless in forced mode) 
    200198      !------------------------------------------------------------------ 
    201       CALL ice_alb( t_su, h_i, h_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 
     199      CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_frac, h_ip, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 
    202200      ! 
    203201      alb_ice(:,:,:) = ( 1._wp - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
    … …  
    246244      CALL iom_put( "vfxlam"     , wfx_lam              )        ! lateral melt  
    247245      CALL iom_put( "vfxice"     , wfx_ice              )        ! total ice growth/melt  
    248       IF ( ln_pnd )   CALL iom_put( "vfxpnd", wfx_pnd   )        ! melt pond water flux 
     246      IF ( ln_pnd_fwb )   CALL iom_put( "vfxpnd", wfx_pnd   )        ! melt pond water flux 
    249247 
    250248      IF ( iom_use( "vfxthin" ) ) THEN   ! ice production for open water + thin ice (<20cm) => comparable to observations   

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90

    r8586 r8637  
    9797      et_i(:,:)  = SUM( SUM( e_i(:,:,:,:), dim=4 ), dim=3 ) 
    9898 
    99       ! MV MP 2016 
    100       IF ( ln_pnd ) THEN                     ! Melt pond 
    101          at_ip(:,:) = SUM( a_ip(:,:,:), dim=3 ) 
    102          vt_ip(:,:) = SUM( v_ip(:,:,:), dim=3 ) 
    103       ENDIF 
    104       ! END MP 2016 
     99      at_ip(:,:) = SUM( a_ip(:,:,:), dim=3 ) ! melt ponds 
     100      vt_ip(:,:) = SUM( v_ip(:,:,:), dim=3 ) 
    105101 
    106102      ato_i(:,:) = 1._wp - at_i(:,:)    ! open water fraction   
    … …  
    134130         END DO 
    135131         ! 
     132         !                           ! put rt0 where there is no ice 
     133         WHERE( at_i(:,:)<=epsi20 ) 
     134            tm_su(:,:) = rt0 
     135            tm_si(:,:) = rt0 
     136            tm_i (:,:) = rt0 
     137         END WHERE 
     138 
    136139         DEALLOCATE( z1_at_i , z1_vt_i ) 
    137140      ENDIF 
    … …  
    160163!!                a criteria for icy area (i.e. a_i > epsi20 and v_i > epsi20 ) 
    161164 
    162       !------------------------------------------------------- 
    163       ! Ice thickness, snow thickness, ice salinity, ice age 
    164       !------------------------------------------------------- 
     165      !--------------------------------------------------------------- 
     166      ! Ice thickness, snow thickness, ice salinity, ice age and ponds 
     167      !--------------------------------------------------------------- 
    165168      !                                            !--- inverse of the ice area 
    166169      WHERE( a_i(:,:,:) > epsi20 )   ;   z1_a_i(:,:,:) = 1._wp / a_i(:,:,:) 
    … …  
    171174      ELSEWHERE                      ;   z1_v_i(:,:,:) = 0._wp 
    172175      END WHERE 
    173       ! 
    174       h_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:)    !--- ice thickness 
     176      !                                           !--- ice thickness 
     177      h_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:) 
    175178 
    176179      zhmax    =          hi_max(jpl) 
    177180      z1_zhmax =  1._wp / hi_max(jpl)                
    178       WHERE( h_i(:,:,jpl) > zhmax )               !--- bound h_i by hi_max (i.e. 99 m) with associated update of ice area 
     181      WHERE( h_i(:,:,jpl) > zhmax )   ! bound h_i by hi_max (i.e. 99 m) with associated update of ice area 
    179182         h_i  (:,:,jpl) = zhmax 
    180183         a_i   (:,:,jpl) = v_i(:,:,jpl) * z1_zhmax  
    181          z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl)          ! NB: v_i always /=0 as h_i > hi_max 
     184         z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl) 
    182185      END WHERE 
    183  
    184       h_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:)    !--- snow thickness 
    185        
    186       o_i(:,:,:)  = oa_i(:,:,:) * z1_a_i(:,:,:)    !--- ice age 
    187  
    188       IF( nn_icesal == 2 ) THEN                    !--- salinity (with a minimum value imposed everywhere) 
     186      !                                           !--- snow thickness 
     187      h_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:) 
     188      !                                           !--- ice age       
     189      o_i(:,:,:) = oa_i(:,:,:) * z1_a_i(:,:,:) 
     190      !                                           !--- pond fraction and thickness       
     191      a_ip_frac(:,:,:) = a_ip(:,:,:) * z1_a_i(:,:,:) 
     192      WHERE( a_ip_frac(:,:,:) > epsi20 )   ;   h_ip(:,:,:) = v_ip(:,:,:) * z1_a_i(:,:,:) / a_ip_frac(:,:,:) 
     193      ELSEWHERE                            ;   h_ip(:,:,:) = 0._wp 
     194      END WHERE 
     195      ! 
     196      !                                           !---  salinity (with a minimum value imposed everywhere)      
     197      IF( nn_icesal == 2 ) THEN 
    189198         WHERE( v_i(:,:,:) > epsi20 )   ;   s_i(:,:,:) = MAX( rn_simin , MIN( rn_simax, sv_i(:,:,:) * z1_v_i(:,:,:) ) ) 
    190199         ELSEWHERE                      ;   s_i(:,:,:) = rn_simin 
    191200         END WHERE 
    192201      ENDIF 
    193  
    194       CALL ice_var_salprof      ! salinity profile 
     202      CALL ice_var_salprof   ! salinity profile 
    195203 
    196204      !------------------- 
    … …  
    235243      vt_s (:,:) = SUM( v_s, dim=3 ) 
    236244      at_i (:,:) = SUM( a_i, dim=3 ) 
    237  
    238 ! MV MP 2016 
    239 ! probably should resum for melt ponds ??? 
    240  
    241245      ! 
    242246   END SUBROUTINE ice_var_glo2eqv 
    … …  
    251255      !!------------------------------------------------------------------- 
    252256      ! 
    253       v_i (:,:,:) = h_i(:,:,:) * a_i(:,:,:) 
    254       v_s (:,:,:) = h_s(:,:,:) * a_i(:,:,:) 
    255       sv_i(:,:,:) = s_i(:,:,:) * v_i(:,:,:) 
     257      v_i (:,:,:) = h_i (:,:,:) * a_i (:,:,:) 
     258      v_s (:,:,:) = h_s (:,:,:) * a_i (:,:,:) 
     259      sv_i(:,:,:) = s_i (:,:,:) * v_i (:,:,:) 
     260      v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 
    256261      ! 
    257262   END SUBROUTINE ice_var_eqv2glo 
    … …  
    471476               wfx_res(ji,jj)  = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl)   * rhosn * r1_rdtice 
    472477               hfx_res(ji,jj)  = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_s (ji,jj,1,jl) * r1_rdtice ! W.m-2 <0 
     478               IF( ln_pnd_fwb ) THEN  
     479                  wfx_res(ji,jj)  = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_ip(ji,jj,jl) * rhofw * r1_rdtice 
     480               ENDIF 
    473481               !----------------------------------------------------------------- 
    474482               ! Zap snow energy  
    … …  
    491499               h_s (ji,jj,jl) = h_s (ji,jj,jl) * zswitch(ji,jj) 
    492500 
    493             END DO 
    494          END DO 
    495  
    496          IF( ln_pnd ) THEN  
    497             DO jj = 1 , jpj 
    498                DO ji = 1 , jpi 
    499                   IF( ln_pnd_fw )   & 
    500                      &   wfx_res(ji,jj)  = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_ip(ji,jj,jl) * rhofw * r1_rdtice 
    501                   a_ip (ji,jj,jl) = a_ip (ji,jj,jl) * zswitch(ji,jj) 
    502                   v_ip (ji,jj,jl) = v_ip (ji,jj,jl) * zswitch(ji,jj) 
    503                END DO 
    504             END DO 
    505          ENDIF 
    506           
     501               a_ip (ji,jj,jl) = a_ip (ji,jj,jl) * zswitch(ji,jj) 
     502               v_ip (ji,jj,jl) = v_ip (ji,jj,jl) * zswitch(ji,jj) 
     503 
     504            END DO 
     505         END DO 
     506         
    507507      END DO  
    508508 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icewri.F90

    r8586 r8637  
    128128      CALL iom_put( "snowvol"     , vt_s  * zswi      )        ! snow volume 
    129129       
    130       IF ( ln_pnd ) THEN 
    131          CALL iom_put( "iceamp"  , at_ip  * zswi        )   ! melt pond total fraction 
    132          CALL iom_put( "icevmp"  , vt_ip  * zswi        )   ! melt pond total volume per unit area 
    133       ENDIF 
     130      CALL iom_put( "iceamp"  , at_ip  * zswi        )   ! melt pond total fraction 
     131      CALL iom_put( "icevmp"  , vt_ip  * zswi        )   ! melt pond total volume per unit area 
    134132 
    135133      !---------------------------------- 
    … …  
    145143      IF ( iom_use('brinevol_cat') )  CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 )          ! brine volume 
    146144 
    147       IF ( ln_pnd ) THEN 
    148          IF ( iom_use('iceamp_cat') )  CALL iom_put( "iceamp_cat"     , a_ip       * zswi2   )       ! melt pond frac for categories 
    149          IF ( iom_use('icevmp_cat') )  CALL iom_put( "icevmp_cat"     , v_ip       * zswi2   )       ! melt pond frac for categories 
    150          IF ( iom_use('icehmp_cat') )  CALL iom_put( "icehmp_cat"     , h_ip       * zswi2   )       ! melt pond frac for categories 
    151          IF ( iom_use('iceafp_cat') )  CALL iom_put( "iceafp_cat"     , a_ip_frac  * zswi2   )       ! melt pond frac for categories 
    152       ENDIF 
     145      IF ( iom_use('iceamp_cat') )  CALL iom_put( "iceamp_cat"     , a_ip       * zswi2   )       ! melt pond frac for categories 
     146      IF ( iom_use('icevmp_cat') )  CALL iom_put( "icevmp_cat"     , v_ip       * zswi2   )       ! melt pond frac for categories 
     147      IF ( iom_use('icehmp_cat') )  CALL iom_put( "icehmp_cat"     , h_ip       * zswi2   )       ! melt pond frac for categories 
     148      IF ( iom_use('iceafp_cat') )  CALL iom_put( "iceafp_cat"     , a_ip_frac  * zswi2   )       ! melt pond frac for categories 
    153149 
    154150      !-------------------------------- 
    … …  
    316312      CALL histdef( kid, "sidive", "Ice divergence"          , "10-8s-1",   & 
    317313      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )  
    318  
    319       ! MV MP 2016 
    320       IF ( ln_pnd ) THEN 
    321          CALL histdef( kid, "si_amp", "Melt pond fraction"      , "%"      ,   & 
    322       &         jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
    323          CALL histdef( kid, "si_vmp", "Melt pond volume"        ,  "m"     ,   & 
    324       &         jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
    325       ENDIF 
    326       ! END MV MP 2016 
    327  
     314      CALL histdef( kid, "si_amp", "Melt pond fraction"      , "%"      ,   & 
     315      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
     316      CALL histdef( kid, "si_vmp", "Melt pond volume"        ,  "m"     ,   & 
     317      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
    328318      CALL histdef( kid, "vfxbog", "Ice bottom production"   , "m/s"    ,   & 
    329319      &      jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 
    … …  
    370360      CALL histwrite( kid, "sidive", kt, divu_i*1.0e8   , jpi*jpj, (/1/) ) 
    371361 
    372       ! MV MP 2016 
    373       IF ( ln_pnd ) THEN 
    374          CALL histwrite( kid, "si_amp", kt, at_ip         , jpi*jpj, (/1/) ) 
    375          CALL histwrite( kid, "si_vmp", kt, vt_ip         , jpi*jpj, (/1/) ) 
    376       ENDIF 
    377       ! END MV MP 2016 
     362      CALL histwrite( kid, "si_amp", kt, at_ip         , jpi*jpj, (/1/) ) 
     363      CALL histwrite( kid, "si_vmp", kt, vt_ip         , jpi*jpj, (/1/) ) 
    378364 
    379365      CALL histwrite( kid, "vfxbog", kt, wfx_bog        , jpi*jpj, (/1/) ) 
    … …  
    384370      CALL histwrite( kid, "vfxbom", kt, wfx_bom        , jpi*jpj, (/1/) ) 
    385371      CALL histwrite( kid, "vfxsum", kt, wfx_sum        , jpi*jpj, (/1/) ) 
    386       IF ( ln_pnd ) & 
    387          CALL histwrite( kid, "vfxpnd", kt, wfx_pnd     , jpi*jpj, (/1/) ) 
     372      CALL histwrite( kid, "vfxpnd", kt, wfx_pnd        , jpi*jpj, (/1/) ) 
    388373 
    389374      CALL histwrite( kid, "sithicat", kt, h_i         , jpi*jpj*jpl, (/1/) )     

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/NST_SRC/agrif_ice.F90

    r8586 r8637  
    1010   !!---------------------------------------------------------------------- 
    1111   !!   'key_agrif'                                              AGRIF zoom 
     12   !!   'key_lim3'                                       LIM3 sea-ice model 
    1213   !!----------------------------------------------------------------------    
    1314   IMPLICIT NONE 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/NST_SRC/agrif_oce.F90

    r8586 r8637  
    4848 
    4949   ! Barotropic arrays used to store open boundary data during time-splitting loop: 
    50    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_w, vbdy_w, hbdy_w 
    51    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_e, vbdy_e, hbdy_e 
    52    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_n, vbdy_n, hbdy_n 
    53    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_s, vbdy_s, hbdy_s 
     50   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_w, vbdy_w, hbdy_w 
     51   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_e, vbdy_e, hbdy_e 
     52   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_n, vbdy_n, hbdy_n 
     53   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::  ubdy_s, vbdy_s, hbdy_s 
    5454 
    55    INTEGER :: tsn_id                                                  ! AGRIF profile for tracers interpolation and update 
    56    INTEGER :: un_interp_id, vn_interp_id                              ! AGRIF profiles for interpolations 
    57    INTEGER :: un_update_id, vn_update_id                              ! AGRIF profiles for udpates 
    58    INTEGER :: tsn_sponge_id, un_sponge_id, vn_sponge_id               ! AGRIF profiles for sponge layers 
     55 
     56   INTEGER, PUBLIC :: tsn_id                                                  ! AGRIF profile for tracers interpolation and update 
     57   INTEGER, PUBLIC :: un_interp_id, vn_interp_id                              ! AGRIF profiles for interpolations 
     58   INTEGER, PUBLIC :: un_update_id, vn_update_id                              ! AGRIF profiles for udpates 
     59   INTEGER, PUBLIC :: tsn_sponge_id, un_sponge_id, vn_sponge_id               ! AGRIF profiles for sponge layers 
    5960# if defined key_top 
    60    INTEGER :: trn_id, trn_sponge_id 
     61   INTEGER, PUBLIC :: trn_id, trn_sponge_id 
    6162# endif   
    62    INTEGER :: unb_id, vnb_id, ub2b_interp_id, vb2b_interp_id 
    63    INTEGER :: ub2b_update_id, vb2b_update_id 
    64    INTEGER :: e3t_id, e1u_id, e2v_id, sshn_id 
    65    INTEGER :: scales_t_id 
    66 # if defined key_zdftke 
    67    INTEGER :: avt_id, avm_id, en_id 
    68 # endif   
    69    INTEGER :: umsk_id, vmsk_id 
    70    INTEGER :: kindic_agr 
    71  
     63   INTEGER, PUBLIC :: unb_id, vnb_id, ub2b_interp_id, vb2b_interp_id 
     64   INTEGER, PUBLIC :: ub2b_update_id, vb2b_update_id 
     65   INTEGER, PUBLIC :: e3t_id, e1u_id, e2v_id, sshn_id 
     66   INTEGER, PUBLIC :: scales_t_id 
     67   INTEGER, PUBLIC :: avt_id, avm_id, en_id                ! TKE related identificators 
     68   INTEGER, PUBLIC :: umsk_id, vmsk_id 
     69   INTEGER, PUBLIC :: kindic_agr 
     70    
    7271   !!---------------------------------------------------------------------- 
    73    !! NEMO/NST 3.3.1 , NEMO Consortium (2011) 
     72   !! NEMO/NST 4.0 , NEMO Consortium (2017) 
    7473   !! $Id$ 
    7574   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90

    r8586 r8637  
    594594   END SUBROUTINE Agrif_ssh_ts 
    595595 
    596 # if defined key_zdftke 
    597596 
    598597   SUBROUTINE Agrif_tke 
    … …  
    615614   END SUBROUTINE Agrif_tke 
    616615    
    617 # endif 
    618616 
    619617   SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice.F90

    r8586 r8637  
    1616   !!---------------------------------------------------------------------- 
    1717   USE oce             ! ocean dynamics and tracers variables 
    18    USE ice             ! LIM_3 ice variables 
    19    USE icevar 
    20    USE icectl 
     18   USE ice             ! sea-ice: variables 
     19   USE icevar          ! sea-ice: operations 
     20   USE iceitd          ! sea-ice: rebining 
     21   USE icectl          ! sea-ice: control prints 
    2122   USE phycst          ! physical constant 
    2223   USE eosbn2          ! equation of state 
    … …  
    3940 
    4041   !!---------------------------------------------------------------------- 
    41    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     42   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    4243   !! $Id: bdyice.F90 8306 2017-07-10 10:18:03Z clem $ 
    4344   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    249250      END DO !jl 
    250251      ! 
     252      ! --- In case categories are out of bounds, do a remapping --- ! 
     253      !     i.e. inputs have not the same ice thickness distribution  
     254      !          (set by rn_himean) than the regional simulation 
     255      IF( jpl > 1 )   CALL ice_itd_reb( kt ) 
    251256      !       
    252257      IF( nn_timing == 1 )   CALL timing_stop('bdy_ice_frs') 
    … …  
    273278      !!------------------------------------------------------------------------------ 
    274279      ! 
    275       IF( nn_timing == 1 ) CALL timing_start('bdy_ice_dyn') 
     280      IF( ln_timing )  CALL timing_start('bdy_ice_dyn') 
    276281      ! 
    277282      DO ib_bdy=1, nb_bdy 
    … …  
    350355      END DO 
    351356      ! 
    352       IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_dyn') 
     357      IF( ln_timing )  CALL timing_stop('bdy_ice_dyn') 
    353358      ! 
    354359    END SUBROUTINE bdy_ice_dyn 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90

    r8586 r8637  
    4040   USE lib_fortran    ! to use key_nosignedzero 
    4141#if defined key_lim3 
    42    USE ice     , ONLY :   u_ice, v_ice, jpl, a_i_b, at_i_b 
     42   USE ice     , ONLY :   u_ice, v_ice, jpl, a_i_b, at_i_b, tm_su 
    4343   USE icethd_dh      ! for CALL ice_thd_snwblow 
    4444#endif 
    … …  
    9292   REAL(wp), PARAMETER ::   Ls     =    2.839e6     ! latent heat of sublimation 
    9393   REAL(wp), PARAMETER ::   Stef   =    5.67e-8     ! Stefan Boltzmann constant 
    94    REAL(wp), PARAMETER ::   Cd_ice =    1.4e-3      ! iovi 1.63e-3     ! transfer coefficient over ice 
     94   REAL(wp), PARAMETER ::   Cd_ice =    1.4e-3      ! transfer coefficient over ice 
    9595   REAL(wp), PARAMETER ::   albo   =    0.066       ! ocean albedo assumed to be constant 
    9696   ! 
    … …  
    108108   REAL(wp) ::   rn_zu          ! z(u)   : height of wind measurements 
    109109!!gm ref namelist initialize it so remove the setting to false below 
    110    LOGICAL  ::   ln_Cd_L12 = .FALSE. !  Modify the drag ice-atm and oce-atm depending on ice concentration (from Lupkes et al. JGR2012) 
     110   LOGICAL  ::   ln_Cd_L12 = .FALSE. !  Modify the drag ice-atm depending on ice concentration (from Lupkes et al. JGR2012) 
     111   LOGICAL  ::   ln_Cd_L15 = .FALSE. !  Modify the drag ice-atm depending on ice concentration (from Lupkes et al. JGR2015) 
    111112   ! 
    112    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_oce   ! air-ocean drag (clem) 
     113!!cr   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_oce   ! air-ocean drag (clem) 
     114   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_atm                    ! transfer coefficient for momentum      (tau) 
     115   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Ch_atm                    ! transfer coefficient for sensible heat (Q_sens) 
     116   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Ce_atm                    ! tansfert coefficient for evaporation   (Q_lat) 
     117   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_zu                      ! air temperature at wind speed height (needed by Lupkes 2015 bulk scheme) 
     118   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   q_zu                      ! air spec. hum.  at wind speed height (needed by Lupkes 2015 bulk scheme) 
     119   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   cdn_oce, chn_oce, cen_oce ! needed by Lupkes 2015 bulk scheme 
    113120 
    114121   INTEGER  ::   nblk           ! choice of the bulk algorithm 
    … …  
    132139      !!             ***  ROUTINE sbc_blk_alloc *** 
    133140      !!------------------------------------------------------------------- 
    134       ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc ) 
     141!!cr      ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc ) 
     142      ALLOCATE( Cd_atm (jpi,jpj), Ch_atm (jpi,jpj), Ce_atm (jpi,jpj), t_zu(jpi,jpj), q_zu(jpi,jpj), & 
     143         &      cdn_oce(jpi,jpj), chn_oce(jpi,jpj), cen_oce(jpi,jpj), STAT=sbc_blk_alloc ) 
    135144      ! 
    136145      IF( lk_mpp             )   CALL mpp_sum ( sbc_blk_alloc ) 
    … …  
    164173         &                 ln_NCAR, ln_COARE_3p0, ln_COARE_3p5, ln_ECMWF,             &   ! bulk algorithm 
    165174         &                 cn_dir , ln_taudif, rn_zqt, rn_zu,                         &  
    166          &                 rn_pfac, rn_efac, rn_vfac, ln_Cd_L12 
     175         &                 rn_pfac, rn_efac, rn_vfac, ln_Cd_L12, ln_Cd_L15 
    167176      !!--------------------------------------------------------------------- 
    168177      ! 
    … …  
    255264         WRITE(numout,*) '      factor applied on ocean/ice velocity                rn_vfac      = ', rn_vfac 
    256265         WRITE(numout,*) '         (form absolute (=0) to relative winds(=1))' 
     266         WRITE(numout,*) '      use ice-atm drag from Lupkes2012                    ln_Cd_L12    = ', ln_Cd_L12 
     267         WRITE(numout,*) '      use ice-atm drag from Lupkes2015                    ln_Cd_L15    = ', ln_Cd_L15 
    257268         ! 
    258269         WRITE(numout,*) 
    … …  
    361372      REAL(wp), DIMENSION(jpi,jpj) ::   zqlw, zqsb        ! long wave and sensible heat fluxes 
    362373      REAL(wp), DIMENSION(jpi,jpj) ::   zqla, zevap       ! latent heat fluxes and evaporation 
    363       REAL(wp), DIMENSION(jpi,jpj) ::   zCd               ! transfer coefficient for momentum      (tau) 
    364       REAL(wp), DIMENSION(jpi,jpj) ::   zCh               ! transfer coefficient for sensible heat (Q_sens) 
    365       REAL(wp), DIMENSION(jpi,jpj) ::   zCe               ! tansfert coefficient for evaporation   (Q_lat) 
    366374      REAL(wp), DIMENSION(jpi,jpj) ::   zst               ! surface temperature in Kelvin 
    367       REAL(wp), DIMENSION(jpi,jpj) ::   zt_zu             ! air temperature at wind speed height 
    368       REAL(wp), DIMENSION(jpi,jpj) ::   zq_zu             ! air spec. hum.  at wind speed height 
    369375      REAL(wp), DIMENSION(jpi,jpj) ::   zU_zu             ! bulk wind speed at height zu  [m/s] 
    370376      REAL(wp), DIMENSION(jpi,jpj) ::   ztpot             ! potential temperature of air at z=rn_zqt [K] 
    … …  
    418424      zqlw(:,:) = (  sf(jp_qlw)%fnow(:,:,1) - Stef * zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:)  ) * tmask(:,:,1)   ! Long  Wave 
    419425 
    420  
    421  
    422426      ! ----------------------------------------------------------------------------- ! 
    423427      !     II    Turbulent FLUXES                                                    ! 
    … …  
    435439      ! 
    436440      CASE( np_NCAR      )   ;   CALL turb_ncar    ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! NCAR-COREv2 
    437          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     441         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    438442      CASE( np_COARE_3p0 )   ;   CALL turb_coare   ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! COARE v3.0 
    439          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     443         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    440444      CASE( np_COARE_3p5 )   ;   CALL turb_coare3p5( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! COARE v3.5 
    441          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     445         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    442446      CASE( np_ECMWF     )   ;   CALL turb_ecmwf   ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! ECMWF 
    443          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     447         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    444448      CASE DEFAULT 
    445449         CALL ctl_stop( 'STOP', 'sbc_oce: non-existing bulk formula selected' ) 
    … …  
    448452      !                          ! Compute true air density : 
    449453      IF( ABS(rn_zu - rn_zqt) > 0.01 ) THEN     ! At zu: (probably useless to remove zrho*grav*rn_zu from SLP...) 
    450          zrhoa(:,:) = rho_air( zt_zu(:,:)             , zq_zu(:,:)             , sf(jp_slp)%fnow(:,:,1) ) 
     454         zrhoa(:,:) = rho_air( t_zu(:,:)              , q_zu(:,:)              , sf(jp_slp)%fnow(:,:,1) ) 
    451455      ELSE                                      ! At zt: 
    452456         zrhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1) ) 
    453457      END IF 
    454458 
    455       Cd_oce(:,:) = zCd(:,:)     ! record value of pure ocean-atm. drag (clem) 
     459!!      CALL iom_put( "Cd_oce", Cd_atm)  ! output value of pure ocean-atm. transfer coef. 
     460!!      CALL iom_put( "Ch_oce", Ch_atm)  ! output value of pure ocean-atm. transfer coef. 
    456461 
    457462      DO jj = 1, jpj             ! tau module, i and j component 
    458463         DO ji = 1, jpi 
    459             zztmp = zrhoa(ji,jj)  * zU_zu(ji,jj) * zCd(ji,jj)   ! using bulk wind speed 
     464            zztmp = zrhoa(ji,jj)  * zU_zu(ji,jj) * Cd_atm(ji,jj)   ! using bulk wind speed 
    460465            taum  (ji,jj) = zztmp * wndm  (ji,jj) 
    461466            zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) 
    … …  
    492497      IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 
    493498         !! q_air and t_air are given at 10m (wind reference height) 
    494          zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*zCe(:,:)*(zsq(:,:) - sf(jp_humi)%fnow(:,:,1)) )  ! Evaporation, using bulk wind speed 
    495          zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*zCh(:,:)*(zst(:,:) - ztpot(:,:)               )   ! Sensible Heat, using bulk wind speed 
     499         zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - sf(jp_humi)%fnow(:,:,1)) ) ! Evaporation, using bulk wind speed 
     500         zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - ztpot(:,:)             )   ! Sensible Heat, using bulk wind speed 
    496501      ELSE 
    497502         !! q_air and t_air are not given at 10m (wind reference height) 
    498503         ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 
    499          zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*zCe(:,:)*(zsq(:,:) - zq_zu(:,:) ) ) ! Evaporation ! using bulk wind speed 
    500          zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*zCh(:,:)*(zst(:,:) - zt_zu(:,:) )   ! Sensible Heat ! using bulk wind speed 
     504         zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - q_zu(:,:) ) ) ! Evaporation, using bulk wind speed 
     505         zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - t_zu(:,:) )   ! Sensible Heat, using bulk wind speed 
    501506      ENDIF 
    502507 
    … …  
    505510 
    506511      IF(ln_ctl) THEN 
    507          CALL prt_ctl( tab2d_1=zqla  , clinfo1=' blk_oce: zqla   : ', tab2d_2=zCe , clinfo2=' Ce  : ' ) 
    508          CALL prt_ctl( tab2d_1=zqsb  , clinfo1=' blk_oce: zqsb   : ', tab2d_2=zCh , clinfo2=' Ch  : ' ) 
     512         CALL prt_ctl( tab2d_1=zqla  , clinfo1=' blk_oce: zqla   : ', tab2d_2=Ce_atm , clinfo2=' Ce_oce  : ' ) 
     513         CALL prt_ctl( tab2d_1=zqsb  , clinfo1=' blk_oce: zqsb   : ', tab2d_2=Ch_atm , clinfo2=' Ch_oce  : ' ) 
    509514         CALL prt_ctl( tab2d_1=zqlw  , clinfo1=' blk_oce: zqlw   : ', tab2d_2=qsr, clinfo2=' qsr : ' ) 
    510515         CALL prt_ctl( tab2d_1=zsq   , clinfo1=' blk_oce: zsq    : ', tab2d_2=zst, clinfo2=' zst : ' ) 
    … …  
    574579      !!--------------------------------------------------------------------- 
    575580      INTEGER  ::   ji, jj    ! dummy loop indices 
    576       REAL(wp) ::   zwndi_f , zwndj_f, zwnorm_f      ! relative wind module and components at F-point 
    577       REAL(wp) ::   zwndi_t , zwndj_t                ! relative wind components at T-point 
    578       REAL(wp), DIMENSION(jpi,jpj) ::   zCd, zrhoa   ! transfer coefficient for momentum      (tau) 
     581      REAL(wp) ::   zwndi_f , zwndj_f, zwnorm_f   ! relative wind module and components at F-point 
     582      REAL(wp) ::   zwndi_t , zwndj_t             ! relative wind components at T-point 
     583      REAL(wp), DIMENSION(jpi,jpj) ::   zrhoa     ! transfer coefficient for momentum      (tau) 
    579584      !!--------------------------------------------------------------------- 
    580585      ! 
    581586      IF( ln_timing )   CALL timing_start('blk_ice_tau') 
    582587      ! 
    583       IF( ln_Cd_L12 ) THEN 
    584          CALL Cdn10_Lupkes2012( zCd )  ! air-ice drag = F(ice concentration) (see Lupkes et al., 2012) 
    585       ELSE 
    586          zCd(:,:) = Cd_ice             ! constant air-ice drag 
    587       ENDIF 
    588  
    589       ! local scalars ( place there for vector optimisation purposes) 
    590       ! Computing density of air! Way denser that 1.2 over sea-ice !!! 
    591       !! 
    592       zrhoa (:,:) =  rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 
    593  
    594       !!gm brutal.... 
    595       utau_ice  (:,:) = 0._wp 
    596       vtau_ice  (:,:) = 0._wp 
    597       wndm_ice  (:,:) = 0._wp 
    598       !!gm end 
    599  
    600       ! ----------------------------------------------------------------------------- ! 
    601       !    Wind components and module relative to the moving ocean ( U10m - U_ice )   ! 
    602       ! ----------------------------------------------------------------------------- ! 
     588      ! set transfer coefficients to default sea-ice values 
     589      Cd_atm(:,:) = Cd_ice 
     590      Ch_atm(:,:) = Cd_ice 
     591      Ce_atm(:,:) = Cd_ice 
     592 
     593      wndm_ice(:,:) = 0._wp      !!gm brutal.... 
     594 
     595      ! ------------------------------------------------------------ ! 
     596      !    Wind module relative to the moving ice ( U10m - U_ice )   ! 
     597      ! ------------------------------------------------------------ ! 
    603598      SELECT CASE( cp_ice_msh ) 
    604599      CASE( 'I' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
    … …  
    606601         DO jj = 2, jpjm1 
    607602            DO ji = 2, jpim1   ! B grid : NO vector opt 
    608                ! ... scalar wind at I-point (fld being at T-point) 
    609                zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)   & 
    610                   &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) - rn_vfac * u_ice(ji,jj) 
    611                zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
    612                   &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
    613                zwnorm_f = zrhoa(ji,jj) * zCd(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
    614                ! ... ice stress at I-point 
    615                utau_ice(ji,jj) = zwnorm_f * zwndi_f 
    616                vtau_ice(ji,jj) = zwnorm_f * zwndj_f 
    617603               ! ... scalar wind at T-point (fld being at T-point) 
    618604               zwndi_t = sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.25 * (  u_ice(ji,jj+1) + u_ice(ji+1,jj+1)   & 
    … …  
    623609            END DO 
    624610         END DO 
    625          CALL lbc_lnk( utau_ice, 'I', -1. ) 
    626          CALL lbc_lnk( vtau_ice, 'I', -1. ) 
    627611         CALL lbc_lnk( wndm_ice, 'T',  1. ) 
    628612         ! 
    629613      CASE( 'C' )                  ! C-grid ice dynamics :   U & V-points (same as ocean) 
    630          DO jj = 2, jpj 
    631             DO ji = fs_2, jpi   ! vect. opt. 
     614         DO jj = 2, jpjm1 
     615            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    632616               zwndi_t = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( u_ice(ji-1,jj  ) + u_ice(ji,jj) )  ) 
    633617               zwndj_t = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( v_ice(ji  ,jj-1) + v_ice(ji,jj) )  ) 
    … …  
    635619            END DO 
    636620         END DO 
     621         CALL lbc_lnk( wndm_ice, 'T',  1. ) 
     622         ! 
     623      END SELECT 
     624 
     625      ! Make ice-atm. drag dependent on ice concentration 
     626      IF    ( ln_Cd_L12 ) THEN   ! calculate new drag from Lupkes(2012) equations 
     627         CALL Cdn10_Lupkes2012( Cd_atm ) 
     628         Ch_atm(:,:) = Cd_atm(:,:)       ! momentum and heat transfer coef. are considered identical 
     629      ELSEIF( ln_Cd_L15 ) THEN   ! calculate new drag from Lupkes(2015) equations 
     630         CALL Cdn10_Lupkes2015( Cd_atm, Ch_atm )  
     631      ENDIF 
     632 
     633!!      CALL iom_put( "Cd_ice", Cd_atm)  ! output value of pure ice-atm. transfer coef. 
     634!!      CALL iom_put( "Ch_ice", Ch_atm)  ! output value of pure ice-atm. transfer coef. 
     635 
     636      ! local scalars ( place there for vector optimisation purposes) 
     637      ! Computing density of air! Way denser that 1.2 over sea-ice !!! 
     638      zrhoa (:,:) =  rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 
     639 
     640      !!gm brutal.... 
     641      utau_ice  (:,:) = 0._wp 
     642      vtau_ice  (:,:) = 0._wp 
     643      !!gm end 
     644 
     645      ! ------------------------------------------------------------ ! 
     646      !    Wind stress relative to the moving ice ( U10m - U_ice )   ! 
     647      ! ------------------------------------------------------------ ! 
     648      SELECT CASE( cp_ice_msh ) 
     649      CASE( 'I' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
     650         DO jj = 2, jpjm1 
     651            DO ji = 2, jpim1   ! B grid : NO vector opt 
     652               ! ... scalar wind at I-point (fld being at T-point) 
     653               zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)   & 
     654                  &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) - rn_vfac * u_ice(ji,jj) 
     655               zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
     656                  &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
     657               ! ... ice stress at I-point 
     658               zwnorm_f = zrhoa(ji,jj) * Cd_atm(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
     659               utau_ice(ji,jj) = zwnorm_f * zwndi_f 
     660               vtau_ice(ji,jj) = zwnorm_f * zwndj_f 
     661            END DO 
     662         END DO 
     663         CALL lbc_lnk( utau_ice, 'I', -1. ) 
     664         CALL lbc_lnk( vtau_ice, 'I', -1. ) 
     665         ! 
     666      CASE( 'C' )                  ! C-grid ice dynamics :   U & V-points (same as ocean) 
    637667         DO jj = 2, jpjm1 
    638668            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    639                utau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * zCd(ji,jj) * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )                          & 
     669               utau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd_atm(ji,jj) * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )            & 
    640670                  &          * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 
    641                vtau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * zCd(ji,jj) * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )                          & 
     671               vtau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd_atm(ji,jj) * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )            & 
    642672                  &          * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 
    643673            END DO 
    … …  
    645675         CALL lbc_lnk( utau_ice, 'U', -1. ) 
    646676         CALL lbc_lnk( vtau_ice, 'V', -1. ) 
    647          CALL lbc_lnk( wndm_ice, 'T',  1. ) 
    648677         ! 
    649678      END SELECT 
    … …  
    684713      REAL(wp), DIMENSION(jpi,jpj)     ::   zevap, zsnw   ! evaporation and snw distribution after wind blowing (LIM3) 
    685714      REAL(wp), DIMENSION(jpi,jpj)     ::   zrhoa 
    686       REAL(wp), DIMENSION(jpi,jpj)     ::   zCd            ! transfer coefficient for momentum      (tau) 
    687715      !!--------------------------------------------------------------------- 
    688716      ! 
    689717      IF( ln_timing )   CALL timing_start('blk_ice_flx') 
    690718      ! 
    691       IF( ln_Cd_L12 ) THEN 
    692          CALL Cdn10_Lupkes2012( zCd )  ! air-ice drag = F(ice concentration) (see Lupkes et al., 2012) 
    693       ELSE 
    694          zCd(:,:) = Cd_ice             ! constant air-ice drag 
    695       ENDIF 
    696  
    697       ! 
    698       ! local scalars ( place there for vector optimisation purposes) 
     719      ! 
     720      ! local scalars 
    699721      zcoef_dqlw   = 4.0 * 0.95 * Stef 
    700722      zcoef_dqla   = -Ls * 11637800. * (-5897.8) 
    … …  
    724746               ! ----------------------------! 
    725747 
    726                ! ... turbulent heat fluxes 
     748               ! ... turbulent heat fluxes with Ch_atm recalculated in blk_ice_tau 
    727749               ! Sensible Heat 
    728                z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * zCd(ji,jj) * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 
     750               z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Ch_atm(ji,jj) * wndm_ice(ji,jj) * (ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1)) 
    729751               ! Latent Heat 
    730                qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls  * zCd(ji,jj) * wndm_ice(ji,jj)   & 
    731                   &                         * (  11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 
     752               qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls  * Ch_atm(ji,jj) * wndm_ice(ji,jj) *  & 
     753                  &                ( 11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 
    732754               ! Latent heat sensitivity for ice (Dqla/Dt) 
    733755               IF( qla_ice(ji,jj,jl) > 0._wp ) THEN 
    734                   dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * zCd(ji,jj) * wndm_ice(ji,jj) / ( zst2 ) * EXP( -5897.8 / ptsu(ji,jj,jl) ) 
     756                  dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * Ch_atm(ji,jj) * wndm_ice(ji,jj) / zst2 * EXP(-5897.8 / ptsu(ji,jj,jl)) 
    735757               ELSE 
    736758                  dqla_ice(ji,jj,jl) = 0._wp 
    … …  
    738760 
    739761               ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 
    740                z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * zCd(ji,jj) * wndm_ice(ji,jj) 
     762               z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Ch_atm(ji,jj) * wndm_ice(ji,jj) 
    741763 
    742764               ! ----------------------------! 
    … …  
    935957#if defined key_lim3 
    936958 
    937    SUBROUTINE Cdn10_Lupkes2012( pCd ) 
     959   SUBROUTINE Cdn10_Lupkes2012( Cd ) 
    938960      !!---------------------------------------------------------------------- 
    939961      !!                      ***  ROUTINE  Cdn10_Lupkes2012  *** 
    … …  
    965987      !! 
    966988      !!---------------------------------------------------------------------- 
    967       REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pCd   ! air-ice drag coefficient 
     989      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Cd 
    968990      REAL(wp), PARAMETER ::   zCe   = 2.23e-03_wp 
    969991      REAL(wp), PARAMETER ::   znu   = 1._wp 
    … …  
    973995      !!---------------------------------------------------------------------- 
    974996      zcoef = znu + 1._wp / ( 10._wp * zbeta ) 
    975       ! 
     997 
    976998      ! generic drag over a cell partly covered by ice 
    977       !!pCd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) +  &                        ! pure ocean drag 
    978       !!   &       Cd_ice      *           at_i_b(:,:)   +  &                        ! pure ice drag 
    979       !!   &       zCe         * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu   ! change due to sea-ice morphology 
     999      !!Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) +  &                        ! pure ocean drag 
     1000      !!   &      Cd_ice      *           at_i_b(:,:)   +  &                        ! pure ice drag 
     1001      !!   &      zCe         * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu   ! change due to sea-ice morphology 
    9801002 
    9811003      ! ice-atm drag 
    982       pCd(:,:) = Cd_ice +  &                                                         ! pure ice drag 
    983          &       zCe    * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp)  ! change due to sea-ice morphology 
    984       ! 
     1004      Cd(:,:) = Cd_ice +  &                                                         ! pure ice drag 
     1005         &      zCe    * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp)  ! change due to sea-ice morphology 
     1006       
    9851007   END SUBROUTINE Cdn10_Lupkes2012 
     1008 
     1009 
     1010   SUBROUTINE Cdn10_Lupkes2015( Cd, Ch ) 
     1011      !!---------------------------------------------------------------------- 
     1012      !!                      ***  ROUTINE  Cdn10_Lupkes2015  *** 
     1013      !! 
     1014      !! ** pUrpose :    1lternative turbulent transfert coefficients formulation 
     1015      !!                 between sea-ice and atmosphere with distinct momentum  
     1016      !!                 and heat coefficients depending on sea-ice concentration  
     1017      !!                 and atmospheric stability (no meltponds effect for now). 
     1018      !!                 
     1019      !! ** Method :     The parameterization is adapted from Lupkes et al. (2015) 
     1020      !!                 and ECHAM6 atmospheric model. Compared to Lupkes2012 scheme, 
     1021      !!                 it considers specific skin and form drags (Andreas et al. 2010) 
     1022      !!                 to compute neutral transfert coefficients for both heat and  
     1023      !!                 momemtum fluxes. Atmospheric stability effect on transfert 
     1024      !!                 coefficient is also taken into account following Louis (1979). 
     1025      !! 
     1026      !! ** References : Lupkes et al. JGR 2015 (theory) 
     1027      !!                 Lupkes et al. ECHAM6 documentation 2015 (implementation) 
     1028      !! 
     1029      !!---------------------------------------------------------------------- 
     1030      ! 
     1031      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Cd 
     1032      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Ch 
     1033      REAL(wp), DIMENSION(jpi,jpj)            ::   zst, zqo_sat, zqi_sat 
     1034      ! 
     1035      ! ECHAM6 constants 
     1036      REAL(wp), PARAMETER ::   z0_skin_ice  = 0.69e-3_wp  ! Eq. 43 [m] 
     1037      REAL(wp), PARAMETER ::   z0_form_ice  = 0.57e-3_wp  ! Eq. 42 [m] 
     1038      REAL(wp), PARAMETER ::   z0_ice       = 1.00e-3_wp  ! Eq. 15 [m] 
     1039      REAL(wp), PARAMETER ::   zce10        = 2.80e-3_wp  ! Eq. 41 
     1040      REAL(wp), PARAMETER ::   zbeta        = 1.1_wp      ! Eq. 41 
     1041      REAL(wp), PARAMETER ::   zc           = 5._wp       ! Eq. 13 
     1042      REAL(wp), PARAMETER ::   zc2          = zc * zc 
     1043      REAL(wp), PARAMETER ::   zam          = 2. * zc     ! Eq. 14 
     1044      REAL(wp), PARAMETER ::   zah          = 3. * zc     ! Eq. 30 
     1045      REAL(wp), PARAMETER ::   z1_alpha     = 1._wp / 0.2_wp  ! Eq. 51 
     1046      REAL(wp), PARAMETER ::   z1_alphaf    = z1_alpha    ! Eq. 56 
     1047      REAL(wp), PARAMETER ::   zbetah       = 1.e-3_wp    ! Eq. 26 
     1048      REAL(wp), PARAMETER ::   zgamma       = 1.25_wp     ! Eq. 26 
     1049      REAL(wp), PARAMETER ::   z1_gamma     = 1._wp / zgamma 
     1050      REAL(wp), PARAMETER ::   r1_3         = 1._wp / 3._wp 
     1051      ! 
     1052      INTEGER  ::   ji, jj         ! dummy loop indices 
     1053      REAL(wp) ::   zthetav_os, zthetav_is, zthetav_zu 
     1054      REAL(wp) ::   zrib_o, zrib_i 
     1055      REAL(wp) ::   zCdn_skin_ice, zCdn_form_ice, zCdn_ice 
     1056      REAL(wp) ::   zChn_skin_ice, zChn_form_ice 
     1057      REAL(wp) ::   z0w, z0i, zfmi, zfmw, zfhi, zfhw 
     1058      REAL(wp) ::   zCdn_form_tmp 
     1059      !!---------------------------------------------------------------------- 
     1060 
     1061      ! Momentum Neutral Transfert Coefficients (should be a constant) 
     1062      zCdn_form_tmp = zce10 * ( LOG( 10._wp / z0_form_ice + 1._wp ) / LOG( rn_zu / z0_form_ice + 1._wp ) )**2   ! Eq. 40 
     1063      zCdn_skin_ice = ( vkarmn                                      / LOG( rn_zu / z0_skin_ice + 1._wp ) )**2   ! Eq. 7 
     1064      zCdn_ice      = zCdn_skin_ice   ! Eq. 7 (cf Lupkes email for details) 
     1065      !zCdn_ice     = 1.89e-3         ! old ECHAM5 value (cf Eq. 32) 
     1066 
     1067      ! Heat Neutral Transfert Coefficients 
     1068      zChn_skin_ice = vkarmn**2 / ( LOG( rn_zu / z0_ice + 1._wp ) * LOG( rn_zu * z1_alpha / z0_skin_ice + 1._wp ) )   ! Eq. 50 + Eq. 52 (cf Lupkes email for details) 
     1069      
     1070      ! Atmospheric and Surface Variables 
     1071      zst(:,:)     = sst_m(:,:) + rt0                                       ! convert SST from Celcius to Kelvin 
     1072      zqo_sat(:,:) = 0.98_wp * q_sat( zst(:,:)  , sf(jp_slp)%fnow(:,:,1) )  ! saturation humidity over ocean [kg/kg] 
     1073      zqi_sat(:,:) = 0.98_wp * q_sat( tm_su(:,:), sf(jp_slp)%fnow(:,:,1) )  ! saturation humidity over ice   [kg/kg] 
     1074      ! 
     1075      DO jj = 2, jpjm1           ! reduced loop is necessary for reproducibility 
     1076         DO ji = fs_2, fs_jpim1 
     1077            ! Virtual potential temperature [K] 
     1078            zthetav_os = zst(ji,jj)   * ( 1._wp + rctv0 * zqo_sat(ji,jj) )   ! over ocean 
     1079            zthetav_is = tm_su(ji,jj) * ( 1._wp + rctv0 * zqi_sat(ji,jj) )   ! ocean ice 
     1080            zthetav_zu = t_zu (ji,jj) * ( 1._wp + rctv0 * q_zu(ji,jj)    )   ! at zu 
     1081             
     1082            ! Bulk Richardson Number (could use Ri_bulk function from aerobulk instead) 
     1083            zrib_o = grav / zthetav_os * ( zthetav_zu - zthetav_os ) * rn_zu / MAX( 0.5, wndm(ji,jj)     )**2   ! over ocean 
     1084            zrib_i = grav / zthetav_is * ( zthetav_zu - zthetav_is ) * rn_zu / MAX( 0.5, wndm_ice(ji,jj) )**2   ! over ice 
     1085             
     1086            ! Momentum and Heat Neutral Transfert Coefficients 
     1087            zCdn_form_ice = zCdn_form_tmp * at_i_b(ji,jj) * ( 1._wp - at_i_b(ji,jj) )**zbeta  ! Eq. 40 
     1088            zChn_form_ice = zCdn_form_ice / ( 1._wp + ( LOG( z1_alphaf ) / vkarmn ) * SQRT( zCdn_form_ice ) )               ! Eq. 53  
     1089                        
     1090            ! Momentum and Heat Stability functions (possibility to use psi_m_ecmwf instead) 
     1091            z0w = rn_zu * EXP( -1._wp * vkarmn / SQRT( Cdn_oce(ji,jj) ) ) ! over water 
     1092            z0i = z0_skin_ice                                             ! over ice (cf Lupkes email for details) 
     1093            IF( zrib_o <= 0._wp ) THEN 
     1094               zfmw = 1._wp - zam * zrib_o / ( 1._wp + 3._wp * zc2 * Cdn_oce(ji,jj) * SQRT( -zrib_o * ( rn_zu / z0w + 1._wp ) ) )  ! Eq. 10 
     1095               zfhw = ( 1._wp + ( zbetah * ( zthetav_os - zthetav_zu )**r1_3 / ( Chn_oce(ji,jj) * MAX(0.01, wndm(ji,jj)) )   &     ! Eq. 26 
     1096                  &             )**zgamma )**z1_gamma 
     1097            ELSE 
     1098               zfmw = 1._wp / ( 1._wp + zam * zrib_o / SQRT( 1._wp + zrib_o ) )   ! Eq. 12 
     1099               zfhw = 1._wp / ( 1._wp + zah * zrib_o / SQRT( 1._wp + zrib_o ) )   ! Eq. 28 
     1100            ENDIF 
     1101             
     1102            IF( zrib_i <= 0._wp ) THEN 
     1103               zfmi = 1._wp - zam * zrib_i / (1._wp + 3._wp * zc2 * zCdn_ice * SQRT( -zrib_i * ( rn_zu / z0i + 1._wp)))   ! Eq.  9 
     1104               zfhi = 1._wp - zah * zrib_i / (1._wp + 3._wp * zc2 * zCdn_ice * SQRT( -zrib_i * ( rn_zu / z0i + 1._wp)))   ! Eq. 25 
     1105            ELSE 
     1106               zfmi = 1._wp / ( 1._wp + zam * zrib_i / SQRT( 1._wp + zrib_i ) )   ! Eq. 11 
     1107               zfhi = 1._wp / ( 1._wp + zah * zrib_i / SQRT( 1._wp + zrib_i ) )   ! Eq. 27 
     1108            ENDIF 
     1109             
     1110            ! Momentum Transfert Coefficients (Eq. 38) 
     1111            Cd(ji,jj) = zCdn_skin_ice *   zfmi +  & 
     1112               &        zCdn_form_ice * ( zfmi * at_i_b(ji,jj) + zfmw * ( 1._wp - at_i_b(ji,jj) ) ) / MAX( 1.e-06, at_i_b(ji,jj) ) 
     1113             
     1114            ! Heat Transfert Coefficients (Eq. 49) 
     1115            Ch(ji,jj) = zChn_skin_ice *   zfhi +  & 
     1116               &        zChn_form_ice * ( zfhi * at_i_b(ji,jj) + zfhw * ( 1._wp - at_i_b(ji,jj) ) ) / MAX( 1.e-06, at_i_b(ji,jj) ) 
     1117            ! 
     1118         END DO 
     1119      END DO 
     1120      CALL lbc_lnk_multi( Cd, 'T',  1., Ch, 'T', 1. ) 
     1121      ! 
     1122   END SUBROUTINE Cdn10_Lupkes2015 
     1123 
    9861124#endif 
    987     
    9881125 
    9891126   !!====================================================================== 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_algo_coare.F90

    r7646 r8637  
    4444   USE lib_fortran    ! to use key_nosignedzero 
    4545 
    46  
    4746   IMPLICIT NONE 
    4847   PRIVATE 
    … …  
    5049   PUBLIC ::   TURB_COARE   ! called by sbcblk.F90 
    5150 
    52    !! COARE own values for given constants: 
    53    REAL(wp), PARAMETER :: & 
    54       &   zi0     = 600.,  &   !: scale height of the atmospheric boundary layer...1 
    55       &  Beta0    = 1.25, &   !: gustiness parameter 
    56       &  rctv0    = 0.608     !: constant to obtain virtual temperature... 
    57  
    58  
     51   !                                           !! COARE own values for given constants: 
     52   REAL(wp), PARAMETER ::   zi0     = 600.      ! scale height of the atmospheric boundary layer...1 
     53   REAL(wp), PARAMETER ::   Beta0   =   1.25    ! gustiness parameter 
     54   REAL(wp), PARAMETER ::   rctv0   =   0.608   ! constant to obtain virtual temperature... 
     55 
     56   !!---------------------------------------------------------------------- 
    5957CONTAINS 
    6058 
    6159   SUBROUTINE turb_coare( zt, zu, sst, t_zt, ssq, q_zt, U_zu, & 
    62       &                   Cd, Ch, Ce, t_zu, q_zu, U_blk ) 
     60      &                   Cd, Ch, Ce, t_zu, q_zu, U_blk,      & 
     61      &                   Cdn, Chn, Cen                       ) 
    6362      !!---------------------------------------------------------------------- 
    6463      !!                      ***  ROUTINE  turb_coare  *** 
    … …  
    106105      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu           [kg/kg] 
    107106      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     107      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    108108      ! 
    109109      INTEGER :: j_itt 
    … …  
    246246      Ce   = ztmp0*q_star/dq_zu 
    247247      ! 
     248      ztmp1 = zu + z0 
     249      Cdn = vkarmn*vkarmn / (log(ztmp1/z0 )*log(ztmp1/z0 )) 
     250      Chn = vkarmn*vkarmn / (log(ztmp1/z0t)*log(ztmp1/z0t)) 
     251      Cen = Chn 
     252      ! 
    248253      CALL wrk_dealloc( jpi,jpj, u_star, t_star, q_star, zeta_u, dt_zu, dq_zu ) 
    249254      CALL wrk_dealloc( jpi,jpj, znu_a, z0, z0t, ztmp0, ztmp1, ztmp2 ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_algo_coare3p5.F90

    r7646 r8637  
    5858CONTAINS 
    5959 
    60    SUBROUTINE turb_coare3p5( zt, zu, sst, t_zt, ssq, q_zt, U_zu, & 
    61       &                   Cd, Ch, Ce, t_zu, q_zu, U_blk ) 
     60   SUBROUTINE turb_coare3p5( zt, zu, sst, t_zt, ssq, q_zt, U_zu,  & 
     61      &                      Cd, Ch, Ce, t_zu, q_zu, U_blk,       & 
     62      &                      Cdn, Chn, Cen                        ) 
    6263      !!---------------------------------------------------------------------------------- 
    6364      !!                      ***  ROUTINE  turb_coare3p5  *** 
    … …  
    105106      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu             [kg/kg] 
    106107      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     108      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    107109      ! 
    108110      INTEGER :: j_itt 
    … …  
    252254      Ce   = ztmp0*q_star/dq_zu 
    253255      ! 
     256      ztmp1 = zu + z0 
     257      Cdn = vkarmn*vkarmn / (log(ztmp1/z0 )*log(ztmp1/z0 )) 
     258      Chn = vkarmn*vkarmn / (log(ztmp1/z0t)*log(ztmp1/z0t)) 
     259      Cen = Chn 
     260      ! 
    254261      CALL wrk_dealloc( jpi,jpj, u_star, t_star, q_star, zeta_u, dt_zu, dq_zu ) 
    255262      CALL wrk_dealloc( jpi,jpj, znu_a, z0, z0t, ztmp0, ztmp1, ztmp2 ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_algo_ecmwf.F90

    r7646 r8637  
    6464 
    6565   SUBROUTINE TURB_ECMWF( zt, zu, sst, t_zt, ssq , q_zt , U_zu,   & 
    66       &                   Cd, Ch, Ce , t_zu, q_zu, U_blk ) 
     66      &                   Cd, Ch, Ce , t_zu, q_zu, U_blk,         & 
     67      &                   Cdn, Chn, Cen                           ) 
    6768      !!---------------------------------------------------------------------------------- 
    6869      !!                      ***  ROUTINE  turb_ecmwf  *** 
    … …  
    112113      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu           [kg/kg] 
    113114      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     115      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    114116      ! 
    115117      INTEGER :: j_itt 
    … …  
    266268            dt_zu = t_zu - sst ;  dt_zu = SIGN( MAX(ABS(dt_zu),1.E-6), dt_zu ) 
    267269            dq_zu = q_zu - ssq ;  dq_zu = SIGN( MAX(ABS(dq_zu),1.E-9), dq_zu ) 
     270 
    268271         END IF 
    269272 
    … …  
    271274         ztmp1 = zu + z0 
    272275         ztmp0 = ztmp1*Linv 
    273          func_m = log(ztmp1) - LOG(z0) - psi_m_ecmwf(ztmp0) + psi_m_ecmwf(z0*Linv) 
     276         func_m = log(ztmp1) - LOG(z0 ) - psi_m_ecmwf(ztmp0) + psi_m_ecmwf(z0 *Linv) 
    274277         func_h = log(ztmp1) - LOG(z0t) - psi_h_ecmwf(ztmp0) + psi_h_ecmwf(z0t*Linv) 
    275278 
    … …  
    280283      ztmp1 = log((zu + z0)/z0q) - psi_h_ecmwf((zu + z0)*Linv) + psi_h_ecmwf(z0q*Linv)   ! func_q 
    281284      Ce = vkarmn*vkarmn/(func_m*ztmp1) 
     285 
     286      ztmp1 = zu + z0 
     287      Cdn = vkarmn*vkarmn / (log(ztmp1/z0 )*log(ztmp1/z0 )) 
     288      Chn = vkarmn*vkarmn / (log(ztmp1/z0t)*log(ztmp1/z0t)) 
     289      Cen = vkarmn*vkarmn / (log(ztmp1/z0q)*log(ztmp1/z0q)) 
    282290 
    283291      CALL wrk_dealloc( jpi,jpj,   u_star, t_star, q_star, func_m, func_h, dt_zu, dq_zu, Linv ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_algo_ncar.F90

    r7753 r8637  
    5353 
    5454   SUBROUTINE turb_ncar( zt, zu, sst, t_zt, ssq, q_zt, U_zu, & 
    55       &                  Cd, Ch, Ce, t_zu, q_zu, U_blk ) 
     55      &                  Cd, Ch, Ce, t_zu, q_zu, U_blk,      & 
     56      &                  Cdn, Chn, Cen                       ) 
    5657      !!---------------------------------------------------------------------------------- 
    5758      !!                      ***  ROUTINE  turb_ncar  *** 
    … …  
    112113      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu           [kg/kg] 
    113114      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     115      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    114116      ! 
    115117      INTEGER ::   j_itt 
    … …  
    199201            ztmp0 = MAX( 0.25 , U_blk/(1. + sqrt_Cd_n10/vkarmn*(LOG(zu/10.) - ztmp2)) ) ! U_n10 (ztmp2 == psi_m(zeta_u)) 
    200202            ztmp0 = cd_neutral_10m(ztmp0)                                               ! Cd_n10 
     203            Cdn(:,:) = ztmp0 
    201204            sqrt_Cd_n10 = sqrt(ztmp0) 
    202205 
    203206            stab    = 0.5 + sign(0.5,zeta_u)                           ! update stability 
    204207            Cx_n10  = 1.e-3*sqrt_Cd_n10*(18.*stab + 32.7*(1. - stab))  ! L&Y 2004 eq. (6c-6d)    (Cx_n10 == Ch_n10) 
     208            Chn(:,:) = Cx_n10 
    205209 
    206210            !! Update of transfer coefficients: 
    … …  
    216220 
    217221         Cx_n10  = 1.e-3 * (34.6 * sqrt_Cd_n10)  ! L&Y 2004 eq. (6b)    ! Cx_n10 == Ce_n10 
     222         Cen(:,:) = Cx_n10 
    218223         ztmp1 = 1. + Cx_n10*ztmp0 
    219224         Ce  = Cx_n10*ztmp2 / ztmp1  ! L&Y 2004 eq. (10c) 
    220  
     225         ! 
    221226      END DO 
    222  
     227      ! 
    223228      CALL wrk_dealloc( jpi,jpj,   Cx_n10, sqrt_Cd_n10, zeta_u, stab ) 
    224229      CALL wrk_dealloc( jpi,jpj,   zpsi_h_u, ztmp0, ztmp1, ztmp2 ) 
    225  
     230      ! 
    226231      IF( nn_timing == 1 )   CALL timing_stop('turb_ncar') 
    227  
     232      ! 
    228233   END SUBROUTINE turb_ncar 
    229234 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r8586 r8637  
    3333   USE geo2ocean      !  
    3434   USE oce   , ONLY : tsn, un, vn, sshn, ub, vb, sshb, fraqsr_1lev 
    35    USE albedooce      !  
     35   USE ocealb         !  
    3636   USE eosbn2         !  
    3737   USE sbcrnf, ONLY : l_rnfcpl 
    … …  
    173173                                         !   -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) 
    174174   TYPE ::   DYNARR      
    175       REAL(wp), POINTER, DIMENSION(:,:,:)    ::   z3    
     175      REAL(wp), POINTER, DIMENSION(:,:,:) ::   z3    
    176176   END TYPE DYNARR 
    177177 
    178178   TYPE( DYNARR ), SAVE, DIMENSION(jprcv) ::   frcv                     ! all fields recieved from the atmosphere 
    179179 
    180    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   albedo_oce_mix    ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
     180   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   alb_oce_mix    ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
    181181 
    182182   REAL(wp) ::   rpref = 101000._wp   ! reference atmospheric pressure[N/m2]  
    … …  
    202202      ierr(:) = 0 
    203203      ! 
    204       ALLOCATE( albedo_oce_mix(jpi,jpj), nrcvinfo(jprcv),  STAT=ierr(1) ) 
     204      ALLOCATE( alb_oce_mix(jpi,jpj), nrcvinfo(jprcv),  STAT=ierr(1) ) 
    205205       
    206206#if ! defined key_lim3 && ! defined key_cice 
    … …  
    736736      !     2. receiving mixed oce-ice solar radiation  
    737737      IF ( TRIM ( sn_snd_alb%cldes ) == 'mixed oce-ice' .OR. TRIM ( sn_rcv_qsr%cldes ) == 'mixed oce-ice' ) THEN 
    738          CALL albedo_oce( zaos, zacs ) 
     738         CALL oce_alb( zaos, zacs ) 
    739739         ! Due to lack of information on nebulosity : mean clear/overcast sky 
    740          albedo_oce_mix(:,:) = ( zacs(:,:) + zaos(:,:) ) * 0.5 
     740         alb_oce_mix(:,:) = ( zacs(:,:) + zaos(:,:) ) * 0.5 
    741741      ENDIF 
    742742 
    … …  
    18851885!       ( see OASIS3 user guide, 5th edition, p39 ) 
    18861886         zqsr_ice(:,:,1) = frcv(jpr_qsrmix)%z3(:,:,1) * ( 1.- palbi(:,:,1) )   & 
    1887             &            / (  1.- ( albedo_oce_mix(:,:  ) * ziceld(:,:)       & 
    1888             &                     + palbi         (:,:,1) * picefr(:,:) ) ) 
     1887            &            / (  1.- ( alb_oce_mix(:,:  ) * ziceld(:,:)       & 
     1888            &                     + palbi      (:,:,1) * picefr(:,:) ) ) 
    18891889      END SELECT 
    18901890      IF( ln_dm2dc .AND. ln_cpl ) THEN   ! modify qsr to include the diurnal cycle 
    … …  
    20522052                   ztmp1(:,:) = SUM( alb_ice (:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) / SUM( a_i(:,:,1:jpl), dim=3 ) 
    20532053                ELSEWHERE 
    2054                    ztmp1(:,:) = albedo_oce_mix(:,:) 
     2054                   ztmp1(:,:) = alb_oce_mix(:,:) 
    20552055                END WHERE 
    20562056             CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%clcat' ) 
    … …  
    20802080 
    20812081      IF( ssnd(jps_albmix)%laction ) THEN                         ! mixed ice-ocean 
    2082          ztmp1(:,:) = albedo_oce_mix(:,:) * zfr_l(:,:) 
     2082         ztmp1(:,:) = alb_oce_mix(:,:) * zfr_l(:,:) 
    20832083         DO jl=1,jpl 
    20842084            ztmp1(:,:) = ztmp1(:,:) + alb_ice(:,:,jl) * a_i(:,:,jl) 
    2085          ENDDO 
     2085         END DO 
    20862086         CALL cpl_snd( jps_albmix, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info ) 
    20872087      ENDIF 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssm.F90

    r8586 r8637  
    236236      ENDIF 
    237237      ! 
    238       IF( .NOT. l_ssm_mean ) THEN   ! default initialisation. needed by lim_istate 
     238      IF( .NOT.l_ssm_mean ) THEN   ! default initialisation. needed by iceistate 
    239239         ! 
    240240         IF(lwp) WRITE(numout,*) '   default initialisation of ss._m arrays' 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfiwm.F90

    r8568 r8637  
    3333   PRIVATE 
    3434 
    35    PUBLIC   zdf_iwm         ! called in step module  
    36    PUBLIC   zdf_iwm_init    ! called in nemogcm module  
    37  
    38    !                       !!* Namelist  namzdf_iwm : internal wave-driven mixing * 
    39    INTEGER  ::  nn_zpyc     ! pycnocline-intensified mixing energy proportional to N (=1) or N^2 (=2) 
    40    LOGICAL  ::  ln_mevar    ! variable (=T) or constant (=F) mixing efficiency 
    41    LOGICAL  ::  ln_tsdiff   ! account for differential T/S wave-driven mixing (=T) or not (=F) 
    42  
    43    REAL(wp) ::  r1_6 = 1._wp / 6._wp 
    44  
    45    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ebot_iwm     ! power available from high-mode wave breaking (W/m2) 
    46    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   epyc_iwm     ! power available from low-mode, pycnocline-intensified wave breaking (W/m2) 
    47    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ecri_iwm     ! power available from low-mode, critical slope wave breaking (W/m2) 
    48    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hbot_iwm     ! WKB decay scale for high-mode energy dissipation (m) 
    49    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hcri_iwm     ! decay scale for low-mode critical slope dissipation (m) 
    50    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   emix_iwm     ! local energy density available for mixing (W/kg) 
    51    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   bflx_iwm     ! buoyancy flux Kz * N^2 (W/kg) 
    52    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   pcmap_iwm    ! vertically integrated buoyancy flux (W/m2) 
    53    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   zav_ratio    ! S/T diffusivity ratio (only for ln_tsdiff=T) 
    54    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   zav_wave     ! Internal wave-induced diffusivity 
     35   PUBLIC   zdf_iwm        ! called in step module  
     36   PUBLIC   zdf_iwm_init   ! called in nemogcm module  
     37 
     38   !                      !!* Namelist  namzdf_iwm : internal wave-driven mixing * 
     39   INTEGER ::  nn_zpyc     ! pycnocline-intensified mixing energy proportional to N (=1) or N^2 (=2) 
     40   LOGICAL ::  ln_mevar    ! variable (=T) or constant (=F) mixing efficiency 
     41   LOGICAL ::  ln_tsdiff   ! account for differential T/S wave-driven mixing (=T) or not (=F) 
     42 
     43   REAL(wp)::  r1_6 = 1._wp / 6._wp 
     44 
     45   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ebot_iwm   ! power available from high-mode wave breaking (W/m2) 
     46   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   epyc_iwm   ! power available from low-mode, pycnocline-intensified wave breaking (W/m2) 
     47   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ecri_iwm   ! power available from low-mode, critical slope wave breaking (W/m2) 
     48   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hbot_iwm   ! WKB decay scale for high-mode energy dissipation (m) 
     49   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hcri_iwm   ! decay scale for low-mode critical slope dissipation (m) 
    5550 
    5651   !! * Substitutions 
    … …  
    6762      !!                ***  FUNCTION zdf_iwm_alloc  *** 
    6863      !!---------------------------------------------------------------------- 
    69       ALLOCATE(     ebot_iwm(jpi,jpj),  epyc_iwm(jpi,jpj),  ecri_iwm(jpi,jpj)    ,   & 
    70       &             hbot_iwm(jpi,jpj),  hcri_iwm(jpi,jpj),  emix_iwm(jpi,jpj,jpk),   & 
    71       &         bflx_iwm(jpi,jpj,jpk), pcmap_iwm(jpi,jpj), zav_ratio(jpi,jpj,jpk),   &  
    72       &         zav_wave(jpi,jpj,jpk), STAT=zdf_iwm_alloc     ) 
     64      ALLOCATE( ebot_iwm(jpi,jpj),  epyc_iwm(jpi,jpj),  ecri_iwm(jpi,jpj) ,     & 
     65      &         hbot_iwm(jpi,jpj),  hcri_iwm(jpi,jpj)                     , STAT=zdf_iwm_alloc ) 
    7366      ! 
    7467      IF( lk_mpp             )   CALL mpp_sum ( zdf_iwm_alloc ) 
    … …  
    8578      !! 
    8679      !! ** Method  : - internal wave-driven vertical mixing is given by: 
    87       !!                  Kz_wave = min(  100 cm2/s, f(  Reb = emix_iwm /( Nu * N^2 )  ) 
    88       !!              where emix_iwm is the 3D space distribution of the wave-breaking  
     80      !!                  Kz_wave = min(  100 cm2/s, f(  Reb = zemx_iwm /( Nu * N^2 )  ) 
     81      !!              where zemx_iwm is the 3D space distribution of the wave-breaking  
    8982      !!              energy and Nu the molecular kinematic viscosity. 
    9083      !!              The function f(Reb) is linear (constant mixing efficiency) 
    9184      !!              if the namelist parameter ln_mevar = F and nonlinear if ln_mevar = T. 
    9285      !! 
    93       !!              - Compute emix_iwm, the 3D power density that allows to compute 
     86      !!              - Compute zemx_iwm, the 3D power density that allows to compute 
    9487      !!              Reb and therefrom the wave-induced vertical diffusivity. 
    9588      !!              This is divided into three components: 
    9689      !!                 1. Bottom-intensified low-mode dissipation at critical slopes 
    97       !!                     emix_iwm(z) = ( ecri_iwm / rau0 ) * EXP( -(H-z)/hcri_iwm ) 
     90      !!                     zemx_iwm(z) = ( ecri_iwm / rau0 ) * EXP( -(H-z)/hcri_iwm ) 
    9891      !!                                   / ( 1. - EXP( - H/hcri_iwm ) ) * hcri_iwm 
    9992      !!              where hcri_iwm is the characteristic length scale of the bottom  
    10093      !!              intensification, ecri_iwm a map of available power, and H the ocean depth. 
    10194      !!                 2. Pycnocline-intensified low-mode dissipation 
    102       !!                     emix_iwm(z) = ( epyc_iwm / rau0 ) * ( sqrt(rn2(z))^nn_zpyc ) 
     95      !!                     zemx_iwm(z) = ( epyc_iwm / rau0 ) * ( sqrt(rn2(z))^nn_zpyc ) 
    10396      !!                                   / SUM( sqrt(rn2(z))^nn_zpyc * e3w(z) ) 
    10497      !!              where epyc_iwm is a map of available power, and nn_zpyc 
    … …  
    10699      !!              energy dissipation. 
    107100      !!                 3. WKB-height dependent high mode dissipation 
    108       !!                     emix_iwm(z) = ( ebot_iwm / rau0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 
     101      !!                     zemx_iwm(z) = ( ebot_iwm / rau0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 
    109102      !!                                   / SUM( rn2(z) * EXP(-z_wkb(z)/hbot_iwm) * e3w(z) ) 
    110103      !!              where hbot_iwm is the characteristic length scale of the WKB bottom  
    … …  
    121114      !!                     avs  = avt  +    av_wave * diffusivity_ratio(Reb) 
    122115      !! 
    123       !! ** Action  : - Define emix_iwm used to compute internal wave-induced mixing 
    124       !!              - avt, avs, avm, increased by internal wave-driven mixing     
     116      !! ** Action  : - avt, avs, avm, increased by tide internal wave-driven mixing     
    125117      !! 
    126118      !! References :  de Lavergne et al. 2015, JPO; 2016, in prep. 
    … …  
    132124      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    133125      REAL(wp) ::   zztmp        ! scalar workspace 
    134       REAL(wp), DIMENSION(jpi,jpj)     ::  zfact     ! Used for vertical structure 
    135       REAL(wp), DIMENSION(jpi,jpj)     ::  zhdep     ! Ocean depth 
    136       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zwkb      ! WKB-stretched height above bottom 
    137       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zweight   ! Weight for high mode vertical distribution 
    138       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  znu_t     ! Molecular kinematic viscosity (T grid) 
    139       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  znu_w     ! Molecular kinematic viscosity (W grid) 
    140       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zReb      ! Turbulence intensity parameter 
     126      REAL(wp), DIMENSION(jpi,jpj)     ::   zfact       ! Used for vertical structure 
     127      REAL(wp), DIMENSION(jpi,jpj)     ::   zhdep       ! Ocean depth 
     128      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwkb        ! WKB-stretched height above bottom 
     129      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zweight     ! Weight for high mode vertical distribution 
     130      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znu_t       ! Molecular kinematic viscosity (T grid) 
     131      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znu_w       ! Molecular kinematic viscosity (W grid) 
     132      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zReb        ! Turbulence intensity parameter 
     133      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zemx_iwm    ! local energy density available for mixing (W/kg) 
     134      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zav_ratio   ! S/T diffusivity ratio (only for ln_tsdiff=T) 
     135      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zav_wave    ! Internal wave-induced diffusivity 
     136      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   z3d  ! 3D workspace used for iom_put  
     137      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   z2d  ! 2D     -      -    -     - 
    141138      !!---------------------------------------------------------------------- 
    142139      ! 
    143140      IF( ln_timing )   CALL timing_start('zdf_iwm') 
    144141      ! 
    145       !                      ! ----------------------------- ! 
    146       !                      !  Internal wave-driven mixing  !  (compute zav_wave) 
    147       !                      ! ----------------------------- ! 
     142      !                       !* Set to zero the 1st and last vertical levels of appropriate variables 
     143      zemx_iwm (:,:,1) = 0._wp   ;   zemx_iwm (:,:,jpk) = 0._wp 
     144      zav_ratio(:,:,1) = 0._wp   ;   zav_ratio(:,:,jpk) = 0._wp 
     145      zav_wave (:,:,1) = 0._wp   ;   zav_wave (:,:,jpk) = 0._wp 
     146      ! 
     147      !                       ! ----------------------------- ! 
     148      !                       !  Internal wave-driven mixing  !  (compute zav_wave) 
     149      !                       ! ----------------------------- ! 
    148150      !                              
    149       !                        !* Critical slope mixing: distribute energy over the time-varying ocean depth, 
     151      !                       !* Critical slope mixing: distribute energy over the time-varying ocean depth, 
    150152      !                                                 using an exponential decay from the seafloor. 
    151153      DO jj = 1, jpj                ! part independent of the level 
    … …  
    156158         END DO 
    157159      END DO 
    158  
     160!!gm gde3w ==>>>  check for ssh taken into account.... seem OK gde3w_n=gdept_n - sshn 
    159161      DO jk = 2, jpkm1              ! complete with the level-dependent part 
    160          emix_iwm(:,:,jk) = zfact(:,:) * (  EXP( ( gde3w_n(:,:,jk  ) - zhdep(:,:) ) / hcri_iwm(:,:) )                      & 
     162         zemx_iwm(:,:,jk) = zfact(:,:) * (  EXP( ( gde3w_n(:,:,jk  ) - zhdep(:,:) ) / hcri_iwm(:,:) )                      & 
    161163            &                             - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_iwm(:,:) )  ) * wmask(:,:,jk)   & 
    162164            &                          / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 
    … …  
    186188         ! 
    187189         DO jk = 2, jpkm1              ! complete with the level-dependent part 
    188             emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zfact(:,:) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
     190            zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zfact(:,:) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
    189191         END DO 
    190192         ! 
    … …  
    203205         ! 
    204206         DO jk = 2, jpkm1              ! complete with the level-dependent part 
    205             emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zfact(:,:) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 
     207            zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zfact(:,:) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 
    206208         END DO 
    207209         ! 
    … …  
    253255      ! 
    254256      DO jk = 2, jpkm1              ! complete with the level-dependent part 
    255          emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk)   & 
     257         zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk)   & 
    256258            &                                / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 
    257259!!gm  use of e3t_n just above? 
    … …  
    269271      ! Calculate turbulence intensity parameter Reb 
    270272      DO jk = 2, jpkm1 
    271          zReb(:,:,jk) = emix_iwm(:,:,jk) / MAX( 1.e-20_wp, znu_w(:,:,jk) * rn2(:,:,jk) ) 
     273         zReb(:,:,jk) = zemx_iwm(:,:,jk) / MAX( 1.e-20_wp, znu_w(:,:,jk) * rn2(:,:,jk) ) 
    272274      END DO 
    273275      ! 
    … …  
    349351      !                             !* output internal wave-driven mixing coefficient 
    350352      CALL iom_put( "av_wave", zav_wave ) 
    351                                     !* output useful diagnostics: N^2, Kz * N^2 (bflx_iwm),  
    352                                     !  vertical integral of rau0 * Kz * N^2 (pcmap_iwm), energy density (emix_iwm) 
     353                                    !* output useful diagnostics: Kz*N^2 ,  
     354!!gm Kz*N2 should take into account the ratio avs/avt if it is used.... (see diaar5) 
     355                                    !  vertical integral of rau0 * Kz * N^2 , energy density (zemx_iwm) 
    353356      IF( iom_use("bflx_iwm") .OR. iom_use("pcmap_iwm") ) THEN 
    354          bflx_iwm(:,:,:) = MAX( 0._wp, rn2(:,:,:) ) * zav_wave(:,:,:) 
    355          pcmap_iwm(:,:) = 0._wp 
     357         ALLOCATE( z2d(jpi,jpj) , z3d(jpi,jpj,jpk) ) 
     358         z3d(:,:,:) = MAX( 0._wp, rn2(:,:,:) ) * zav_wave(:,:,:) 
     359         z2d(:,:) = 0._wp 
    356360         DO jk = 2, jpkm1 
    357             pcmap_iwm(:,:) = pcmap_iwm(:,:) + e3w_n(:,:,jk) * bflx_iwm(:,:,jk) * wmask(:,:,jk) 
    358          END DO 
    359          pcmap_iwm(:,:) = rau0 * pcmap_iwm(:,:) 
    360          CALL iom_put( "bflx_iwm", bflx_iwm ) 
    361          CALL iom_put( "pcmap_iwm", pcmap_iwm ) 
    362       ENDIF 
    363       CALL iom_put( "bn2", rn2 ) 
    364       CALL iom_put( "emix_iwm", emix_iwm ) 
     361            z2d(:,:) = z2d(:,:) + e3w_n(:,:,jk) * z3d(:,:,jk) * wmask(:,:,jk) 
     362         END DO 
     363         z2d(:,:) = rau0 * z2d(:,:) 
     364         CALL iom_put( "bflx_iwm", z3d ) 
     365         CALL iom_put( "pcmap_iwm", z2d ) 
     366         DEALLOCATE( z2d , z3d ) 
     367      ENDIF 
     368      CALL iom_put( "emix_iwm", zemx_iwm ) 
    365369       
    366370      IF(ln_ctl)   CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm - av_wave: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk) 
    … …  
    466470      ecri_iwm(:,:) = ecri_iwm(:,:) * ssmask(:,:) 
    467471 
    468       ! Set once for all to zero the first and last vertical levels of appropriate variables 
    469       emix_iwm (:,:, 1 ) = 0._wp 
    470       emix_iwm (:,:,jpk) = 0._wp 
    471       zav_ratio(:,:, 1 ) = 0._wp 
    472       zav_ratio(:,:,jpk) = 0._wp 
    473       zav_wave (:,:, 1 ) = 0._wp 
    474       zav_wave (:,:,jpk) = 0._wp 
    475  
    476472      zbot = glob_sum( e1e2t(:,:) * ebot_iwm(:,:) ) 
    477473      zpyc = glob_sum( e1e2t(:,:) * epyc_iwm(:,:) ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/step.F90

    r8586 r8637  
    6262      !!                     ***  ROUTINE stp  *** 
    6363      !! 
    64       !! ** Purpose : - Time stepping of OPA (momentum and active tracer eqs.) 
    65       !!              - Time stepping of LIM (dynamic and thermodynamic eqs.) 
    66       !!              - Tme stepping  of TRC (passive tracer eqs.) 
     64      !! ** Purpose : - Time stepping of OPA  (momentum and active tracer eqs.) 
     65      !!              - Time stepping of ESIM (dynamic and thermodynamic eqs.) 
     66      !!              - Time stepping of TRC (passive tracer eqs.) 
    6767      !! 
    6868      !! ** Method  : -1- Update forcings and data 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/RK3_SRC/BDY/bdyice.F90

    r8586 r8637  
    1616   !!---------------------------------------------------------------------- 
    1717   USE oce             ! ocean dynamics and tracers variables 
    18    USE ice             ! LIM_3 ice variables 
    19    USE icevar 
    20    USE icectl 
     18   USE ice             ! sea-ice: variables 
     19   USE icevar          ! sea-ice: operations 
     20   USE iceitd          ! sea-ice: rebining 
     21   USE icectl          ! sea-ice: control prints 
    2122   USE phycst          ! physical constant 
    2223   USE eosbn2          ! equation of state 
    … …  
    3940 
    4041   !!---------------------------------------------------------------------- 
    41    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     42   !! NEMO/OPA 4.0 , NEMO Consortium (2017) 
    4243   !! $Id: bdyice.F90 8306 2017-07-10 10:18:03Z clem $ 
    4344   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    249250      END DO !jl 
    250251      ! 
     252      ! --- In case categories are out of bounds, do a remapping --- ! 
     253      !     i.e. inputs have not the same ice thickness distribution  
     254      !          (set by rn_himean) than the regional simulation 
     255      IF( jpl > 1 )   CALL ice_itd_reb( kt ) 
    251256      !       
    252257      IF( nn_timing == 1 )   CALL timing_stop('bdy_ice_frs') 
    … …  
    273278      !!------------------------------------------------------------------------------ 
    274279      ! 
    275       IF( nn_timing == 1 ) CALL timing_start('bdy_ice_dyn') 
     280      IF( ln_timing )  CALL timing_start('bdy_ice_dyn') 
    276281      ! 
    277282      DO ib_bdy=1, nb_bdy 
    … …  
    350355      END DO 
    351356      ! 
    352       IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_dyn') 
     357      IF( ln_timing )  CALL timing_stop('bdy_ice_dyn') 
    353358      ! 
    354359    END SUBROUTINE bdy_ice_dyn 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/RK3_SRC/SBC/sbcblk.F90

    r8586 r8637  
    4040   USE lib_fortran    ! to use key_nosignedzero 
    4141#if defined key_lim3 
    42    USE ice     , ONLY :   u_ice, v_ice, jpl, a_i_b, at_i_b 
     42   USE ice     , ONLY :   u_ice, v_ice, jpl, a_i_b, at_i_b, tm_su 
    4343   USE icethd_dh      ! for CALL ice_thd_snwblow 
    4444#endif 
    … …  
    9292   REAL(wp), PARAMETER ::   Ls     =    2.839e6     ! latent heat of sublimation 
    9393   REAL(wp), PARAMETER ::   Stef   =    5.67e-8     ! Stefan Boltzmann constant 
    94    REAL(wp), PARAMETER ::   Cd_ice =    1.4e-3      ! iovi 1.63e-3     ! transfer coefficient over ice 
     94   REAL(wp), PARAMETER ::   Cd_ice =    1.4e-3      ! transfer coefficient over ice 
    9595   REAL(wp), PARAMETER ::   albo   =    0.066       ! ocean albedo assumed to be constant 
    9696   ! 
    … …  
    108108   REAL(wp) ::   rn_zu          ! z(u)   : height of wind measurements 
    109109!!gm ref namelist initialize it so remove the setting to false below 
    110    LOGICAL  ::   ln_Cd_L12 = .FALSE. !  Modify the drag ice-atm and oce-atm depending on ice concentration (from Lupkes et al. JGR2012) 
     110   LOGICAL  ::   ln_Cd_L12 = .FALSE. !  Modify the drag ice-atm depending on ice concentration (from Lupkes et al. JGR2012) 
     111   LOGICAL  ::   ln_Cd_L15 = .FALSE. !  Modify the drag ice-atm depending on ice concentration (from Lupkes et al. JGR2015) 
    111112   ! 
    112    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_oce   ! air-ocean drag (clem) 
     113!!cr   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_oce   ! air-ocean drag (clem) 
     114   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Cd_atm                    ! transfer coefficient for momentum      (tau) 
     115   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Ch_atm                    ! transfer coefficient for sensible heat (Q_sens) 
     116   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   Ce_atm                    ! tansfert coefficient for evaporation   (Q_lat) 
     117   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_zu                      ! air temperature at wind speed height (needed by Lupkes 2015 bulk scheme) 
     118   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   q_zu                      ! air spec. hum.  at wind speed height (needed by Lupkes 2015 bulk scheme) 
     119   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   cdn_oce, chn_oce, cen_oce ! needed by Lupkes 2015 bulk scheme 
    113120 
    114121   INTEGER  ::   nblk           ! choice of the bulk algorithm 
    … …  
    132139      !!             ***  ROUTINE sbc_blk_alloc *** 
    133140      !!------------------------------------------------------------------- 
    134       ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc ) 
     141!!cr      ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc ) 
     142      ALLOCATE( Cd_atm (jpi,jpj), Ch_atm (jpi,jpj), Ce_atm (jpi,jpj), t_zu(jpi,jpj), q_zu(jpi,jpj), & 
     143         &      cdn_oce(jpi,jpj), chn_oce(jpi,jpj), cen_oce(jpi,jpj), STAT=sbc_blk_alloc ) 
    135144      ! 
    136145      IF( lk_mpp             )   CALL mpp_sum ( sbc_blk_alloc ) 
    … …  
    164173         &                 ln_NCAR, ln_COARE_3p0, ln_COARE_3p5, ln_ECMWF,             &   ! bulk algorithm 
    165174         &                 cn_dir , ln_taudif, rn_zqt, rn_zu,                         &  
    166          &                 rn_pfac, rn_efac, rn_vfac, ln_Cd_L12 
     175         &                 rn_pfac, rn_efac, rn_vfac, ln_Cd_L12, ln_Cd_L15 
    167176      !!--------------------------------------------------------------------- 
    168177      ! 
    … …  
    255264         WRITE(numout,*) '      factor applied on ocean/ice velocity                rn_vfac      = ', rn_vfac 
    256265         WRITE(numout,*) '         (form absolute (=0) to relative winds(=1))' 
     266         WRITE(numout,*) '      use ice-atm drag from Lupkes2012                    ln_Cd_L12    = ', ln_Cd_L12 
     267         WRITE(numout,*) '      use ice-atm drag from Lupkes2015                    ln_Cd_L15    = ', ln_Cd_L15 
    257268         ! 
    258269         WRITE(numout,*) 
    … …  
    361372      REAL(wp), DIMENSION(jpi,jpj) ::   zqlw, zqsb        ! long wave and sensible heat fluxes 
    362373      REAL(wp), DIMENSION(jpi,jpj) ::   zqla, zevap       ! latent heat fluxes and evaporation 
    363       REAL(wp), DIMENSION(jpi,jpj) ::   zCd               ! transfer coefficient for momentum      (tau) 
    364       REAL(wp), DIMENSION(jpi,jpj) ::   zCh               ! transfer coefficient for sensible heat (Q_sens) 
    365       REAL(wp), DIMENSION(jpi,jpj) ::   zCe               ! tansfert coefficient for evaporation   (Q_lat) 
    366374      REAL(wp), DIMENSION(jpi,jpj) ::   zst               ! surface temperature in Kelvin 
    367       REAL(wp), DIMENSION(jpi,jpj) ::   zt_zu             ! air temperature at wind speed height 
    368       REAL(wp), DIMENSION(jpi,jpj) ::   zq_zu             ! air spec. hum.  at wind speed height 
    369375      REAL(wp), DIMENSION(jpi,jpj) ::   zU_zu             ! bulk wind speed at height zu  [m/s] 
    370376      REAL(wp), DIMENSION(jpi,jpj) ::   ztpot             ! potential temperature of air at z=rn_zqt [K] 
    … …  
    418424      zqlw(:,:) = (  sf(jp_qlw)%fnow(:,:,1) - Stef * zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:)  ) * tmask(:,:,1)   ! Long  Wave 
    419425 
    420  
    421  
    422426      ! ----------------------------------------------------------------------------- ! 
    423427      !     II    Turbulent FLUXES                                                    ! 
    … …  
    435439      ! 
    436440      CASE( np_NCAR      )   ;   CALL turb_ncar    ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! NCAR-COREv2 
    437          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     441         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    438442      CASE( np_COARE_3p0 )   ;   CALL turb_coare   ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! COARE v3.0 
    439          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     443         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    440444      CASE( np_COARE_3p5 )   ;   CALL turb_coare3p5( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! COARE v3.5 
    441          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     445         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    442446      CASE( np_ECMWF     )   ;   CALL turb_ecmwf   ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm,   &  ! ECMWF 
    443          &                                              zCd, zCh, zCe, zt_zu, zq_zu, zU_zu ) 
     447         &                                           Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 
    444448      CASE DEFAULT 
    445449         CALL ctl_stop( 'STOP', 'sbc_oce: non-existing bulk formula selected' ) 
    … …  
    448452      !                          ! Compute true air density : 
    449453      IF( ABS(rn_zu - rn_zqt) > 0.01 ) THEN     ! At zu: (probably useless to remove zrho*grav*rn_zu from SLP...) 
    450          zrhoa(:,:) = rho_air( zt_zu(:,:)             , zq_zu(:,:)             , sf(jp_slp)%fnow(:,:,1) ) 
     454         zrhoa(:,:) = rho_air( t_zu(:,:)              , q_zu(:,:)              , sf(jp_slp)%fnow(:,:,1) ) 
    451455      ELSE                                      ! At zt: 
    452456         zrhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1) ) 
    453457      END IF 
    454458 
    455       Cd_oce(:,:) = zCd(:,:)     ! record value of pure ocean-atm. drag (clem) 
     459!!      CALL iom_put( "Cd_oce", Cd_atm)  ! output value of pure ocean-atm. transfer coef. 
     460!!      CALL iom_put( "Ch_oce", Ch_atm)  ! output value of pure ocean-atm. transfer coef. 
    456461 
    457462      DO jj = 1, jpj             ! tau module, i and j component 
    458463         DO ji = 1, jpi 
    459             zztmp = zrhoa(ji,jj)  * zU_zu(ji,jj) * zCd(ji,jj)   ! using bulk wind speed 
     464            zztmp = zrhoa(ji,jj)  * zU_zu(ji,jj) * Cd_atm(ji,jj)   ! using bulk wind speed 
    460465            taum  (ji,jj) = zztmp * wndm  (ji,jj) 
    461466            zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) 
    … …  
    492497      IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 
    493498         !! q_air and t_air are given at 10m (wind reference height) 
    494          zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*zCe(:,:)*(zsq(:,:) - sf(jp_humi)%fnow(:,:,1)) )  ! Evaporation, using bulk wind speed 
    495          zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*zCh(:,:)*(zst(:,:) - ztpot(:,:)               )   ! Sensible Heat, using bulk wind speed 
     499         zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - sf(jp_humi)%fnow(:,:,1)) ) ! Evaporation, using bulk wind speed 
     500         zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - ztpot(:,:)             )   ! Sensible Heat, using bulk wind speed 
    496501      ELSE 
    497502         !! q_air and t_air are not given at 10m (wind reference height) 
    498503         ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 
    499          zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*zCe(:,:)*(zsq(:,:) - zq_zu(:,:) ) ) ! Evaporation ! using bulk wind speed 
    500          zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*zCh(:,:)*(zst(:,:) - zt_zu(:,:) )   ! Sensible Heat ! using bulk wind speed 
     504         zevap(:,:) = rn_efac*MAX( 0._wp,             zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - q_zu(:,:) ) ) ! Evaporation, using bulk wind speed 
     505         zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - t_zu(:,:) )   ! Sensible Heat, using bulk wind speed 
    501506      ENDIF 
    502507 
    … …  
    505510 
    506511      IF(ln_ctl) THEN 
    507          CALL prt_ctl( tab2d_1=zqla  , clinfo1=' blk_oce: zqla   : ', tab2d_2=zCe , clinfo2=' Ce  : ' ) 
    508          CALL prt_ctl( tab2d_1=zqsb  , clinfo1=' blk_oce: zqsb   : ', tab2d_2=zCh , clinfo2=' Ch  : ' ) 
     512         CALL prt_ctl( tab2d_1=zqla  , clinfo1=' blk_oce: zqla   : ', tab2d_2=Ce_atm , clinfo2=' Ce_oce  : ' ) 
     513         CALL prt_ctl( tab2d_1=zqsb  , clinfo1=' blk_oce: zqsb   : ', tab2d_2=Ch_atm , clinfo2=' Ch_oce  : ' ) 
    509514         CALL prt_ctl( tab2d_1=zqlw  , clinfo1=' blk_oce: zqlw   : ', tab2d_2=qsr, clinfo2=' qsr : ' ) 
    510515         CALL prt_ctl( tab2d_1=zsq   , clinfo1=' blk_oce: zsq    : ', tab2d_2=zst, clinfo2=' zst : ' ) 
    … …  
    574579      !!--------------------------------------------------------------------- 
    575580      INTEGER  ::   ji, jj    ! dummy loop indices 
    576       REAL(wp) ::   zwndi_f , zwndj_f, zwnorm_f      ! relative wind module and components at F-point 
    577       REAL(wp) ::   zwndi_t , zwndj_t                ! relative wind components at T-point 
    578       REAL(wp), DIMENSION(jpi,jpj) ::   zCd, zrhoa   ! transfer coefficient for momentum      (tau) 
     581      REAL(wp) ::   zwndi_f , zwndj_f, zwnorm_f   ! relative wind module and components at F-point 
     582      REAL(wp) ::   zwndi_t , zwndj_t             ! relative wind components at T-point 
     583      REAL(wp), DIMENSION(jpi,jpj) ::   zrhoa     ! transfer coefficient for momentum      (tau) 
    579584      !!--------------------------------------------------------------------- 
    580585      ! 
    581586      IF( ln_timing )   CALL timing_start('blk_ice_tau') 
    582587      ! 
    583       IF( ln_Cd_L12 ) THEN 
    584          CALL Cdn10_Lupkes2012( zCd )  ! air-ice drag = F(ice concentration) (see Lupkes et al., 2012) 
    585       ELSE 
    586          zCd(:,:) = Cd_ice             ! constant air-ice drag 
    587       ENDIF 
    588  
    589       ! local scalars ( place there for vector optimisation purposes) 
    590       ! Computing density of air! Way denser that 1.2 over sea-ice !!! 
    591       !! 
    592       zrhoa (:,:) =  rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 
    593  
    594       !!gm brutal.... 
    595       utau_ice  (:,:) = 0._wp 
    596       vtau_ice  (:,:) = 0._wp 
    597       wndm_ice  (:,:) = 0._wp 
    598       !!gm end 
    599  
    600       ! ----------------------------------------------------------------------------- ! 
    601       !    Wind components and module relative to the moving ocean ( U10m - U_ice )   ! 
    602       ! ----------------------------------------------------------------------------- ! 
     588      ! set transfer coefficients to default sea-ice values 
     589      Cd_atm(:,:) = Cd_ice 
     590      Ch_atm(:,:) = Cd_ice 
     591      Ce_atm(:,:) = Cd_ice 
     592 
     593      wndm_ice(:,:) = 0._wp      !!gm brutal.... 
     594 
     595      ! ------------------------------------------------------------ ! 
     596      !    Wind module relative to the moving ice ( U10m - U_ice )   ! 
     597      ! ------------------------------------------------------------ ! 
    603598      SELECT CASE( cp_ice_msh ) 
    604599      CASE( 'I' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
    … …  
    606601         DO jj = 2, jpjm1 
    607602            DO ji = 2, jpim1   ! B grid : NO vector opt 
    608                ! ... scalar wind at I-point (fld being at T-point) 
    609                zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)   & 
    610                   &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) - rn_vfac * u_ice(ji,jj) 
    611                zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
    612                   &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
    613                zwnorm_f = zrhoa(ji,jj) * zCd(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
    614                ! ... ice stress at I-point 
    615                utau_ice(ji,jj) = zwnorm_f * zwndi_f 
    616                vtau_ice(ji,jj) = zwnorm_f * zwndj_f 
    617603               ! ... scalar wind at T-point (fld being at T-point) 
    618604               zwndi_t = sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.25 * (  u_ice(ji,jj+1) + u_ice(ji+1,jj+1)   & 
    … …  
    623609            END DO 
    624610         END DO 
    625          CALL lbc_lnk( utau_ice, 'I', -1. ) 
    626          CALL lbc_lnk( vtau_ice, 'I', -1. ) 
    627611         CALL lbc_lnk( wndm_ice, 'T',  1. ) 
    628612         ! 
    629613      CASE( 'C' )                  ! C-grid ice dynamics :   U & V-points (same as ocean) 
    630          DO jj = 2, jpj 
    631             DO ji = fs_2, jpi   ! vect. opt. 
     614         DO jj = 2, jpjm1 
     615            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    632616               zwndi_t = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( u_ice(ji-1,jj  ) + u_ice(ji,jj) )  ) 
    633617               zwndj_t = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( v_ice(ji  ,jj-1) + v_ice(ji,jj) )  ) 
    … …  
    635619            END DO 
    636620         END DO 
     621         CALL lbc_lnk( wndm_ice, 'T',  1. ) 
     622         ! 
     623      END SELECT 
     624 
     625      ! Make ice-atm. drag dependent on ice concentration 
     626      IF    ( ln_Cd_L12 ) THEN   ! calculate new drag from Lupkes(2012) equations 
     627         CALL Cdn10_Lupkes2012( Cd_atm ) 
     628         Ch_atm(:,:) = Cd_atm(:,:)       ! momentum and heat transfer coef. are considered identical 
     629      ELSEIF( ln_Cd_L15 ) THEN   ! calculate new drag from Lupkes(2015) equations 
     630         CALL Cdn10_Lupkes2015( Cd_atm, Ch_atm )  
     631      ENDIF 
     632 
     633!!      CALL iom_put( "Cd_ice", Cd_atm)  ! output value of pure ice-atm. transfer coef. 
     634!!      CALL iom_put( "Ch_ice", Ch_atm)  ! output value of pure ice-atm. transfer coef. 
     635 
     636      ! local scalars ( place there for vector optimisation purposes) 
     637      ! Computing density of air! Way denser that 1.2 over sea-ice !!! 
     638      zrhoa (:,:) =  rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 
     639 
     640      !!gm brutal.... 
     641      utau_ice  (:,:) = 0._wp 
     642      vtau_ice  (:,:) = 0._wp 
     643      !!gm end 
     644 
     645      ! ------------------------------------------------------------ ! 
     646      !    Wind stress relative to the moving ice ( U10m - U_ice )   ! 
     647      ! ------------------------------------------------------------ ! 
     648      SELECT CASE( cp_ice_msh ) 
     649      CASE( 'I' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
     650         DO jj = 2, jpjm1 
     651            DO ji = 2, jpim1   ! B grid : NO vector opt 
     652               ! ... scalar wind at I-point (fld being at T-point) 
     653               zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)   & 
     654                  &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) - rn_vfac * u_ice(ji,jj) 
     655               zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
     656                  &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
     657               ! ... ice stress at I-point 
     658               zwnorm_f = zrhoa(ji,jj) * Cd_atm(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
     659               utau_ice(ji,jj) = zwnorm_f * zwndi_f 
     660               vtau_ice(ji,jj) = zwnorm_f * zwndj_f 
     661            END DO 
     662         END DO 
     663         CALL lbc_lnk( utau_ice, 'I', -1. ) 
     664         CALL lbc_lnk( vtau_ice, 'I', -1. ) 
     665         ! 
     666      CASE( 'C' )                  ! C-grid ice dynamics :   U & V-points (same as ocean) 
    637667         DO jj = 2, jpjm1 
    638668            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    639                utau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * zCd(ji,jj) * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )                          & 
     669               utau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd_atm(ji,jj) * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )            & 
    640670                  &          * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 
    641                vtau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * zCd(ji,jj) * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )                          & 
     671               vtau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd_atm(ji,jj) * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )            & 
    642672                  &          * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 
    643673            END DO 
    … …  
    645675         CALL lbc_lnk( utau_ice, 'U', -1. ) 
    646676         CALL lbc_lnk( vtau_ice, 'V', -1. ) 
    647          CALL lbc_lnk( wndm_ice, 'T',  1. ) 
    648677         ! 
    649678      END SELECT 
    … …  
    684713      REAL(wp), DIMENSION(jpi,jpj)     ::   zevap, zsnw   ! evaporation and snw distribution after wind blowing (LIM3) 
    685714      REAL(wp), DIMENSION(jpi,jpj)     ::   zrhoa 
    686       REAL(wp), DIMENSION(jpi,jpj)     ::   zCd            ! transfer coefficient for momentum      (tau) 
    687715      !!--------------------------------------------------------------------- 
    688716      ! 
    689717      IF( ln_timing )   CALL timing_start('blk_ice_flx') 
    690718      ! 
    691       IF( ln_Cd_L12 ) THEN 
    692          CALL Cdn10_Lupkes2012( zCd )  ! air-ice drag = F(ice concentration) (see Lupkes et al., 2012) 
    693       ELSE 
    694          zCd(:,:) = Cd_ice             ! constant air-ice drag 
    695       ENDIF 
    696  
    697       ! 
    698       ! local scalars ( place there for vector optimisation purposes) 
     719      ! 
     720      ! local scalars 
    699721      zcoef_dqlw   = 4.0 * 0.95 * Stef 
    700722      zcoef_dqla   = -Ls * 11637800. * (-5897.8) 
    … …  
    724746               ! ----------------------------! 
    725747 
    726                ! ... turbulent heat fluxes 
     748               ! ... turbulent heat fluxes with Ch_atm recalculated in blk_ice_tau 
    727749               ! Sensible Heat 
    728                z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * zCd(ji,jj) * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 
     750               z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Ch_atm(ji,jj) * wndm_ice(ji,jj) * (ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1)) 
    729751               ! Latent Heat 
    730                qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls  * zCd(ji,jj) * wndm_ice(ji,jj)   & 
    731                   &                         * (  11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 
     752               qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls  * Ch_atm(ji,jj) * wndm_ice(ji,jj) *  & 
     753                  &                ( 11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 
    732754               ! Latent heat sensitivity for ice (Dqla/Dt) 
    733755               IF( qla_ice(ji,jj,jl) > 0._wp ) THEN 
    734                   dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * zCd(ji,jj) * wndm_ice(ji,jj) / ( zst2 ) * EXP( -5897.8 / ptsu(ji,jj,jl) ) 
     756                  dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * Ch_atm(ji,jj) * wndm_ice(ji,jj) / zst2 * EXP(-5897.8 / ptsu(ji,jj,jl)) 
    735757               ELSE 
    736758                  dqla_ice(ji,jj,jl) = 0._wp 
    … …  
    738760 
    739761               ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 
    740                z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * zCd(ji,jj) * wndm_ice(ji,jj) 
     762               z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Ch_atm(ji,jj) * wndm_ice(ji,jj) 
    741763 
    742764               ! ----------------------------! 
    … …  
    935957#if defined key_lim3 
    936958 
    937    SUBROUTINE Cdn10_Lupkes2012( pCd ) 
     959   SUBROUTINE Cdn10_Lupkes2012( Cd ) 
    938960      !!---------------------------------------------------------------------- 
    939961      !!                      ***  ROUTINE  Cdn10_Lupkes2012  *** 
    … …  
    965987      !! 
    966988      !!---------------------------------------------------------------------- 
    967       REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pCd   ! air-ice drag coefficient 
     989      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Cd 
    968990      REAL(wp), PARAMETER ::   zCe   = 2.23e-03_wp 
    969991      REAL(wp), PARAMETER ::   znu   = 1._wp 
    … …  
    973995      !!---------------------------------------------------------------------- 
    974996      zcoef = znu + 1._wp / ( 10._wp * zbeta ) 
    975       ! 
     997 
    976998      ! generic drag over a cell partly covered by ice 
    977       !!pCd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) +  &                        ! pure ocean drag 
    978       !!   &       Cd_ice      *           at_i_b(:,:)   +  &                        ! pure ice drag 
    979       !!   &       zCe         * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu   ! change due to sea-ice morphology 
     999      !!Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) +  &                        ! pure ocean drag 
     1000      !!   &      Cd_ice      *           at_i_b(:,:)   +  &                        ! pure ice drag 
     1001      !!   &      zCe         * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu   ! change due to sea-ice morphology 
    9801002 
    9811003      ! ice-atm drag 
    982       pCd(:,:) = Cd_ice +  &                                                         ! pure ice drag 
    983          &       zCe    * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp)  ! change due to sea-ice morphology 
    984       ! 
     1004      Cd(:,:) = Cd_ice +  &                                                         ! pure ice drag 
     1005         &      zCe    * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp)  ! change due to sea-ice morphology 
     1006       
    9851007   END SUBROUTINE Cdn10_Lupkes2012 
     1008 
     1009 
     1010   SUBROUTINE Cdn10_Lupkes2015( Cd, Ch ) 
     1011      !!---------------------------------------------------------------------- 
     1012      !!                      ***  ROUTINE  Cdn10_Lupkes2015  *** 
     1013      !! 
     1014      !! ** pUrpose :    1lternative turbulent transfert coefficients formulation 
     1015      !!                 between sea-ice and atmosphere with distinct momentum  
     1016      !!                 and heat coefficients depending on sea-ice concentration  
     1017      !!                 and atmospheric stability (no meltponds effect for now). 
     1018      !!                 
     1019      !! ** Method :     The parameterization is adapted from Lupkes et al. (2015) 
     1020      !!                 and ECHAM6 atmospheric model. Compared to Lupkes2012 scheme, 
     1021      !!                 it considers specific skin and form drags (Andreas et al. 2010) 
     1022      !!                 to compute neutral transfert coefficients for both heat and  
     1023      !!                 momemtum fluxes. Atmospheric stability effect on transfert 
     1024      !!                 coefficient is also taken into account following Louis (1979). 
     1025      !! 
     1026      !! ** References : Lupkes et al. JGR 2015 (theory) 
     1027      !!                 Lupkes et al. ECHAM6 documentation 2015 (implementation) 
     1028      !! 
     1029      !!---------------------------------------------------------------------- 
     1030      ! 
     1031      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Cd 
     1032      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   Ch 
     1033      REAL(wp), DIMENSION(jpi,jpj)            ::   zst, zqo_sat, zqi_sat 
     1034      ! 
     1035      ! ECHAM6 constants 
     1036      REAL(wp), PARAMETER ::   z0_skin_ice  = 0.69e-3_wp  ! Eq. 43 [m] 
     1037      REAL(wp), PARAMETER ::   z0_form_ice  = 0.57e-3_wp  ! Eq. 42 [m] 
     1038      REAL(wp), PARAMETER ::   z0_ice       = 1.00e-3_wp  ! Eq. 15 [m] 
     1039      REAL(wp), PARAMETER ::   zce10        = 2.80e-3_wp  ! Eq. 41 
     1040      REAL(wp), PARAMETER ::   zbeta        = 1.1_wp      ! Eq. 41 
     1041      REAL(wp), PARAMETER ::   zc           = 5._wp       ! Eq. 13 
     1042      REAL(wp), PARAMETER ::   zc2          = zc * zc 
     1043      REAL(wp), PARAMETER ::   zam          = 2. * zc     ! Eq. 14 
     1044      REAL(wp), PARAMETER ::   zah          = 3. * zc     ! Eq. 30 
     1045      REAL(wp), PARAMETER ::   z1_alpha     = 1._wp / 0.2_wp  ! Eq. 51 
     1046      REAL(wp), PARAMETER ::   z1_alphaf    = z1_alpha    ! Eq. 56 
     1047      REAL(wp), PARAMETER ::   zbetah       = 1.e-3_wp    ! Eq. 26 
     1048      REAL(wp), PARAMETER ::   zgamma       = 1.25_wp     ! Eq. 26 
     1049      REAL(wp), PARAMETER ::   z1_gamma     = 1._wp / zgamma 
     1050      REAL(wp), PARAMETER ::   r1_3         = 1._wp / 3._wp 
     1051      ! 
     1052      INTEGER  ::   ji, jj         ! dummy loop indices 
     1053      REAL(wp) ::   zthetav_os, zthetav_is, zthetav_zu 
     1054      REAL(wp) ::   zrib_o, zrib_i 
     1055      REAL(wp) ::   zCdn_skin_ice, zCdn_form_ice, zCdn_ice 
     1056      REAL(wp) ::   zChn_skin_ice, zChn_form_ice 
     1057      REAL(wp) ::   z0w, z0i, zfmi, zfmw, zfhi, zfhw 
     1058      REAL(wp) ::   zCdn_form_tmp 
     1059      !!---------------------------------------------------------------------- 
     1060 
     1061      ! Momentum Neutral Transfert Coefficients (should be a constant) 
     1062      zCdn_form_tmp = zce10 * ( LOG( 10._wp / z0_form_ice + 1._wp ) / LOG( rn_zu / z0_form_ice + 1._wp ) )**2   ! Eq. 40 
     1063      zCdn_skin_ice = ( vkarmn                                      / LOG( rn_zu / z0_skin_ice + 1._wp ) )**2   ! Eq. 7 
     1064      zCdn_ice      = zCdn_skin_ice   ! Eq. 7 (cf Lupkes email for details) 
     1065      !zCdn_ice     = 1.89e-3         ! old ECHAM5 value (cf Eq. 32) 
     1066 
     1067      ! Heat Neutral Transfert Coefficients 
     1068      zChn_skin_ice = vkarmn**2 / ( LOG( rn_zu / z0_ice + 1._wp ) * LOG( rn_zu * z1_alpha / z0_skin_ice + 1._wp ) )   ! Eq. 50 + Eq. 52 (cf Lupkes email for details) 
     1069      
     1070      ! Atmospheric and Surface Variables 
     1071      zst(:,:)     = sst_m(:,:) + rt0                                       ! convert SST from Celcius to Kelvin 
     1072      zqo_sat(:,:) = 0.98_wp * q_sat( zst(:,:)  , sf(jp_slp)%fnow(:,:,1) )  ! saturation humidity over ocean [kg/kg] 
     1073      zqi_sat(:,:) = 0.98_wp * q_sat( tm_su(:,:), sf(jp_slp)%fnow(:,:,1) )  ! saturation humidity over ice   [kg/kg] 
     1074      ! 
     1075      DO jj = 2, jpjm1           ! reduced loop is necessary for reproducibility 
     1076         DO ji = fs_2, fs_jpim1 
     1077            ! Virtual potential temperature [K] 
     1078            zthetav_os = zst(ji,jj)   * ( 1._wp + rctv0 * zqo_sat(ji,jj) )   ! over ocean 
     1079            zthetav_is = tm_su(ji,jj) * ( 1._wp + rctv0 * zqi_sat(ji,jj) )   ! ocean ice 
     1080            zthetav_zu = t_zu (ji,jj) * ( 1._wp + rctv0 * q_zu(ji,jj)    )   ! at zu 
     1081             
     1082            ! Bulk Richardson Number (could use Ri_bulk function from aerobulk instead) 
     1083            zrib_o = grav / zthetav_os * ( zthetav_zu - zthetav_os ) * rn_zu / MAX( 0.5, wndm(ji,jj)     )**2   ! over ocean 
     1084            zrib_i = grav / zthetav_is * ( zthetav_zu - zthetav_is ) * rn_zu / MAX( 0.5, wndm_ice(ji,jj) )**2   ! over ice 
     1085             
     1086            ! Momentum and Heat Neutral Transfert Coefficients 
     1087            zCdn_form_ice = zCdn_form_tmp * at_i_b(ji,jj) * ( 1._wp - at_i_b(ji,jj) )**zbeta  ! Eq. 40 
     1088            zChn_form_ice = zCdn_form_ice / ( 1._wp + ( LOG( z1_alphaf ) / vkarmn ) * SQRT( zCdn_form_ice ) )               ! Eq. 53  
     1089                        
     1090            ! Momentum and Heat Stability functions (possibility to use psi_m_ecmwf instead) 
     1091            z0w = rn_zu * EXP( -1._wp * vkarmn / SQRT( Cdn_oce(ji,jj) ) ) ! over water 
     1092            z0i = z0_skin_ice                                             ! over ice (cf Lupkes email for details) 
     1093            IF( zrib_o <= 0._wp ) THEN 
     1094               zfmw = 1._wp - zam * zrib_o / ( 1._wp + 3._wp * zc2 * Cdn_oce(ji,jj) * SQRT( -zrib_o * ( rn_zu / z0w + 1._wp ) ) )  ! Eq. 10 
     1095               zfhw = ( 1._wp + ( zbetah * ( zthetav_os - zthetav_zu )**r1_3 / ( Chn_oce(ji,jj) * MAX(0.01, wndm(ji,jj)) )   &     ! Eq. 26 
     1096                  &             )**zgamma )**z1_gamma 
     1097            ELSE 
     1098               zfmw = 1._wp / ( 1._wp + zam * zrib_o / SQRT( 1._wp + zrib_o ) )   ! Eq. 12 
     1099               zfhw = 1._wp / ( 1._wp + zah * zrib_o / SQRT( 1._wp + zrib_o ) )   ! Eq. 28 
     1100            ENDIF 
     1101             
     1102            IF( zrib_i <= 0._wp ) THEN 
     1103               zfmi = 1._wp - zam * zrib_i / (1._wp + 3._wp * zc2 * zCdn_ice * SQRT( -zrib_i * ( rn_zu / z0i + 1._wp)))   ! Eq.  9 
     1104               zfhi = 1._wp - zah * zrib_i / (1._wp + 3._wp * zc2 * zCdn_ice * SQRT( -zrib_i * ( rn_zu / z0i + 1._wp)))   ! Eq. 25 
     1105            ELSE 
     1106               zfmi = 1._wp / ( 1._wp + zam * zrib_i / SQRT( 1._wp + zrib_i ) )   ! Eq. 11 
     1107               zfhi = 1._wp / ( 1._wp + zah * zrib_i / SQRT( 1._wp + zrib_i ) )   ! Eq. 27 
     1108            ENDIF 
     1109             
     1110            ! Momentum Transfert Coefficients (Eq. 38) 
     1111            Cd(ji,jj) = zCdn_skin_ice *   zfmi +  & 
     1112               &        zCdn_form_ice * ( zfmi * at_i_b(ji,jj) + zfmw * ( 1._wp - at_i_b(ji,jj) ) ) / MAX( 1.e-06, at_i_b(ji,jj) ) 
     1113             
     1114            ! Heat Transfert Coefficients (Eq. 49) 
     1115            Ch(ji,jj) = zChn_skin_ice *   zfhi +  & 
     1116               &        zChn_form_ice * ( zfhi * at_i_b(ji,jj) + zfhw * ( 1._wp - at_i_b(ji,jj) ) ) / MAX( 1.e-06, at_i_b(ji,jj) ) 
     1117            ! 
     1118         END DO 
     1119      END DO 
     1120      CALL lbc_lnk_multi( Cd, 'T',  1., Ch, 'T', 1. ) 
     1121      ! 
     1122   END SUBROUTINE Cdn10_Lupkes2015 
     1123 
    9861124#endif 
    987     
    9881125 
    9891126   !!====================================================================== 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/RK3_SRC/SBC/sbcblk_algo_coare.F90

    r7646 r8637  
    4444   USE lib_fortran    ! to use key_nosignedzero 
    4545 
    46  
    4746   IMPLICIT NONE 
    4847   PRIVATE 
    … …  
    5049   PUBLIC ::   TURB_COARE   ! called by sbcblk.F90 
    5150 
    52    !! COARE own values for given constants: 
    53    REAL(wp), PARAMETER :: & 
    54       &   zi0     = 600.,  &   !: scale height of the atmospheric boundary layer...1 
    55       &  Beta0    = 1.25, &   !: gustiness parameter 
    56       &  rctv0    = 0.608     !: constant to obtain virtual temperature... 
    57  
    58  
     51   !                                           !! COARE own values for given constants: 
     52   REAL(wp), PARAMETER ::   zi0     = 600.      ! scale height of the atmospheric boundary layer...1 
     53   REAL(wp), PARAMETER ::   Beta0   =   1.25    ! gustiness parameter 
     54   REAL(wp), PARAMETER ::   rctv0   =   0.608   ! constant to obtain virtual temperature... 
     55 
     56   !!---------------------------------------------------------------------- 
    5957CONTAINS 
    6058 
    6159   SUBROUTINE turb_coare( zt, zu, sst, t_zt, ssq, q_zt, U_zu, & 
    62       &                   Cd, Ch, Ce, t_zu, q_zu, U_blk ) 
     60      &                   Cd, Ch, Ce, t_zu, q_zu, U_blk,      & 
     61      &                   Cdn, Chn, Cen                       ) 
    6362      !!---------------------------------------------------------------------- 
    6463      !!                      ***  ROUTINE  turb_coare  *** 
    … …  
    106105      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu           [kg/kg] 
    107106      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     107      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    108108      ! 
    109109      INTEGER :: j_itt 
    … …  
    246246      Ce   = ztmp0*q_star/dq_zu 
    247247      ! 
     248      ztmp1 = zu + z0 
     249      Cdn = vkarmn*vkarmn / (log(ztmp1/z0 )*log(ztmp1/z0 )) 
     250      Chn = vkarmn*vkarmn / (log(ztmp1/z0t)*log(ztmp1/z0t)) 
     251      Cen = Chn 
     252      ! 
    248253      CALL wrk_dealloc( jpi,jpj, u_star, t_star, q_star, zeta_u, dt_zu, dq_zu ) 
    249254      CALL wrk_dealloc( jpi,jpj, znu_a, z0, z0t, ztmp0, ztmp1, ztmp2 ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/RK3_SRC/SBC/sbcblk_algo_coare3p5.F90

    r7646 r8637  
    5858CONTAINS 
    5959 
    60    SUBROUTINE turb_coare3p5( zt, zu, sst, t_zt, ssq, q_zt, U_zu, & 
    61       &                   Cd, Ch, Ce, t_zu, q_zu, U_blk ) 
     60   SUBROUTINE turb_coare3p5( zt, zu, sst, t_zt, ssq, q_zt, U_zu,  & 
     61      &                      Cd, Ch, Ce, t_zu, q_zu, U_blk,       & 
     62      &                      Cdn, Chn, Cen                        ) 
    6263      !!---------------------------------------------------------------------------------- 
    6364      !!                      ***  ROUTINE  turb_coare3p5  *** 
    … …  
    105106      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu             [kg/kg] 
    106107      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     108      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    107109      ! 
    108110      INTEGER :: j_itt 
    … …  
    252254      Ce   = ztmp0*q_star/dq_zu 
    253255      ! 
     256      ztmp1 = zu + z0 
     257      Cdn = vkarmn*vkarmn / (log(ztmp1/z0 )*log(ztmp1/z0 )) 
     258      Chn = vkarmn*vkarmn / (log(ztmp1/z0t)*log(ztmp1/z0t)) 
     259      Cen = Chn 
     260      ! 
    254261      CALL wrk_dealloc( jpi,jpj, u_star, t_star, q_star, zeta_u, dt_zu, dq_zu ) 
    255262      CALL wrk_dealloc( jpi,jpj, znu_a, z0, z0t, ztmp0, ztmp1, ztmp2 ) 

  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/RK3_SRC/SBC/sbcblk_algo_ecmwf.F90

    r7646 r8637  
    6464 
    6565   SUBROUTINE TURB_ECMWF( zt, zu, sst, t_zt, ssq , q_zt , U_zu,   & 
    66       &                   Cd, Ch, Ce , t_zu, q_zu, U_blk ) 
     66      &                   Cd, Ch, Ce , t_zu, q_zu, U_blk,         & 
     67      &                   Cdn, Chn, Cen                           ) 
    6768      !!---------------------------------------------------------------------------------- 
    6869      !!                      ***  ROUTINE  turb_ecmwf  *** 
    … …  
    112113      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   q_zu     ! spec. humidity adjusted at zu           [kg/kg] 
    113114      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   U_blk    ! bulk wind at 10m                          [m/s] 
     115      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   Cdn, Chn, Cen ! neutral transfer coefficients 
    114116      ! 
    115117      INTEGER :: j_itt 
    … …  
    266268            dt_zu = t_zu - sst ;  dt_zu = SIGN( MAX(ABS(dt_zu),1.E-6), dt_zu ) 
    267269            dq_zu = q_zu - ssq ;  dq_zu = SIGN( MAX(ABS(dq_zu),1.E-9), dq_zu ) 
     270 
    268271         END IF 
    269272 
    … …  
    271274         ztmp1 = zu + z0 
    272275         ztmp0 = ztmp1*Linv 
    273          func_m = log(ztmp1) - LOG(z0) - psi_m_ecmwf(ztmp0) + psi_m_ecmwf(z0*Linv) 
     276         func_m = log(ztmp1) - LOG(z0 ) - psi_m_ecmwf(ztmp0) + psi_m_ecmwf(z0 *Linv) 
    274277         func_h = log(ztmp1) - LOG(z0t) - psi_h_ecmwf(ztmp0) + psi_h_ecmwf(z0t*Linv) 
    275278 
    … …  
    280283      ztmp1 = log((zu + z0)/z0q) -