Changeset 3625
- Timestamp:
- 2012-11-21T14:19:18+01:00 (10 years ago)
- Location:
- branches/2012/dev_NOC_2012_rev3555
- Files:
-
- 106 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2012/dev_NOC_2012_rev3555/DOC/TexFiles/Chapters/Chap_SBC.tex
r3609 r3625 1146 1146 \label{SBC_cice} 1147 1147 1148 It is now possible to couple a global NEMO configuration (without AGRIF) to the CICE sea-ice1148 It is now possible to couple a regional or global NEMO configuration (without AGRIF) to the CICE sea-ice 1149 1149 model by using \key{cice}. The CICE code can be obtained from 1150 1150 \href{http://oceans11.lanl.gov/trac/CICE/}{LANL} and the additional 'hadgem3' drivers will be required, -
branches/2012/dev_NOC_2012_rev3555/DOC/TexFiles/Chapters/Introduction.tex
r3308 r3625 63 63 \citep{OASIS2006}. Two-way nesting is also available through an interface to the 64 64 AGRIF package (Adaptative Grid Refinement in \textsc{Fortran}) \citep{Debreu_al_CG2008}. 65 The interface code for coupling to an alternative sea ice model (CICE, \citet{Hunke2008}) is now 66 available although this is currently only designed for global domains, without the use of AGRIF. 65 The interface code for coupling to an alternative sea ice model (CICE, \citet{Hunke2008}) 66 has now been upgraded so that it works for both global and regional domains, although AGRIF 67 is still not available. 67 68 68 69 Other model characteristics are the lateral boundary conditions (chapter~\ref{LBC}). -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/AMM12/EXP00/namelist
r3609 r3625 137 137 ! =1 use observed ice-cover , 138 138 ! =2 ice-model used ("key_lim3" or "key_lim2) 139 nn_ice_embd = 0 ! =0 levitating ice (no mass exchange, concentration/dilution effect) 140 ! =1 levitating ice with mass and salt exchange but no presure effect 141 ! =2 embedded sea-ice (full salt and mass exchanges and pressure) 139 142 ln_dm2dc = .false. ! daily mean to diurnal cycle on short wave 140 143 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/AMM12_PISCES/EXP00/namelist
r3609 r3625 137 137 ! =1 use observed ice-cover , 138 138 ! =2 ice-model used ("key_lim3" or "key_lim2) 139 nn_ice_embd = 0 ! =0 levitating ice (no mass exchange, concentration/dilution effect) 140 ! =1 levitating ice with mass and salt exchange but no presure effect 141 ! =2 embedded sea-ice (full salt and mass exchanges and pressure) 139 142 ln_dm2dc = .false. ! daily mean to diurnal cycle on short wave 140 143 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/GYRE/EXP00/namelist
r3614 r3625 137 137 ! =1 use observed ice-cover , 138 138 ! =2 ice-model used ("key_lim3" or "key_lim2) 139 nn_ice_embd = 0 ! =0 levitating ice (no mass exchange, concentration/dilution effect) 140 ! =1 levitating ice with mass and salt exchange but no presure effect 141 ! =2 embedded sea-ice (full salt and mass exchanges and pressure) 139 142 ln_dm2dc = .false. ! daily mean to diurnal cycle on short wave 140 143 ln_rnf = .false. ! runoffs (T => fill namsbc_rnf) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/iodef.xml
r3294 r3625 124 124 125 125 <field id="empmr" description="Net Upward Water Flux" unit="kg/m2/s" /> 126 <field id=" empsmr" description="concentration/dilution water flux" unit="kg/m2/s"/>126 <field id="saltflx" description="Downward Salt Flux" unit="PSU/m2/s" /> 127 127 <field id="snowpre" description="Snow precipitation" unit="kg/m2/s" /> 128 128 <field id="runoffs" description="River Runoffs" unit="Kg/m2/s" /> … … 145 145 <field id="qsb_oce" description="Sensible Downward Heat Flux over open ocean" unit="W/m2" /> 146 146 <field id="qla_oce" description="Latent Downward Heat Flux over open ocean" unit="W/m2" /> 147 <field id="qhc_oce" description="Downward Heat Content of E-P over open ocean" unit="W/m2" /> 147 148 <field id="taum_oce" description="wind stress module over open ocean" unit="N/m2" /> 148 149 … … 173 174 <field id="v_imasstr" description="Sea-ice mass transport along j-axis" unit="kg/s" /> 174 175 176 <!-- available if not defined key_vvl --> 177 <field id="emp_x_sst" description="Concentration/Dilution term on SST" unit="kgC/m2/s" /> 178 <field id="emp_x_sss" description="Concentration/Dilution term on SSS" unit="kgPSU/m2/s" /> 175 179 <!-- available key_coupled --> 176 180 <field id="snow_ao_cea" description="Snow over ice-free ocean (cell average)" unit="kg/m2/s" /> … … 1016 1020 <field ref="empmr" name="sowaflup" /> 1017 1021 <field ref="qsr" name="soshfldo" /> 1018 <field ref=" empsmr" name="sowaflcd" />1022 <field ref="saltflx" name="sosfldow" /> 1019 1023 <field ref="qt" name="sohefldo" /> 1020 1024 <field ref="mldr10_1" name="somxl010" /> -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
r3609 r3625 137 137 ! =1 use observed ice-cover , 138 138 ! =2 ice-model used ("key_lim3" or "key_lim2) 139 nn_ice_embd = 0 ! =0 levitating ice (no mass exchange, concentration/dilution effect) 140 ! =1 levitating ice with mass and salt exchange but no presure effect 141 ! =2 embedded sea-ice (full salt and mass exchanges and pressure) 139 142 ln_dm2dc = .false. ! daily mean to diurnal cycle on short wave 140 143 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist
r3614 r3625 137 137 ! =1 use observed ice-cover , 138 138 ! =2 ice-model used ("key_lim3" or "key_lim2) 139 nn_ice_embd = 0 ! =0 levitating ice (no mass exchange, concentration/dilution effect) 140 ! =1 levitating ice with mass and salt exchange but no presure effect 141 ! =2 embedded sea-ice (full salt and mass exchanges and pressure) 139 142 ln_dm2dc = .false. ! daily mean to diurnal cycle on short wave 140 143 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) … … 650 653 sn_mld = 'dyna_grid_T' , 120 , 'somixhgt' , .true. , .true. , 'yearly' , '' , '' 651 654 sn_emp = 'dyna_grid_T' , 120 , 'sowaflcd' , .true. , .true. , 'yearly' , '' , '' 655 ! sn_emp = 'dyna_grid_T' , 120 , 'sowaflup' , .true. , .true. , 'yearly' , '' , '' ! v3.5+ 656 ! sn_sfx = 'dyna_grid_T' , 120 , 'sosfldow' , .true. , .true. , 'yearly' , '' , '' ! v3.5+ 652 657 sn_ice = 'dyna_grid_T' , 120 , 'soicecov' , .true. , .true. , 'yearly' , '' , '' 653 658 sn_qsr = 'dyna_grid_T' , 120 , 'soshfldo' , .true. , .true. , 'yearly' , '' , '' -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/ice_2.F90
r2715 r3625 19 19 PUBLIC ice_alloc_2 ! Called in iceini_2.F90 20 20 21 INTEGER , PUBLIC :: numit !: ice iteration index22 REAL(wp), PUBLIC :: rdt_ice !: ice time step21 INTEGER , PUBLIC :: numit !: ice iteration index 22 REAL(wp), PUBLIC :: rdt_ice !: ice time step 23 23 24 24 ! !!* namicerun read in iceini * … … 98 98 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qstoif !: Energy stored in the brine pockets 99 99 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fbif !: Heat flux at the ice base 100 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdmsnif !: Variation of snow mass 101 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdmicif !: Variation of ice mass 100 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_snw !: Variation of snow mass over 1 time step [Kg/m2] 101 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_snw !: Heat content associated with rdm_snw [J/m2] 102 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_ice !: Variation of ice mass over 1 time step [Kg/m2] 103 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_ice !: Heat content associated with rdm_ice [J/m2] 102 104 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qldif !: heat balance of the lead (or of the open ocean) 103 105 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qcmif !: Energy needed to freeze the ocean surface layer … … 153 155 154 156 ALLOCATE(phicif(jpi,jpj) , pfrld (jpi,jpj) , qstoif (jpi,jpj) , & 155 & fbif (jpi,jpj) , rdmsnif(jpi,jpj) , rdmicif(jpi,jpj) , & 157 & fbif (jpi,jpj) , rdm_snw(jpi,jpj) , rdq_snw(jpi,jpj) , & 158 & rdm_ice(jpi,jpj) , rdq_ice(jpi,jpj) , & 156 159 & qldif (jpi,jpj) , qcmif (jpi,jpj) , fdtcn (jpi,jpj) , & 157 160 & qdtcn (jpi,jpj) , thcm (jpi,jpj) , STAT=ierr(4) ) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/iceini_2.F90
r3294 r3625 13 13 !! 'key_lim2' : LIM 2.0 sea-ice model 14 14 !!---------------------------------------------------------------------- 15 !! ice_init_2 16 !! ice_run_2 15 !! ice_init_2 : sea-ice model initialization 16 !! ice_run_2 : Definition some run parameter for ice model 17 17 !!---------------------------------------------------------------------- 18 USE phycst ! physical constants 19 USE dom_oce ! ocean domain 20 USE sbc_oce ! surface boundary condition: ocean 21 USE sbc_ice ! LIM2 surface boundary condition 22 USE dom_ice_2 ! LIM2 ice domain 23 USE par_ice_2 ! LIM2 parameters 24 USE thd_ice_2 ! LIM2 thermodynamical variables 25 USE ice_2 ! LIM2 ice variable 26 USE limmsh_2 ! LIM2 mesh 27 USE limistate_2 ! LIM2 initial state 28 USE limrst_2 ! LIM2 restart 29 USE limsbc_2 ! LIM2 surface boundary condition 30 USE in_out_manager ! I/O manager 31 USE lib_mpp ! MPP library 18 USE phycst ! physical constants 19 USE dom_oce ! ocean domain 20 USE sbc_oce ! surface boundary condition: ocean 21 USE sbc_ice ! LIM2 surface boundary condition 22 USE dom_ice_2 ! LIM2 ice domain 23 USE par_ice_2 ! LIM2 parameters 24 USE thd_ice_2 ! LIM2 thermodynamical variables 25 USE ice_2 ! LIM2 ice variable 26 USE limmsh_2 ! LIM2 mesh 27 USE limistate_2 ! LIM2 initial state 28 USE limrst_2 ! LIM2 restart 29 USE limsbc_2 ! LIM2 surface boundary condition 30 USE in_out_manager ! I/O manager 31 USE lib_mpp ! MPP library 32 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 32 33 33 34 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limadv_2.F90
r3294 r3625 14 14 !! 'key_lim2' LIM 2.0 sea-ice model 15 15 !!---------------------------------------------------------------------- 16 !! lim_adv_x_2 : advection of sea ice on x axis17 !! lim_adv_y_2 : advection of sea ice on y axis16 !! lim_adv_x_2 : advection of sea ice on x axis 17 !! lim_adv_y_2 : advection of sea ice on y axis 18 18 !!---------------------------------------------------------------------- 19 19 USE dom_oce … … 21 21 USE ice_2 22 22 USE lbclnk 23 USE in_out_manager ! I/O manager 24 USE lib_mpp ! MPP library 25 USE wrk_nemo ! work arrays 26 USE prtctl ! Print control 23 USE in_out_manager ! I/O manager 24 USE lib_mpp ! MPP library 25 USE wrk_nemo ! work arrays 26 USE prtctl ! Print control 27 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 27 28 28 29 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limdia_2.F90
r2715 r3625 24 24 USE in_out_manager ! I/O manager 25 25 USE lib_mpp ! MPP library 26 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 26 27 27 28 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limdmp_2.F90
r2715 r3625 19 19 USE in_out_manager ! I/O manager 20 20 USE lib_mpp ! MPP library 21 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 21 22 22 23 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limdyn_2.F90
r3294 r3625 31 31 USE in_out_manager ! I/O manager 32 32 USE prtctl ! Print control 33 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 33 34 34 35 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limhdf_2.F90
r3294 r3625 21 21 USE prtctl ! Print control 22 22 USE in_out_manager ! I/O manager 23 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 23 24 24 25 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limistate_2.F90
r3294 r3625 27 27 USE iom 28 28 USE in_out_manager 29 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 29 30 30 31 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limmsh_2.F90
r3294 r3625 23 23 USE wrk_nemo ! work arrays 24 24 #endif 25 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 25 26 26 27 IMPLICIT NONE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limrhg_2.F90
r3294 r3625 30 30 USE in_out_manager ! I/O manager 31 31 USE prtctl ! Print control 32 USE oce , ONLY : snwice_mass, snwice_mass_b 33 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 32 34 33 35 IMPLICIT NONE … … 80 82 REAL(wp) :: zs21_11, zs21_12, zs21_21, zs21_22 81 83 REAL(wp) :: zs22_11, zs22_12, zs22_21, zs22_22 84 REAL(wp) :: zintb, zintn 82 85 REAL(wp), POINTER, DIMENSION(:,:) :: zfrld, zmass, zcorl 83 86 REAL(wp), POINTER, DIMENSION(:,:) :: za1ct, za2ct, zresr 84 87 REAL(wp), POINTER, DIMENSION(:,:) :: zc1u, zc1v, zc2u, zc2v 85 REAL(wp), POINTER, DIMENSION(:,:) :: zsang 88 REAL(wp), POINTER, DIMENSION(:,:) :: zsang, zpice 86 89 REAL(wp), POINTER, DIMENSION(:,:) :: zu0, zv0 87 90 REAL(wp), POINTER, DIMENSION(:,:) :: zu_n, zv_n … … 93 96 94 97 CALL wrk_alloc( jpi,jpj, zfrld, zmass, zcorl, za1ct, za2ct, zresr ) 95 CALL wrk_alloc( jpi,jpj, zc1u , zc1v , zc2u , zc2v , zsang )98 CALL wrk_alloc( jpi,jpj, zc1u , zc1v , zc2u , zc2v , zsang, zpice ) 96 99 CALL wrk_alloc( jpi,jpj+2, zu0, zv0, zu_n, zv_n, zu_a, zv_a, zviszeta, zviseta, kjstart = 0 ) 97 100 CALL wrk_alloc( jpi,jpj+2, zzfrld, zztms, zi1, zi2, zmasst, zpresh, kjstart = 0 ) … … 129 132 !i zviszeta(:,jpj+1) = 0._wp ; zviseta(:,jpj+1) = 0._wp 130 133 134 IF( nn_ice_embd == 2 ) THEN !== embedded sea ice: compute representative ice top surface ==! 135 ! 136 ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[n/nn_fsbc], n=0,nn_fsbc-1} 137 ! = (1/nn_fsbc)^2 * {SUM[n], n=0,nn_fsbc-1} 138 zintn = REAL( nn_fsbc - 1 ) / REAL( nn_fsbc ) * 0.5_wp 139 ! 140 ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[1-n/nn_fsbc], n=0,nn_fsbc-1} 141 ! = (1/nn_fsbc)^2 * (nn_fsbc^2 - {SUM[n], n=0,nn_fsbc-1}) 142 zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp 143 ! 144 zpice(:,:) = ssh_m(:,:) + ( zintn * snwice_mass(:,:) + zintb * snwice_mass_b(:,:) ) * r1_rau0 145 ! 146 ! 147 ELSE !== non-embedded sea ice: use ocean surface for slope calculation ==! 148 zpice(:,:) = ssh_m(:,:) 149 ENDIF 131 150 132 151 ! Ice mass, ice strength, and wind stress at the center | … … 196 215 197 216 ! Gradient of the sea surface height 198 zgsshx = ( ( ssh_m(ji ,jj ) - ssh_m(ji-1,jj ))/e1u(ji-1,jj ) &199 & + ( ssh_m(ji ,jj-1) - ssh_m(ji-1,jj-1))/e1u(ji-1,jj-1) ) * 0.5_wp200 zgsshy = ( ( ssh_m(ji ,jj ) - ssh_m(ji ,jj-1))/e2v(ji ,jj-1) &201 & + ( ssh_m(ji-1,jj ) - ssh_m(ji-1,jj-1))/e2v(ji-1,jj-1) ) * 0.5_wp217 zgsshx = ( (zpice(ji ,jj ) - zpice(ji-1,jj ))/e1u(ji-1,jj ) & 218 & + (zpice(ji ,jj-1) - zpice(ji-1,jj-1))/e1u(ji-1,jj-1) ) * 0.5_wp 219 zgsshy = ( (zpice(ji ,jj ) - zpice(ji ,jj-1))/e2v(ji ,jj-1) & 220 & + (zpice(ji-1,jj ) - zpice(ji-1,jj-1))/e2v(ji-1,jj-1) ) * 0.5_wp 202 221 203 222 ! Computation of the velocity field taking into account the ice-ice interaction. … … 575 594 576 595 CALL wrk_dealloc( jpi,jpj, zfrld, zmass, zcorl, za1ct, za2ct, zresr ) 577 CALL wrk_dealloc( jpi,jpj, zc1u , zc1v , zc2u , zc2v , zsang )596 CALL wrk_dealloc( jpi,jpj, zc1u , zc1v , zc2u , zc2v , zsang, zpice ) 578 597 CALL wrk_dealloc( jpi,jpj+2, zu0, zv0, zu_n, zv_n, zu_a, zv_a, zviszeta, zviseta, kjstart = 0 ) 579 598 CALL wrk_dealloc( jpi,jpj+2, zzfrld, zztms, zi1, zi2, zmasst, zpresh, kjstart = 0 ) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90
r3294 r3625 9 9 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 10 10 !! - ! 2010-11 (G. Madec) ice-ocean stress computed at each ocean time-step 11 !! 4.0 ! 2011-01 (A. R. Porter, STFC Daresbury) dynamical allocation 11 !! 3.3.1 ! 2011-01 (A. R. Porter, STFC Daresbury) dynamical allocation 12 !! 3.5 ! 2012-11 ((G. Madec, Y. Aksenov, A. Coward) salt and heat fluxes associated with e-p 12 13 !!---------------------------------------------------------------------- 13 14 #if defined key_lim2 … … 28 29 USE sbc_oce ! surface boundary condition: ocean 29 30 USE sbccpl 30 31 USE cpl_oasis3, ONLY : lk_cpl 32 USE oce , ONLY : sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass 31 33 USE albedo ! albedo parameters 32 34 USE lbclnk ! ocean lateral boundary condition - MPP exchanges … … 37 39 USE iom ! I/O library 38 40 USE prtctl ! Print control 39 USE cpl_oasis3, ONLY : lk_cpl41 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 40 42 41 43 IMPLICIT NONE … … 88 90 !! - Update the fluxes provided to the ocean 89 91 !! 90 !! ** Outputs : - qsr : sea heat flux :solar91 !! - qns : sea heat flux : non solar92 !! - emp : freshwater budget: volumeflux93 !! - emps : freshwater budget: concentration/dillution92 !! ** Outputs : - qsr : sea heat flux : solar 93 !! - qns : sea heat flux : non solar (including heat content of the mass flux) 94 !! - emp : freshwater budget: mass flux 95 !! - sfx : freshwater budget: salt flux due to Freezing/Melting 94 96 !! - utau : sea surface i-stress (ocean referential) 95 97 !! - vtau : sea surface j-stress (ocean referential) … … 107 109 INTEGER :: ifvt, i1mfr, idfr, iflt ! - - 108 110 INTEGER :: ial, iadv, ifral, ifrdv ! - - 109 REAL(wp) :: zqsr, zqns, zfm ! local scalars 110 REAL(wp) :: zinda, zfons, zemp ! - - 111 REAL(wp) :: zqsr, zqns, zfmm ! local scalars 112 REAL(wp) :: zinda, zfsalt, zemp ! - - 113 REAL(wp) :: zemp_snw, zqhc, zcd ! - - 114 REAL(wp) :: zswitch ! - - 111 115 REAL(wp), POINTER, DIMENSION(:,:) :: zqnsoce ! 2D workspace 112 116 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace … … 115 119 CALL wrk_alloc( jpi, jpj, zqnsoce ) 116 120 CALL wrk_alloc( jpi, jpj, 1, zalb, zalbp ) 121 122 SELECT CASE( nn_ice_embd ) ! levitating or embedded sea-ice option 123 CASE( 0 ) ; zswitch = 1 ! (0) standard levitating sea-ice : salt exchange only 124 CASE( 1, 2 ) ; zswitch = 0 ! (1) levitating sea-ice: salt and volume exchange but no pressure effect 125 ! (2) embedded sea-ice : salt and volume fluxes and pressure 126 END SELECT ! 117 127 118 128 !------------------------------------------! … … 133 143 ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv 134 144 135 !!$ zinda = 1.0 - AINT( pfrld(ji,jj) ) ! = 0. if pure ocean else 1. (at previous time)136 !!$ 137 !!$ i1mfr = 1.0 - AINT( frld(ji,jj) ) ! = 0. if pure ocean else 1. (at current time)138 !!$ 139 !!$ IF( phicif(ji,jj) <= 0. ) THEN ; ifvt = zinda ! = 1. if (snow and no ice at previous time) else 0.???140 !!$ ELSE ; ifvt = 0. 145 !!$ attempt to explain the tricky flags set above.... 146 !!$ zinda = 1.0 - AINT( pfrld(ji,jj) ) ! = 0. if ice-free ocean else 1. (after ice adv, but before ice thermo) 147 !!$ i1mfr = 1.0 - AINT( frld(ji,jj) ) ! = 0. if ice-free ocean else 1. (after ice thermo) 148 !!$ 149 !!$ IF( phicif(ji,jj) <= 0. ) THEN ; ifvt = zinda ! = zinda if previous thermodynamic step overmelted the ice??? 150 !!$ ELSE ; ifvt = 0. ! 141 151 !!$ ENDIF 142 152 !!$ 143 !!$ IF( frld(ji,jj) >= pfrld(ji,jj) ) THEN ; idfr = 0. ! = 0. if lead fraction increases from previous to current153 !!$ IF( frld(ji,jj) >= pfrld(ji,jj) ) THEN ; idfr = 0. ! = 0. if lead fraction increases due to ice thermodynamics 144 154 !!$ ELSE ; idfr = 1. 145 155 !!$ ENDIF 146 156 !!$ 147 !!$ iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ! = 1. if ice (not only snow) at previous and pure ocean at current157 !!$ iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ! = 1. if ice (not only snow) at previous time and ice-free ocean currently 148 158 !!$ 149 159 !!$ ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr 160 !!$ = i1mfr if ifvt = 1 i.e. 161 !!$ = idfr if ifvt = 0 150 162 !!$! snow no ice ice ice or nothing lead fraction increases 151 163 !!$! at previous now at previous 152 !!$! -> ice a era increases ??? -> ice aera decreases ???164 !!$! -> ice area increases ??? -> ice area decreases ??? 153 165 !!$ 154 166 !!$ iadv = ( 1 - i1mfr ) * zinda … … 174 186 #endif 175 187 ! computation the non solar heat flux at ocean surface 176 zqns = - ( 1. - thcm(ji,jj) ) * zqsr & ! part of the solar energy used in leads 177 & + iflt * ( fscmbq(ji,jj) + ffltbif(ji,jj) ) & 178 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice & 179 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice 180 181 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! ??? 182 ! 188 zqns = - ( 1. - thcm(ji,jj) ) * zqsr & ! part of the solar energy used in leads 189 & + iflt * ( fscmbq(ji,jj) + ffltbif(ji,jj) ) & 190 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice & 191 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice 192 193 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! store residual heat flux (to put into the ocean at the next time-step) 194 zqhc = ( rdq_snw(ji,jj) & 195 & + rdq_ice(ji,jj) * ( 1.- zswitch) ) * r1_rdtice ! heat flux due to snow ( & ice heat content, 196 ! ! if ice/ocean mass exchange active) 183 197 qsr (ji,jj) = zqsr ! solar heat flux 184 qns (ji,jj) = zqns - fdtcn(ji,jj) ! non solar heat flux 198 qns (ji,jj) = zqns - fdtcn(ji,jj) + zqhc ! non solar heat flux 199 ! 200 ! !------------------------------------------! 201 ! ! mass and salt flux at the ocean surface ! 202 ! !------------------------------------------! 203 ! 204 ! mass flux at the ocean-atmosphere interface (open ocean fraction = leads area) 205 #if defined key_coupled 206 ! ! coupled mode: 207 zemp = + emp_tot(ji,jj) & ! net mass flux over the grid cell (ice+ocean area) 208 & - emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) ! minus the mass flux intercepted by sea-ice 209 #else 210 ! ! forced mode: 211 zemp = + emp(ji,jj) * frld(ji,jj) & ! mass flux over open ocean fraction 212 & - tprecip(ji,jj) * ( 1. - frld(ji,jj) ) & ! liquid precip. over ice reaches directly the ocean 213 & + sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) ! snow is intercepted by sea-ice (previous frld) 214 #endif 215 ! 216 ! mass flux at the ocean/ice interface (sea ice fraction) 217 zemp_snw = rdm_snw(ji,jj) * r1_rdtice ! snow melting = pure water that enters the ocean 218 zfmm = rdm_ice(ji,jj) * r1_rdtice ! Freezing minus Melting (F-M) 219 220 ! salt flux at the ice/ocean interface (sea ice fraction) [PSU*kg/m2/s] 221 zfsalt = - sice_0(ji,jj) * zfmm ! F-M salt exchange 222 zcd = soce_0(ji,jj) * zfmm ! concentration/dilution term due to F-M 223 ! 224 ! salt flux only : add concentration dilution term in salt flux and no F-M term in volume flux 225 ! salt and mass fluxes : non concentration dilution term in salt flux and add F-M term in volume flux 226 sfx (ji,jj) = zfsalt + zswitch * zcd ! salt flux (+ C/D if no ice/ocean mass exchange) 227 emp (ji,jj) = zemp + zemp_snw + ( 1.- zswitch) * zfmm ! mass flux (+ F/M mass flux if ice/ocean mass exchange) 228 ! 185 229 END DO 186 230 END DO 231 ! !------------------------------------------! 232 ! ! mass of snow and ice per unit area ! 233 ! !------------------------------------------! 234 IF( nn_ice_embd /= 0 ) THEN ! embedded sea-ice (mass required) 235 snwice_mass_b(:,:) = snwice_mass(:,:) ! save mass from the previous ice time step 236 ! ! new mass per unit area 237 snwice_mass (:,:) = tms(:,:) * ( rhosn * hsnif(:,:) + rhoic * hicif(:,:) ) * ( 1.0 - frld(:,:) ) 238 ! ! time evolution of snow+ice mass 239 snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) / rdt_ice 240 ENDIF 187 241 188 242 CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) ) … … 190 244 CALL iom_put( 'qsr_io_cea', fstric(:,:) * (1.e0 - pfrld(:,:)) ) 191 245 192 !------------------------------------------!193 ! mass flux at the ocean surface !194 !------------------------------------------!195 DO jj = 1, jpj196 DO ji = 1, jpi197 !198 #if defined key_coupled199 ! freshwater exchanges at the ice-atmosphere / ocean interface (coupled mode)200 zemp = emp_tot(ji,jj) - emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) & !201 & + rdmsnif(ji,jj) * r1_rdtice ! freshwaterflux due to snow melting202 #else203 ! computing freshwater exchanges at the ice/ocean interface204 zemp = + emp(ji,jj) * frld(ji,jj) & ! e-p budget over open ocean fraction205 & - tprecip(ji,jj) * ( 1. - frld(ji,jj) ) & ! liquid precipitation reaches directly the ocean206 & + sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! change in ice cover within the time step207 & + rdmsnif(ji,jj) * r1_rdtice ! freshwater flux due to snow melting208 #endif209 !210 ! computing salt exchanges at the ice/ocean interface211 zfons = ( soce_0(ji,jj) - sice_0(ji,jj) ) * ( rdmicif(ji,jj) * r1_rdtice )212 !213 ! converting the salt flux from ice to a freshwater flux from ocean214 zfm = zfons / ( sss_m(ji,jj) + epsi16 )215 !216 emps(ji,jj) = zemp + zfm ! surface ocean concentration/dilution effect (use on SSS evolution)217 emp (ji,jj) = zemp ! surface ocean volume flux (use on sea-surface height evolution)218 !219 END DO220 END DO221 222 246 IF( lk_diaar5 ) THEN ! AR5 diagnostics 223 CALL iom_put( 'isnwmlt_cea' , rdm snif(:,:) * r1_rdtice )224 CALL iom_put( 'fsal_virt_cea', soce_0(:,:) * rdm icif(:,:) * r1_rdtice )225 CALL iom_put( 'fsal_real_cea', - sice_0(:,:) * rdm icif(:,:) * r1_rdtice )247 CALL iom_put( 'isnwmlt_cea' , rdm_snw(:,:) * r1_rdtice ) 248 CALL iom_put( 'fsal_virt_cea', soce_0(:,:) * rdm_ice(:,:) * r1_rdtice ) 249 CALL iom_put( 'fsal_real_cea', - sice_0(:,:) * rdm_ice(:,:) * r1_rdtice ) 226 250 ENDIF 227 251 … … 243 267 IF(ln_ctl) THEN ! control print 244 268 CALL prt_ctl(tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ') 245 CALL prt_ctl(tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2= emps , clinfo2=' emps: ')269 CALL prt_ctl(tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=sfx , clinfo2=' sfx : ') 246 270 CALL prt_ctl(tab2d_1=utau , clinfo1=' lim_sbc: utau : ', mask1=umask, & 247 271 & tab2d_2=vtau , clinfo2=' vtau : ' , mask2=vmask ) … … 439 463 END WHERE 440 464 ENDIF 465 ! ! embedded sea ice 466 IF( nn_ice_embd /= 0 ) THEN ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass 467 snwice_mass (:,:) = tms(:,:) * ( rhosn * hsnif(:,:) + rhoic * hicif(:,:) ) * ( 1.0 - frld(:,:) ) 468 snwice_mass_b(:,:) = snwice_mass(:,:) 469 ELSE 470 snwice_mass (:,:) = 0.e0 ! no mass exchanges 471 snwice_mass_b(:,:) = 0.e0 ! no mass exchanges 472 ENDIF 473 IF( nn_ice_embd == 2 .AND. & ! full embedment (case 2) & no restart : 474 & .NOT.ln_rstart ) THEN ! deplete the initial ssh below sea-ice area 475 sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0 476 sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0 477 ENDIF 441 478 ! 442 479 END SUBROUTINE lim_sbc_init_2 -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limthd_2.F90
r3294 r3625 13 13 !! 'key_lim2' : LIM 2.0 sea-ice model 14 14 !!---------------------------------------------------------------------- 15 !! lim_thd_2 : thermodynamic of sea ice16 !! lim_thd_init_2 : initialisation of sea-ice thermodynamic15 !! lim_thd_2 : thermodynamic of sea ice 16 !! lim_thd_init_2 : initialisation of sea-ice thermodynamic 17 17 !!---------------------------------------------------------------------- 18 USE phycst ! physical constants19 USE dom_oce ! ocean space and time domain variables18 USE phycst ! physical constants 19 USE dom_oce ! ocean space and time domain variables 20 20 USE domvvl 21 21 USE lbclnk 22 USE in_out_manager ! I/O manager22 USE in_out_manager ! I/O manager 23 23 USE lib_mpp 24 USE wrk_nemo ! work arrays25 USE iom ! IOM library26 USE ice_2 ! LIM sea-ice variables27 USE sbc_oce !28 USE sbc_ice !29 USE thd_ice_2 ! LIM thermodynamic sea-ice variables30 USE dom_ice_2 ! LIM sea-ice domain24 USE wrk_nemo ! work arrays 25 USE iom ! IOM library 26 USE ice_2 ! LIM sea-ice variables 27 USE sbc_oce ! 28 USE sbc_ice ! 29 USE thd_ice_2 ! LIM thermodynamic sea-ice variables 30 USE dom_ice_2 ! LIM sea-ice domain 31 31 USE limthd_zdf_2 32 32 USE limthd_lac_2 33 33 USE limtab_2 34 USE prtctl ! Print control 35 USE cpl_oasis3, ONLY : lk_cpl 36 USE diaar5, ONLY : lk_diaar5 37 34 USE prtctl ! Print control 35 USE cpl_oasis3, ONLY : lk_cpl 36 USE diaar5 , ONLY : lk_diaar5 37 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 38 38 39 IMPLICIT NONE 39 40 PRIVATE … … 55 56 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 56 57 !!---------------------------------------------------------------------- 57 58 58 CONTAINS 59 59 … … 89 89 REAL(wp) :: za , zh, zthsnice ! 90 90 REAL(wp) :: zfric_u ! friction velocity 91 REAL(wp) :: zfnsol ! total non solar heat92 REAL(wp) :: zfontn ! heat flux from snow thickness93 91 REAL(wp) :: zfntlat, zpareff ! test. the val. of lead heat budget 94 92 … … 129 127 zdvolif(:,:) = 0.e0 ! total variation of ice volume 130 128 zdvonif(:,:) = 0.e0 ! transformation of snow to sea-ice volume 131 ! zdvonif(:,:) = 0.e0 ! lateral variation of ice volume132 129 zlicegr(:,:) = 0.e0 ! lateral variation of ice volume 133 130 zdvomif(:,:) = 0.e0 ! variation of ice volume at bottom due to melting only … … 137 134 ffltbif(:,:) = 0.e0 ! linked with fstric 138 135 qfvbq (:,:) = 0.e0 ! linked with fstric 139 rdmsnif(:,:) = 0.e0 ! variation of snow mass per unit area 140 rdmicif(:,:) = 0.e0 ! variation of ice mass per unit area 136 rdm_snw(:,:) = 0.e0 ! variation of snow mass over 1 time step 137 rdq_snw(:,:) = 0.e0 ! heat content associated with rdm_snw 138 rdm_ice(:,:) = 0.e0 ! variation of ice mass over 1 time step 139 rdq_ice(:,:) = 0.e0 ! heat content associated with rdm_ice 141 140 zmsk (:,:,:) = 0.e0 142 141 … … 199 198 !-------------------------------------------------------------------------- 200 199 201 sst_m(:,:) = sst_m(:,:) + rt0 202 203 !CDIR NOVERRCHK 204 DO jj = 1, jpj 205 !CDIR NOVERRCHK 200 !CDIR NOVERRCHK 201 DO jj = 1, jpj 202 !CDIR NOVERRCHK 206 203 DO ji = 1, jpi 207 204 zthsnice = hsnif(ji,jj) + hicif(ji,jj) … … 217 214 ! temperature and turbulent mixing (McPhee, 1992) 218 215 zfric_u = MAX ( MIN( SQRT( ust2s(ji,jj) ) , zfric_umax ) , zfric_umin ) ! friction velocity 219 fdtcn(ji,jj) = zindb * rau0 * rcp * 0.006 * zfric_u * ( sst_m(ji,jj) - tfu(ji,jj) )216 fdtcn(ji,jj) = zindb * rau0 * rcp * 0.006 * zfric_u * ( sst_m(ji,jj) + rt0 - tfu(ji,jj) ) 220 217 qdtcn(ji,jj) = zindb * fdtcn(ji,jj) * frld(ji,jj) * rdt_ice 221 218 222 219 ! partial computation of the lead energy budget (qldif) 223 220 #if defined key_coupled 224 qldif(ji,jj) = tms(ji,jj) * rdt_ice &221 qldif(ji,jj) = tms(ji,jj) * rdt_ice & 225 222 & * ( ( qsr_tot(ji,jj) - qsr_ice(ji,jj,1) * zfricp ) * ( 1.0 - thcm(ji,jj) ) & 226 223 & + ( qns_tot(ji,jj) - qns_ice(ji,jj,1) * zfricp ) & 227 224 & + frld(ji,jj) * ( fdtcn(ji,jj) + ( 1.0 - zindb ) * fsbbq(ji,jj) ) ) 228 225 #else 229 zfontn = ( sprecip(ji,jj) / rhosn ) * xlsn ! energy for melting solid precipitation 230 zfnsol = qns(ji,jj) ! total non solar flux over the ocean 231 qldif(ji,jj) = tms(ji,jj) * ( qsr(ji,jj) * ( 1.0 - thcm(ji,jj) ) & 232 & + zfnsol + fdtcn(ji,jj) - zfontn & 233 & + ( 1.0 - zindb ) * fsbbq(ji,jj) ) & 234 & * frld(ji,jj) * rdt_ice 235 !!$ qldif(ji,jj) = tms(ji,jj) * rdt_ice * frld(ji,jj) 236 !!$ & * ( qsr(ji,jj) * ( 1.0 - thcm(ji,jj) ) & 237 !!$ & + qns(ji,jj) + fdtcn(ji,jj) - zfontn & 238 !!$ & + ( 1.0 - zindb ) * fsbbq(ji,jj) ) & 226 qldif(ji,jj) = tms(ji,jj) * rdt_ice * frld(ji,jj) & 227 & * ( qsr(ji,jj) * ( 1.0 - thcm(ji,jj) ) & 228 & + qns(ji,jj) + fdtcn(ji,jj) & 229 & + ( 1.0 - zindb ) * fsbbq(ji,jj) ) 239 230 #endif 240 231 ! parlat : percentage of energy used for lateral ablation (0.0) … … 246 237 247 238 ! energy needed to bring ocean surface layer until its freezing 248 qcmif (ji,jj) = rau0 * rcp * fse3t_m(ji,jj,1) & 249 & * ( tfu(ji,jj) - sst_m(ji,jj) ) * ( 1 - zinda ) 239 qcmif (ji,jj) = rau0 * rcp * fse3t_m(ji,jj,1) * ( tfu(ji,jj) - sst_m(ji,jj) - rt0 ) * ( 1 - zinda ) 250 240 251 241 ! calculate oceanic heat flux. … … 257 247 END DO 258 248 259 sst_m(:,:) = sst_m(:,:) - rt0260 261 249 ! Select icy points and fulfill arrays for the vectorial grid. 262 250 !---------------------------------------------------------------------- … … 312 300 CALL tab_2d_1d_2( nbpb, qldif_1d (1:nbpb) , qldif , jpi, jpj, npb(1:nbpb) ) 313 301 CALL tab_2d_1d_2( nbpb, qstbif_1d (1:nbpb) , qstoif , jpi, jpj, npb(1:nbpb) ) 314 CALL tab_2d_1d_2( nbpb, rdmicif_1d (1:nbpb) , rdmicif , jpi, jpj, npb(1:nbpb) ) 302 CALL tab_2d_1d_2( nbpb, rdm_ice_1d (1:nbpb) , rdm_ice , jpi, jpj, npb(1:nbpb) ) 303 CALL tab_2d_1d_2( nbpb, rdq_ice_1d (1:nbpb) , rdq_ice , jpi, jpj, npb(1:nbpb) ) 315 304 CALL tab_2d_1d_2( nbpb, dmgwi_1d (1:nbpb) , dmgwi , jpi, jpj, npb(1:nbpb) ) 305 CALL tab_2d_1d_2( nbpb, rdm_snw_1d (1:nbpb) , rdm_snw , jpi, jpj, npb(1:nbpb) ) 306 CALL tab_2d_1d_2( nbpb, rdq_snw_1d (1:nbpb) , rdq_snw , jpi, jpj, npb(1:nbpb) ) 316 307 CALL tab_2d_1d_2( nbpb, qlbbq_1d (1:nbpb) , zqlbsbq , jpi, jpj, npb(1:nbpb) ) 317 308 ! … … 332 323 CALL tab_1d_2d_2( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb) , jpi, jpj ) 333 324 CALL tab_1d_2d_2( nbpb, qstoif , npb, qstbif_1d (1:nbpb) , jpi, jpj ) 334 CALL tab_1d_2d_2( nbpb, rdmicif , npb, rdmicif_1d(1:nbpb) , jpi, jpj ) 325 CALL tab_1d_2d_2( nbpb, rdm_ice , npb, rdm_ice_1d(1:nbpb) , jpi, jpj ) 326 CALL tab_1d_2d_2( nbpb, rdq_ice , npb, rdq_ice_1d(1:nbpb) , jpi, jpj ) 335 327 CALL tab_1d_2d_2( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb) , jpi, jpj ) 336 CALL tab_1d_2d_2( nbpb, rdmsnif , npb, rdmsnif_1d(1:nbpb) , jpi, jpj ) 328 CALL tab_1d_2d_2( nbpb, rdm_snw , npb, rdm_snw_1d(1:nbpb) , jpi, jpj ) 329 CALL tab_1d_2d_2( nbpb, rdq_snw , npb, rdq_snw_1d(1:nbpb) , jpi, jpj ) 337 330 CALL tab_1d_2d_2( nbpb, zdvosif , npb, dvsbq_1d (1:nbpb) , jpi, jpj ) 338 331 CALL tab_1d_2d_2( nbpb, zdvobif , npb, dvbbq_1d (1:nbpb) , jpi, jpj ) … … 393 386 IF( nbpac > 0 ) THEN 394 387 ! 395 zlicegr(:,:) = rdm icif(:,:) ! to output the lateral sea-ice growth388 zlicegr(:,:) = rdm_ice(:,:) ! to output the lateral sea-ice growth 396 389 !...Put the variable in a 1-D array for lateral accretion 397 390 CALL tab_2d_1d_2( nbpac, frld_1d (1:nbpac) , frld , jpi, jpj, npac(1:nbpac) ) … … 404 397 CALL tab_2d_1d_2( nbpac, qcmif_1d (1:nbpac) , qcmif , jpi, jpj, npac(1:nbpac) ) 405 398 CALL tab_2d_1d_2( nbpac, qstbif_1d (1:nbpac) , qstoif , jpi, jpj, npac(1:nbpac) ) 406 CALL tab_2d_1d_2( nbpac, rdmicif_1d(1:nbpac) , rdmicif , jpi, jpj, npac(1:nbpac) ) 399 CALL tab_2d_1d_2( nbpac, rdm_ice_1d(1:nbpac) , rdm_ice , jpi, jpj, npac(1:nbpac) ) 400 CALL tab_2d_1d_2( nbpac, rdq_ice_1d(1:nbpac) , rdq_ice , jpi, jpj, npac(1:nbpac) ) 407 401 CALL tab_2d_1d_2( nbpac, dvlbq_1d (1:nbpac) , zdvolif , jpi, jpj, npac(1:nbpac) ) 408 402 CALL tab_2d_1d_2( nbpac, tfu_1d (1:nbpac) , tfu , jpi, jpj, npac(1:nbpac) ) … … 418 412 CALL tab_1d_2d_2( nbpac, tbif(:,:,3), npac(1:nbpac), tbif_1d (1:nbpac , 3 ), jpi, jpj ) 419 413 CALL tab_1d_2d_2( nbpac, qstoif , npac(1:nbpac), qstbif_1d (1:nbpac) , jpi, jpj ) 420 CALL tab_1d_2d_2( nbpac, rdmicif , npac(1:nbpac), rdmicif_1d(1:nbpac) , jpi, jpj ) 414 CALL tab_1d_2d_2( nbpac, rdm_ice , npac(1:nbpac), rdm_ice_1d(1:nbpac) , jpi, jpj ) 415 CALL tab_1d_2d_2( nbpac, rdq_ice , npac(1:nbpac), rdq_ice_1d(1:nbpac) , jpi, jpj ) 421 416 CALL tab_1d_2d_2( nbpac, zdvolif , npac(1:nbpac), dvlbq_1d (1:nbpac) , jpi, jpj ) 422 417 ! … … 449 444 CALL iom_put( 'iceprod_cea' , hicifp (:,:) * zztmp ) ! Ice produced [m/s] 450 445 IF( lk_diaar5 ) THEN 451 CALL iom_put( 'snowmel_cea' , rdm snif(:,:) * zztmp ) ! Snow melt [kg/m2/s]446 CALL iom_put( 'snowmel_cea' , rdm_snw(:,:) * zztmp ) ! Snow melt [kg/m2/s] 452 447 zztmp = rhoic / rdt_ice 453 448 CALL iom_put( 'sntoice_cea' , zdvonif(:,:) * zztmp ) ! Snow to Ice transformation [kg/m2/s] 454 449 CALL iom_put( 'ticemel_cea' , zdvosif(:,:) * zztmp ) ! Melt at Sea Ice top [kg/m2/s] 455 450 CALL iom_put( 'bicemel_cea' , zdvomif(:,:) * zztmp ) ! Melt at Sea Ice bottom [kg/m2/s] 456 zlicegr(:,:) = MAX( 0.e0, rdm icif(:,:)-zlicegr(:,:) )457 CALL iom_put( 'licepro_cea' , zlicegr(:,:) * zztmp ) ! Later eal sea ice growth[kg/m2/s]451 zlicegr(:,:) = MAX( 0.e0, rdm_ice(:,:)-zlicegr(:,:) ) 452 CALL iom_put( 'licepro_cea' , zlicegr(:,:) * zztmp ) ! Lateral sea ice growth [kg/m2/s] 458 453 ENDIF 459 454 ! -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limthd_lac_2.F90
r3294 r3625 7 7 8 8 !!---------------------------------------------------------------------- 9 !! lim_lat_acr_2 10 !!---------------------------------------------------------------------- 11 USE par_oce 9 !! lim_lat_acr_2 : lateral accretion of ice 10 !!---------------------------------------------------------------------- 11 USE par_oce ! ocean parameters 12 12 USE phycst 13 13 USE thd_ice_2 14 14 USE ice_2 15 15 USE limistate_2 16 USE lib_mpp ! MPP library 17 USE wrk_nemo ! work arrays 16 USE lib_mpp ! MPP library 17 USE wrk_nemo ! work arrays 18 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 18 19 19 20 IMPLICIT NONE … … 145 146 frld_1d (ji) = MAX( zfrlnew , zfrlmin(ji) ) 146 147 !--computation of the remaining part of ice thickness which has been already used 147 zdhicbot(ji) = ( frld_1d(ji) - zfrlnew ) * zhice0(ji) / ( 1.0 - zfrlmin(ji) ) &148 148 zdhicbot(ji) = ( frld_1d(ji) - zfrlnew ) * zhice0(ji) / ( 1.0 - zfrlmin(ji) ) & 149 & - ( ( 1.0 - zfrrate ) / ( 1.0 - frld_1d(ji) ) ) * ( zqbgow(ji) / xlic ) 149 150 END DO 150 151 … … 196 197 & ) / zah 197 198 198 tbif_1d(ji,3) = ( iiceform * ( zhnews2 - zdh3 )* zta1 &199 tbif_1d(ji,3) = ( iiceform * ( zhnews2 - zdh3 ) * zta1 & 199 200 & + ( iiceform * zdh3 + ( 1 - iiceform ) * zdh1 ) * zta2 & 200 201 & + ( iiceform * ( zhnews2 - zdh5 ) + ( 1 - iiceform ) * ( zhnews2 - zdh1 ) ) * zta3 & … … 217 218 DO ji = kideb , kiut 218 219 dvlbq_1d (ji) = ( 1. - frld_1d(ji) ) * h_ice_1d(ji) - ( 1. - zfrl_old(ji) ) * zhice_old(ji) 219 rdmicif_1d(ji) = rdmicif_1d(ji) + rhoic * dvlbq_1d(ji) 220 rdm_ice_1d(ji) = rdm_ice_1d(ji) + rhoic * dvlbq_1d(ji) 221 rdq_ice_1d(ji) = rdq_ice_1d(ji) + rcpic * dvlbq_1d(ji) * ( tfu_1d(ji) - rt0 ) ! heat content relative to rt0 220 222 END DO 221 223 -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limthd_zdf_2.F90
r3294 r3625 18 18 USE ice_2 19 19 USE limistate_2 20 USE cpl_oasis3, ONLY : lk_cpl 20 21 USE in_out_manager 21 22 USE lib_mpp ! MPP library 22 23 USE wrk_nemo ! work arrays 23 USE cpl_oasis3, ONLY : lk_cpl24 24 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 25 25 26 IMPLICIT NONE 26 27 PRIVATE … … 86 87 REAL(wp), POINTER, DIMENSION(:) :: zrcpdt ! h_su*rho_su*cp_su/dt(h_su being the thick. of surf. layer) 87 88 REAL(wp), POINTER, DIMENSION(:) :: zts_old ! previous surface temperature 88 REAL(wp), POINTER, DIMENSION(:) :: zidsn , z1midsn , zidsnic ! tempor y variables89 REAL(wp), POINTER, DIMENSION(:) :: zidsn , z1midsn , zidsnic ! temporary variables 89 90 REAL(wp), POINTER, DIMENSION(:) :: zfnet ! net heat flux at the top surface( incl. conductive heat flux) 90 91 REAL(wp), POINTER, DIMENSION(:) :: zsprecip ! snow accumulation … … 98 99 REAL(wp), POINTER, DIMENSION(:) :: zep ! internal temperature of the 2nd layer of the snow/ice system 99 100 REAL(wp), DIMENSION(3) :: & 100 zplediag & ! principle diagonal, subdiag. and supdiag. of the101 zplediag & ! principle diagonal, subdiag. and supdiag. of the 101 102 , zsubdiag & ! tri-diagonal matrix coming from the computation 102 103 , zsupdiag & ! of the temperatures inside the snow-ice system 103 104 , zsmbr ! second member 104 REAL(wp) :: & 105 zhsu & ! thickness of surface layer 106 , zhe & ! effective thickness for compu. of equ. thermal conductivity 107 , zheshth & ! = zhe / thth 108 , zghe & ! correction factor of the thermal conductivity 109 , zumsb & ! parameter for numerical method to solve heat-diffusion eq. 110 , zkhsn & ! conductivity at the snow layer 111 , zkhic & ! conductivity at the ice layers 112 , zkint & ! equivalent conductivity at the snow-ice interface 113 , zkhsnint & ! = zkint*dt / (hsn*rhosn*cpsn) 114 , zkhicint & ! = 2*zkint*dt / (hic*rhoic*cpic) 115 , zpiv1 , zpiv2 & ! tempory scalars used to solve the tri-diagonal system 116 , zb2 , zd2 , zb3 , zd3 & 105 REAL(wp) :: & 106 zhsu & ! thickness of surface layer 107 , zhe & ! effective thickness for compu. of equ. thermal conductivity 108 , zheshth & ! = zhe / thth 109 , zghe & ! correction factor of the thermal conductivity 110 , zumsb & ! parameter for numerical method to solve heat-diffusion eq. 111 , zkhsn & ! conductivity at the snow layer 112 , zkhic & ! conductivity at the ice layers 113 , zkint & ! equivalent conductivity at the snow-ice interface 114 , zkhsnint & ! = zkint*dt / (hsn*rhosn*cpsn) 115 , zkhicint & ! = 2*zkint*dt / (hic*rhoic*cpic) 116 , zpiv1, zpiv2 & ! temporary scalars used to solve the tri-diagonal system 117 , zb2, zd2 & ! temporary scalars used to solve the tri-diagonal system 118 , zb3, zd3 & ! temporary scalars used to solve the tri-diagonal system 117 119 , ztint ! equivalent temperature at the snow-ice interface 118 REAL(wp) :: &119 zexp &! exponential function of the ice thickness120 , zfsab & 121 , zfts & 122 , zdfts & 123 , zdts & 124 , zqsnw_mlt & 125 , zdhsmlt & 126 , zhsn & 127 , zqsn_mlt_rem & 128 , zqice_top_mlt & 129 , zdhssub &! change in snow thick. due to sublimation or evaporation130 , zdhisub &! change in ice thick. due to sublimation or evaporation131 , zdhsn &! snow ice thickness increment132 , zdtsn &! snow internal temp. increment133 , zdtic &! ice internal temp. increment120 REAL(wp) :: & 121 zexp & ! exponential function of the ice thickness 122 , zfsab & ! part of solar radiation stored in brine pockets 123 , zfts & ! value of energy balance function when the temp. equal surf. temp. 124 , zdfts & ! value of derivative of ztfs when the temp. equal surf. temp. 125 , zdts & ! surface temperature increment 126 , zqsnw_mlt & ! energy needed to melt snow 127 , zdhsmlt & ! change in snow thickness due to melt 128 , zhsn & ! snow thickness (previous+accumulation-melt) 129 , zqsn_mlt_rem & ! remaining heat coming from snow melting 130 , zqice_top_mlt &! energy used to melt ice at top surface 131 , zdhssub & ! change in snow thick. due to sublimation or evaporation 132 , zdhisub & ! change in ice thick. due to sublimation or evaporation 133 , zdhsn & ! snow ice thickness increment 134 , zdtsn & ! snow internal temp. increment 135 , zdtic & ! ice internal temp. increment 134 136 , zqnes ! conductive energy due to ice melting in the first ice layer 135 REAL(wp) :: & 136 ztbot & ! temperature at the bottom surface 137 , zfcbot & ! conductive heat flux at bottom surface 138 , zqice_bot & ! energy used for bottom melting/growing 139 , zqice_bot_mlt & ! energy used for bottom melting 140 , zqstbif_bot & ! part of energy stored in brine pockets used for bottom melting 141 , zqstbif_old & ! tempory var. for zqstbif_bot 142 , zdhicmlt & ! change in ice thickness due to bottom melting 143 , zdhicm & ! change in ice thickness var. 144 , zdhsnm & ! change in snow thickness var. 145 , zhsnfi & ! snow thickness var. 146 , zc1, zpc1, zc2, zpc2, zp1, zp2 & ! tempory variables 147 , ztb2, ztb3 148 REAL(wp) :: & 149 zdrmh & ! change in snow/ice thick. after snow-ice formation 150 , zhicnew & ! new ice thickness 151 , zhsnnew & ! new snow thickness 152 , zquot , ztneq & ! tempory temp. variables 153 , zqice, zqicetot & ! total heat inside the snow/ice system 154 , zdfrl & ! change in ice concentration 155 , zdvsnvol & ! change in snow volume 156 , zdrfrl1, zdrfrl2 & ! tempory scalars 157 , zihsn, zidhb, zihic, zihe, zihq, ziexp, ziqf, zihnf, zibmlt, ziqr, zihgnew, zind 137 REAL(wp) :: & 138 ztbot & ! temperature at the bottom surface 139 , zfcbot & ! conductive heat flux at bottom surface 140 , zqice_bot & ! energy used for bottom melting/growing 141 , zqice_bot_mlt &! energy used for bottom melting 142 , zqstbif_bot & ! part of energy stored in brine pockets used for bottom melting 143 , zqstbif_old & ! temporary var. for zqstbif_bot 144 , zdhicmlt & ! change in ice thickness due to bottom melting 145 , zdhicm & ! change in ice thickness var. 146 , zdhsnm & ! change in snow thickness var. 147 , zhsnfi & ! snow thickness var. 148 , zc1, zpc1 & ! temporary variables 149 , zc2, zpc2 & ! temporary variables 150 , zp1, zp2 & ! temporary variables 151 , ztb2, ztb3 ! temporary variables 152 REAL(wp) :: & 153 zdrmh & ! change in snow/ice thick. after snow-ice formation 154 , zhicnew & ! new ice thickness 155 , zhsnnew & ! new snow thickness 156 , zquot & 157 , ztneq & ! temporary temp. variables 158 , zqice & 159 , zqicetot & ! total heat inside the snow/ice system 160 , zdfrl & ! change in ice concentration 161 , zdvsnvol & ! change in snow volume 162 , zdrfrl1, zdrfrl2, zihsn, zidhb, zihic & ! temporary scalars 163 , zihe, zihq, ziexp, ziqf, zihnf & ! temporary scalars 164 , zibmlt, ziqr, zihgnew, zind, ztmp ! temporary scalars 158 165 !!---------------------------------------------------------------------- 159 166 CALL wrk_alloc( jpij, ztsmlt, ztbif , zksn , zkic , zksndh , zfcsu , zfcsudt , zi0 , z1mi0 , zqmax ) … … 169 176 170 177 DO ji = kideb , kiut 178 ! do nothing if the snow (ice) thickness falls below its minimum thickness 171 179 zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 172 180 zihic = MAX( zzero , SIGN( zone , hicdif - h_ice_1d(ji) ) ) 173 !--computation of energy due to surface melting 174 zqcmlts(ji) = ( MAX ( zzero , & 175 & rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 176 !--computation of energy due to bottom melting 177 zqcmltb(ji) = ( MAX( zzero , & 178 & rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 179 & + MAX( zzero , & 180 & rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 181 & ) * ( 1.0 - zihic ) 182 !--limitation of snow/ice system internal temperature 181 !--energy required to bring snow to its melting point (rt0_snow) 182 zqcmlts(ji) = ( MAX ( zzero , rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 183 !--energy required to bring ice to its melting point (rt0_ice) 184 zqcmltb(ji) = ( MAX( zzero , rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 185 & + MAX( zzero , rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 186 & ) * ( 1.0 - zihic ) 187 !--limitation of snow/ice system internal temperature 183 188 tbif_1d(ji,1) = MIN( rt0_snow, tbif_1d(ji,1) ) 184 189 tbif_1d(ji,2) = MIN( rt0_ice , tbif_1d(ji,2) ) … … 480 485 dvsbq_1d(ji) = ( 1.0 - frld_1d(ji) ) * ( h_snow_1d(ji) - zhsnw_old(ji) - zsprecip(ji) ) 481 486 dvsbq_1d(ji) = MIN( zzero , dvsbq_1d(ji) ) 482 rdmsnif_1d(ji) = rhosn * dvsbq_1d(ji) 487 ztmp = rhosn * dvsbq_1d(ji) 488 rdm_snw_1d(ji) = ztmp 489 !--heat content of the water provided to the ocean (referenced to rt0) 490 rdq_snw_1d(ji) = cpic * ztmp * ( rt0_snow - rt0 ) 483 491 !-- If the snow is completely melted the remaining heat is used to melt ice 484 492 zqsn_mlt_rem = MAX( zzero , -zhsn ) * xlsn … … 623 631 !---updating new ice thickness and computing the newly formed ice mass 624 632 zhicnew = zihgnew * zhicnew 625 rdmicif_1d(ji) = rdmicif_1d(ji) + ( 1.0 - frld_1d(ji) ) * ( zhicnew - h_ice_1d(ji) ) * rhoic 633 ztmp = ( 1.0 - frld_1d(ji) ) * ( zhicnew - h_ice_1d(ji) ) * rhoic 634 rdm_ice_1d(ji) = rdm_ice_1d(ji) + ztmp 635 !---heat content of the water provided to the ocean (referenced to rt0) 636 ! use of rt0_ice is OK for melting ice; in the case of freezing, tfu_1d should be used. 637 ! This is done in 9.5 section (see below) 638 rdq_ice_1d(ji) = cpic * ztmp * ( rt0_ice - rt0 ) 626 639 !---updating new snow thickness and computing the newly formed snow mass 627 640 zhsnfi = zhsn + zdhsnm 628 641 h_snow_1d(ji) = MAX( zzero , zhsnfi ) 629 rdmsnif_1d(ji) = rdmsnif_1d(ji) + ( 1.0 - frld_1d(ji) ) * ( h_snow_1d(ji) - zhsn ) * rhosn 642 ztmp = ( 1.0 - frld_1d(ji) ) * ( h_snow_1d(ji) - zhsn ) * rhosn 643 rdm_snw_1d(ji) = rdm_snw_1d(ji) + ztmp 644 !---updating the heat content of the water provided to the ocean (referenced to rt0) 645 rdq_snw_1d(ji) = rdq_snw_1d(ji) + cpic * ztmp * ( rt0_snow - rt0 ) 630 646 !--remaining energy in case of total ablation 631 647 zqocea(ji) = - ( zihsn * xlic * zdhicm + xlsn * ( zhsnfi - h_snow_1d(ji) ) ) * ( 1.0 - frld_1d(ji) ) … … 659 675 tbif_1d(ji,3) = zihgnew * ztb3 + ( 1.0 - zihgnew ) * tfu_1d(ji) 660 676 h_ice_1d(ji) = zhicnew 677 ! update the ice heat content given to the ocean in freezing case 678 ! (part due to difference between rt0_ice and tfu_1d) 679 ztmp = ( 1. - zidhb ) * rhoic * dvbbq_1d(ji) 680 rdq_ice_1d(ji) = rdq_ice_1d(ji) + cpic * ztmp * ( tfu_1d(ji) - rt0_ice ) 661 681 END DO 662 682 … … 700 720 dmgwi_1d(ji) = dmgwi_1d(ji) + ( 1.0 -frld_1d(ji) ) * ( h_snow_1d(ji) - zhsnnew ) * rhosn 701 721 !--- volume change of ice and snow (used for ocean-ice freshwater flux computation) 702 rdmicif_1d(ji) = rdmicif_1d(ji) + ( 1.0 - frld_1d(ji) ) * ( zhicnew - h_ice_1d (ji) ) * rhoic 703 rdmsnif_1d(ji) = rdmsnif_1d(ji) + ( 1.0 - frld_1d(ji) ) * ( zhsnnew - h_snow_1d(ji) ) * rhosn 722 ztmp = ( 1.0 - frld_1d(ji) ) * ( zhicnew - h_ice_1d (ji) ) * rhoic 723 rdm_ice_1d(ji) = rdm_ice_1d(ji) + ztmp 724 rdq_ice_1d(ji) = rdq_ice_1d(ji) + cpic * ztmp * ( tfu_1d(ji) - rt0 ) 725 !!gm BUG ?? snow ==> only needed for nn_ice_embd == 0 (standard levitating sea-ice) 726 ztmp = ( 1.0 - frld_1d(ji) ) * ( zhsnnew - h_snow_1d(ji) ) * rhosn 727 rdm_snw_1d(ji) = rdm_snw_1d(ji) + ztmp 728 rdq_snw_1d(ji) = rdq_snw_1d(ji) + cpic * ztmp * ( rt0_snow - rt0 ) 704 729 705 730 !--- Actualize new snow and ice thickness. … … 748 773 !--variation of ice volume and ice mass 749 774 dvlbq_1d(ji) = zihic * ( zfrl_old(ji) - frld_1d(ji) ) * h_ice_1d(ji) 750 rdmicif_1d(ji) = rdmicif_1d(ji) + dvlbq_1d(ji) * rhoic 775 ztmp = dvlbq_1d(ji) * rhoic 776 rdm_ice_1d(ji) = rdm_ice_1d(ji) + ztmp 777 !!gm 778 !!gm This should be split in two parts: 779 !!gm 1- heat required to bring sea-ice to tfu : this part should be added to the heat flux taken from the ocean 780 !!gm cpic * ztmp * 0.5 * ( tbif_1d(ji,2) + tbif_1d(ji,3) - 2.* rt0_ice ) 781 !!gm 2- heat content of lateral ablation referenced to rt0 : this part only put in rdq_ice_1d 782 !!gm cpic * ztmp * ( rt0_ice - rt0 ) 783 !!gm Currently we put all the heat in rdq_ice_1d 784 rdq_ice_1d(ji) = rdq_ice_1d(ji) + cpic * ztmp * 0.5 * ( tbif_1d(ji,2) + tbif_1d(ji,3) - 2.* rt0 ) 785 ! 751 786 !--variation of snow volume and snow mass 752 zdvsnvol = zihsn * ( zfrl_old(ji) - frld_1d(ji) ) * h_snow_1d(ji) 753 rdmsnif_1d(ji) = rdmsnif_1d(ji) + zdvsnvol * rhosn 787 zdvsnvol = zihsn * ( zfrl_old(ji) - frld_1d(ji) ) * h_snow_1d(ji) 788 ztmp = zdvsnvol * rhosn 789 rdm_snw_1d(ji) = rdm_snw_1d(ji) + ztmp 790 !!gm 791 !!gm This should be split in two parts: 792 !!gm 1- heat required to bring snow to tfu : this part should be added to the heat flux taken from the ocean 793 !!gm cpic * ztmp * ( tbif_1d(ji,1) - rt0_snow ) 794 !!gm 2- heat content of lateral ablation referenced to rt0 : this part only put in rdq_snw_1d 795 !!gm cpic * ztmp * ( rt0_snow - rt0 ) 796 !!gm Currently we put all the heat in rdq_snw_1d 797 rdq_snw_1d(ji) = rdq_snw_1d(ji) + cpic * ztmp * ( tbif_1d(ji,1) - rt0 ) 798 754 799 h_snow_1d(ji) = ziqf * h_snow_1d(ji) 755 800 -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limwri_2.F90
r3294 r3625 13 13 !!---------------------------------------------------------------------- 14 14 !!---------------------------------------------------------------------- 15 !! lim_wri_2 : write of the diagnostics variables in ouput file16 !! lim_wri_init_2 : initialization and namelist read15 !! lim_wri_2 : write of the diagnostics variables in ouput file 16 !! lim_wri_init_2 : initialization and namelist read 17 17 !! lim_wri_state_2 : write for initial state or/and abandon: 18 18 !! > output.init.nc (if ninist = 1 in namelist) … … 26 26 USE ice_2 27 27 28 USE dianam ! build name of file (routine)28 USE dianam ! build name of file (routine) 29 29 USE lbclnk 30 30 USE in_out_manager 31 USE lib_mpp ! MPP library32 USE wrk_nemo ! work arrays31 USE lib_mpp ! MPP library 32 USE wrk_nemo ! work arrays 33 33 USE iom 34 34 USE ioipsl 35 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 35 36 36 37 IMPLICIT NONE … … 173 174 zcmo(ji,jj,13) = qns(ji,jj) 174 175 ! See thersf for the coefficient 175 zcmo(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce !!gm ???176 zcmo(ji,jj,14) = - sfx (ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce !!gm ??? 176 177 zcmo(ji,jj,15) = utau_ice(ji,jj) 177 178 zcmo(ji,jj,16) = vtau_ice(ji,jj) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/limwri_dimg_2.h90
r3294 r3625 118 118 zcmo(ji,jj,13) = qns(ji,jj) 119 119 ! See thersf for the coefficient 120 zcmo(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce120 zcmo(ji,jj,14) = - sfx (ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce 121 121 zcmo(ji,jj,15) = utau_ice(ji,jj) 122 122 zcmo(ji,jj,16) = vtau_ice(ji,jj) … … 161 161 rcmoy(ji,jj,13) = qns(ji,jj) 162 162 ! See thersf for the coefficient 163 rcmoy(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce163 rcmoy(ji,jj,14) = - sfx (ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce 164 164 rcmoy(ji,jj,15) = utau_ice(ji,jj) 165 165 rcmoy(ji,jj,16) = vtau_ice(ji,jj) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_2/thd_ice_2.F90
r2715 r3625 68 68 qstbif_1d , & !: " " qstoif 69 69 fbif_1d , & !: " " fbif 70 rdmicif_1d , & !: " " rdmicif 71 rdmsnif_1d , & !: " " rdmsnif 70 rdm_ice_1d , & !: " " rdm_ice 71 rdq_ice_1d , & !: " " rdq_ice 72 rdm_snw_1d , & !: " " rdm_snw 73 rdq_snw_1d , & !: " " rdq_snw 72 74 qlbbq_1d , & !: " " qlbsbq 73 75 dmgwi_1d , & !: " " dmgwi … … 108 110 & qstbif_1d(jpij), fbif_1d(jpij), Stat=ierr(2)) 109 111 ! 110 ALLOCATE( rdmicif_1d(jpij), rdmsnif_1d(jpij), qlbbq_1d(jpij), & 112 ALLOCATE( rdm_ice_1d(jpij), rdq_ice_1d(jpij) , & 113 & rdm_snw_1d(jpij), rdq_snw_1d(jpij), qlbbq_1d(jpij) , & 111 114 & dmgwi_1d(jpij) , dvsbq_1d(jpij) , rdvomif_1d(jpij), & 112 115 & dvbbq_1d(jpij) , dvlbq_1d(jpij) , dvnbq_1d(jpij) , & -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/dom_ice.F90
r2777 r3625 9 9 USE par_ice ! LIM-3 parameter 10 10 USE in_out_manager ! I/O manager 11 USE lib_mpp ! MPP library 11 USE lib_mpp ! MPP library 12 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 12 13 13 14 IMPLICIT NONE … … 30 31 31 32 !!---------------------------------------------------------------------- 32 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)33 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 33 34 !! $Id$ 34 35 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
r2777 r3625 9 9 #if defined key_lim3 10 10 !!---------------------------------------------------------------------- 11 !! 'key_lim3' : LIM3 sea-ice model12 !!---------------------------------------------------------------------- 13 USE par_ice 14 USE in_out_manager 15 USE lib_mpp 11 !! 'key_lim3' LIM-3 sea-ice model 12 !!---------------------------------------------------------------------- 13 USE par_ice ! LIM sea-ice parameters 14 USE in_out_manager ! I/O manager 15 USE lib_mpp ! MPP library 16 16 17 17 IMPLICIT NONE … … 158 158 !! * Share Module variables 159 159 !!-------------------------------------------------------------------------- 160 INTEGER , PUBLIC :: nstart !: iteration number of the begining of the run 161 INTEGER , PUBLIC :: nlast !: iteration number of the end of the run 162 INTEGER , PUBLIC :: nitrun !: number of iteration 163 INTEGER , PUBLIC :: numit !: iteration number 164 REAL(wp), PUBLIC :: rdt_ice !: ice time step 160 INTEGER , PUBLIC :: nstart !: iteration number of the begining of the run 161 INTEGER , PUBLIC :: nlast !: iteration number of the end of the run 162 INTEGER , PUBLIC :: nitrun !: number of iteration 163 INTEGER , PUBLIC :: numit !: iteration number 164 REAL(wp), PUBLIC :: rdt_ice !: ice time step 165 REAL(wp), PUBLIC :: r1_rdtice !: = 1. / rdt_ice 165 166 166 167 ! !!** ice-dynamic namelist (namicedyn) ** … … 201 202 ! !!** ice-salinity namelist (namicesal) ** 202 203 INTEGER , PUBLIC :: num_sal = 1 !: salinity configuration used in the model 203 ! ! 1 - s constant inspace and time204 ! ! 1 - constant salinity in both space and time 204 205 ! ! 2 - prognostic salinity (s(z,t)) 205 206 ! ! 3 - salinity profile, constant in time 206 ! ! 4 - salinity variations affect only ice thermodynamics207 207 INTEGER , PUBLIC :: sal_prof = 1 !: salinity profile or not 208 208 INTEGER , PUBLIC :: thcon_i_swi = 1 !: thermal conductivity: =1 Untersteiner (1964) ; =2 Pringle et al (2007) … … 264 264 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: phicif !: Old ice thickness 265 265 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fbif !: Heat flux at the ice base 266 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdmsnif !: Variation of snow mass 267 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdmicif !: Variation of ice mass 266 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_snw !: Variation of snow mass over 1 time step [Kg/m2] 267 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_snw !: Heat content associated with rdm_snw [J/m2] 268 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_ice !: Variation of ice mass over 1 time step [Kg/m2] 269 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_ice !: Heat content associated with rdm_ice [J/m2] 268 270 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qldif !: heat balance of the lead (or of the open ocean) 269 271 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qcmif !: Energy needed to bring the ocean to freezing … … 276 278 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qfvbq !: store energy in case of total lateral ablation (?) 277 279 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: dmgwi !: Variation of the mass of snow ice 278 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fsalt_res !: Residual salt flux due to correction of ice thickness279 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fsbri !: Salt flux due to brine rejection280 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fsalt_rpo !: Salt flux associated with porous ridged ice formation281 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fheat_rpo !: Heat flux associated with porous ridged ice formation280 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_thd !: salt flux due to ice growth/melt [PSU/m2/s] 281 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_bri !: salt flux due to brine rejection [PSU/m2/s] 282 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_mec !: salt flux due to porous ridged ice formation [PSU/m2/s] 283 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_res !: residual salt flux due to correction of ice thickness [PSU/m2/s] 282 284 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fhbri !: heat flux due to brine rejection 283 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: f mmec !: Mass flux due to snow loss during compression284 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: f seqv !: Equivalent salt flux due to ice growth/melt285 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: f hmec !: Heat flux due to snow loss during compression286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fh eat_res !: Residual heat flux due to correction of ice thickness285 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fheat_mec !: heat flux associated with porous ridged ice formation [???] 286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fheat_res !: residual heat flux due to correction of ice thickness 287 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fmmec !: mass flux due to snow loss during compression [Kg/m2/s] 288 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fhmec !: heat flux due to snow loss during compression 287 289 288 290 ! temporary arrays for dummy version of the code … … 415 417 416 418 !!---------------------------------------------------------------------- 417 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2010)419 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2010) 418 420 !! $Id$ 419 421 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 444 446 445 447 ii = ii + 1 446 ALLOCATE( firic (jpi,jpj) , fcsic (jpi,jpj) , fleic (jpi,jpj) , qlatic (jpi,jpj) , & 447 & rdvosif (jpi,jpj) , rdvobif(jpi,jpj) , fdvolif (jpi,jpj) , rdvonif (jpi,jpj) , & 448 & sist (jpi,jpj) , icethi (jpi,jpj) , t_bo (jpi,jpj) , hicifp (jpi,jpj) , & 449 & frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , fbif (jpi,jpj) , & 450 & rdmsnif (jpi,jpj) , rdmicif(jpi,jpj) , qldif (jpi,jpj) , qcmif (jpi,jpj) , & 451 & fdtcn (jpi,jpj) , qdtcn (jpi,jpj) , fstric (jpi,jpj) , fscmbq (jpi,jpj) , & 452 & ffltbif (jpi,jpj) , fsbbq (jpi,jpj) , qfvbq (jpi,jpj) , dmgwi (jpi,jpj) , & 453 & fsalt_res(jpi,jpj) , fsbri (jpi,jpj) , fsalt_rpo(jpi,jpj) , fheat_rpo(jpi,jpj) , & 454 & fhbri (jpi,jpj) , fmmec (jpi,jpj) , fseqv (jpi,jpj) , fhmec (jpi,jpj) , & 455 & fheat_res(jpi,jpj) , STAT=ierr(ii) ) 448 ALLOCATE( firic (jpi,jpj) , fcsic (jpi,jpj) , fleic (jpi,jpj) , qlatic (jpi,jpj) , & 449 & rdvosif (jpi,jpj) , rdvobif(jpi,jpj) , fdvolif(jpi,jpj) , rdvonif (jpi,jpj) , & 450 & sist (jpi,jpj) , icethi (jpi,jpj) , t_bo (jpi,jpj) , hicifp (jpi,jpj) , & 451 & frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , fbif (jpi,jpj) , & 452 & rdm_snw (jpi,jpj) , rdq_snw(jpi,jpj) , rdm_ice(jpi,jpj) , rdq_ice (jpi,jpj) , & 453 & qldif (jpi,jpj) , qcmif (jpi,jpj) , & 454 & fdtcn (jpi,jpj) , qdtcn (jpi,jpj) , fstric (jpi,jpj) , fscmbq (jpi,jpj) , & 455 & ffltbif (jpi,jpj) , fsbbq (jpi,jpj) , qfvbq (jpi,jpj) , dmgwi (jpi,jpj) , & 456 & sfx_res (jpi,jpj) , sfx_bri(jpi,jpj) , sfx_mec(jpi,jpj) , fheat_mec(jpi,jpj) , & 457 & fhbri (jpi,jpj) , fmmec (jpi,jpj) , sfx_thd(jpi,jpj) , fhmec (jpi,jpj) , & 458 & fheat_res(jpi,jpj) , STAT=ierr(ii) ) 456 459 457 460 ii = ii + 1 -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/iceini.F90
r3294 r3625 10 10 #if defined key_lim3 11 11 !!---------------------------------------------------------------------- 12 !! 'key_lim3' : LIM sea-ice model 13 !!---------------------------------------------------------------------- 14 !! ice_init : sea-ice model initialization 15 !!---------------------------------------------------------------------- 16 USE phycst ! physical constants 17 USE dom_oce ! ocean domain 18 USE sbc_oce ! Surface boundary condition: ocean fields 19 USE sbc_ice ! Surface boundary condition: ice fields 20 USE ice ! LIM variables 21 USE par_ice ! LIM parameters 22 USE dom_ice ! LIM domain 23 USE thd_ice ! LIM thermodynamical variables 24 USE limitd_me ! LIM ice thickness distribution 25 USE limmsh ! LIM mesh 26 USE limistate ! LIM initial state 27 USE limrst ! LIM restart 28 USE limthd ! LIM ice thermodynamics 29 USE limthd_sal ! LIM ice thermodynamics: salinity 30 USE limvar ! LIM variables 31 USE limsbc ! LIM surface boundary condition 32 USE in_out_manager ! I/O manager 33 USE lib_mpp ! MPP library 12 !! 'key_lim3' LIM sea-ice model 13 !!---------------------------------------------------------------------- 14 !! ice_init : sea-ice model initialization 15 !!---------------------------------------------------------------------- 16 USE phycst ! physical constants 17 USE dom_oce ! ocean domain 18 USE sbc_oce ! Surface boundary condition: ocean fields 19 USE sbc_ice ! Surface boundary condition: ice fields 20 USE ice ! LIM variables 21 USE par_ice ! LIM parameters 22 USE dom_ice ! LIM domain 23 USE thd_ice ! LIM thermodynamical variables 24 USE limitd_me ! LIM ice thickness distribution 25 USE limmsh ! LIM mesh 26 USE limistate ! LIM initial state 27 USE limrst ! LIM restart 28 USE limthd ! LIM ice thermodynamics 29 USE limthd_sal ! LIM ice thermodynamics: salinity 30 USE limvar ! LIM variables 31 USE limsbc ! LIM surface boundary condition 32 USE in_out_manager ! I/O manager 33 USE lib_mpp ! MPP library 34 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 34 35 35 36 IMPLICIT NONE … … 39 40 40 41 !!---------------------------------------------------------------------- 41 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)42 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 42 43 !! $Id$ 43 44 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 79 80 CALL lim_thd_sal_init ! set ice salinity parameters 80 81 ! 81 rdt_ice = nn_fsbc * rdttra(1) ! sea-ice timestep 82 rdt_ice = nn_fsbc * rdttra(1) ! sea-ice timestep 83 r1_rdtice = 1._wp / rdt_ice ! sea-ice timestep inverse 82 84 ! 83 85 CALL lim_msh ! ice mesh initialization -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limadv.F90
r3294 r3625 15 15 !! lim_adv_y : advection of sea ice on y axis 16 16 !!---------------------------------------------------------------------- 17 USE dom_oce ! ocean domain 18 USE dom_ice ! LIM-3 domain 19 USE ice ! LIM-3 variables 20 USE lbclnk ! lateral boundary condition - MPP exchanges 21 USE in_out_manager ! I/O manager 22 USE prtctl ! Print control 23 USE lib_mpp ! MPP library 24 USE wrk_nemo ! work arrays 17 USE dom_oce ! ocean domain 18 USE ice ! LIM-3 variables 19 USE dom_ice ! LIM-3 domain 20 USE lbclnk ! lateral boundary condition - MPP exchanges 21 USE in_out_manager ! I/O manager 22 USE prtctl ! Print control 23 USE lib_mpp ! MPP library 24 USE wrk_nemo ! work arrays 25 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 25 26 26 27 IMPLICIT NONE … … 37 38 # include "vectopt_loop_substitute.h90" 38 39 !!---------------------------------------------------------------------- 39 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)40 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 40 41 !! $Id$ 41 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 88 89 zs2new = MIN( 2.0 * zslpmax - 0.3334 * ABS( zs1new ), & 89 90 & MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj) ) ) 90 zin0 = ( 1.0 - MAX( rzero, sign( rone, -zslpmax) ) ) * tms(ji,jj) ! Case of empty boxes & Apply mask91 zin0 = ( 1.0 - MAX( rzero, SIGN( rone, -zslpmax) ) ) * tms(ji,jj) ! Case of empty boxes & Apply mask 91 92 92 93 ps0 (ji,jj) = zslpmax … … 273 274 zs2new = MIN( ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ), & 274 275 & MAX( ABS( zs1new )-zslpmax, psyy(ji,jj) ) ) 275 zin0 = ( 1.0 - MAX( rzero, sign( rone, -zslpmax) ) ) * tms(ji,jj) ! Case of empty boxes & Apply mask276 zin0 = ( 1.0 - MAX( rzero, SIGN( rone, -zslpmax) ) ) * tms(ji,jj) ! Case of empty boxes & Apply mask 276 277 ! 277 278 ps0 (ji,jj) = zslpmax -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90
r2777 r3625 10 10 #if defined key_lim3 11 11 !!---------------------------------------------------------------------- 12 !! 'key_lim3' : LIM3 sea-ice model12 !! 'key_lim3' LIM-3 sea-ice model 13 13 !!---------------------------------------------------------------------- 14 !! lim_cons : checks whether energy, mass and salt are conserved14 !! lim_cons : checks whether energy, mass and salt are conserved 15 15 !!---------------------------------------------------------------------- 16 USE par_ice ! LIM-3 parameter 17 USE ice ! LIM-3 variables 18 USE dom_ice ! LIM-3 domain 19 USE dom_oce ! ocean domain 20 USE in_out_manager ! I/O manager 21 USE lib_mpp ! MPP library 16 USE par_ice ! LIM-3 parameter 17 USE ice ! LIM-3 variables 18 USE dom_ice ! LIM-3 domain 19 USE dom_oce ! ocean domain 20 USE in_out_manager ! I/O manager 21 USE lib_mpp ! MPP library 22 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 22 23 23 24 IMPLICIT NONE … … 29 30 30 31 !!---------------------------------------------------------------------- 31 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)32 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 32 33 !! $Id$ 33 34 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limdia.F90
r2715 r3625 11 11 !! 'key_lim3' LIM3 sea-ice model 12 12 !!---------------------------------------------------------------------- 13 !! lim_dia : computation and output of the time evolution of keys variables 14 !! lim_dia_init : initialization and namelist read 15 !!---------------------------------------------------------------------- 16 USE ice ! LIM-3: sea-ice variable 17 USE par_ice ! LIM-3: ice parameters 18 USE dom_ice ! LIM-3: sea-ice domain 19 USE dom_oce ! ocean domain 20 USE sbc_oce ! surface boundary condition: ocean fields 21 USE daymod ! model calendar 22 USE phycst ! physical constant 23 USE in_out_manager ! I/O manager 24 USE lib_mpp ! MPP library 25 13 !! lim_dia : computation and output of the time evolution of keys variables 14 !! lim_dia_init : initialization and namelist read 15 !!---------------------------------------------------------------------- 16 USE ice ! LIM-3: sea-ice variable 17 USE par_ice ! LIM-3: ice parameters 18 USE dom_ice ! LIM-3: sea-ice domain 19 USE dom_oce ! ocean domain 20 USE sbc_oce ! surface boundary condition: ocean fields 21 USE daymod ! model calendar 22 USE phycst ! physical constant 23 USE in_out_manager ! I/O manager 24 USE lib_mpp ! MPP library 25 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 26 26 27 IMPLICIT NONE 27 28 PRIVATE … … 70 71 !! the temporal evolution of some key variables 71 72 !!------------------------------------------------------------------- 72 INTEGER :: jv, ji, jj, jl ! dummy loop indices 73 REAL(wp) :: zshift_date ! date from the minimum ice extent 74 REAL(wp) :: zday, zday_min ! current day, day of minimum extent 75 REAL(wp) :: zafy, zamy ! temporary area of fy and my ice 73 INTEGER :: jv, ji, jj, jl ! dummy loop indices 74 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integer 75 REAL(wp) :: zshift_date ! date from the minimum ice extent 76 REAL(wp) :: zday, zday_min ! current day, day of minimum extent 77 REAL(wp) :: zafy, zamy ! temporary area of fy and my ice 76 78 REAL(wp) :: zindb 77 REAL(wp), DIMENSION(jpinfmx) :: vinfor ! temporary workingspace79 REAL(wp), DIMENSION(jpinfmx) :: vinfor ! 1D workspace 78 80 !!------------------------------------------------------------------- 79 81 … … 105 107 IF( tms(ji,jj) == 1 ) THEN 106 108 vinfor(3) = vinfor(3) + at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice area 107 IF ( at_i(ji,jj).GT.0.15)vinfor(5) = vinfor(5) + aire(ji,jj) * 1.e-12_wp !ice extent109 IF ( at_i(ji,jj) > 0.15 ) vinfor(5) = vinfor(5) + aire(ji,jj) * 1.e-12_wp !ice extent 108 110 vinfor(7) = vinfor(7) + vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice volume 109 111 vinfor(9) = vinfor(9) + vt_s(ji,jj)*aire(ji,jj) * 1.e-12_wp !snow volume … … 111 113 vinfor(29) = vinfor(29) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !mean salinity 112 114 ! the computation of this diagnostic is not reliable 113 vinfor(31) = vinfor(31) + vt_i(ji,jj) *( u_ice(ji,jj)*u_ice(ji,jj) +&114 v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12115 vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) * 1.e-12_wp !salt flux116 vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) * 1.e-12_wp !brine drainage flux117 vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) * 1.e-12_wp !equivalent salt flux115 vinfor(31) = vinfor(31) + vt_i(ji,jj) * ( u_ice(ji,jj)*u_ice(ji,jj) & 116 & + v_ice(ji,jj)*v_ice(ji,jj) ) * aire(ji,jj) * 1.e-12 117 vinfor(53) = vinfor(53) + sfx (ji,jj)*aire(ji,jj) * 1.e-12_wp !salt flux 118 vinfor(55) = vinfor(55) + sfx_bri(ji,jj)*aire(ji,jj) * 1.e-12_wp !brine drainage flux 119 vinfor(57) = vinfor(57) + sfx_thd(ji,jj)*aire(ji,jj) * 1.e-12_wp !equivalent salt flux 118 120 vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !SST 119 121 vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !SSS … … 180 182 vinfor(43) = vinfor(43) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp 181 183 vinfor(45) = vinfor(45) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12_wp 182 vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp / rdt_ice! volume acc in OW184 vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp * r1_rdtice ! volume acc in OW 183 185 ENDIF 184 186 END DO … … 231 233 vinfor(51) = zindb*vinfor(51) / MAX(vinfor(27),epsi06) 232 234 233 !! Fram Strait Export 234 !! 83 = area export 235 !! 84 = volume export 236 !! Fram strait in ORCA2 = 5 points 237 !! export = -v_ice*e1t*ddtb*at_i or -v_ice*e1t*ddtb*at_i*h_i 238 jj = 136 ! C grid 239 vinfor(83) = 0.0 240 vinfor(84) = 0.0 241 DO ji = 134, 138 242 vinfor(83) = vinfor(83) - v_ice(ji,jj) * & 243 e1t(ji,jj)*at_i(ji,jj)*rdt_ice * 1.e-12_wp 244 vinfor(84) = vinfor(84) - v_ice(ji,jj) * & 245 e1t(ji,jj)*vt_i(ji,jj)*rdt_ice * 1.e-12_wp 246 END DO 235 IF( cp_cfg == "orca" ) THEN !* ORCA configuration : Fram Strait Export 236 SELECT CASE ( jp_cfg ) 237 CASE ( 2 ) ! ORCA_R2 238 ij0 = 136 ; ij1 = 136 ! Fram strait : 83 = area export 239 ii0 = 134 ; ii1 = 138 ! 84 = volume export 240 DO jj = mj0(ij0),mj1(ij1) 241 DO ji = mi0(ii0),mi1(ii1) 242 vinfor(83) = vinfor(83) - v_ice(ji,jj) * e1t(ji,jj)*at_i(ji,jj)*rdt_ice * 1.e-12_wp 243 vinfor(84) = vinfor(84) - v_ice(ji,jj) * e1t(ji,jj)*vt_i(ji,jj)*rdt_ice * 1.e-12_wp 244 END DO 245 END DO 246 END SELECT 247 !!gm just above, this is NOT the correct way of evaluating the transport ! 248 !!gm mass of snow is missing and v_ice should be the mean between jj and jj+1 249 !!gm Other ORCA configurations should be added 250 ENDIF 247 251 248 252 !!------------------------------------------------------------------- … … 264 268 vinfor(32) = vinfor(32) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + & 265 269 v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12 !ice vel 266 vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Total salt flux267 vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Brine drainage salt flux268 vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Equivalent salt flux270 vinfor(54) = vinfor(54) + sfx (ji,jj)*aire(ji,jj) * 1.e-12_wp ! Total salt flux 271 vinfor(56) = vinfor(56) + sfx_bri(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Brine drainage salt flux 272 vinfor(58) = vinfor(58) + sfx_thd(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Equivalent salt flux 269 273 vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !SST 270 274 vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !SSS … … 331 335 vinfor(44) = vinfor(44) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp 332 336 vinfor(46) = vinfor(46) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12_wp 333 vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp / rdt_ice! volume acc in OW337 vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp * r1_rdtice ! volume acc in OW 334 338 ENDIF 335 339 END DO … … 345 349 END DO 346 350 END DO 347 zindb = 1. 0 - MAX(0.0,SIGN(1.0,-vinfor(4)))!348 vinfor(64) = zindb * vinfor(64) / MAX( vinfor(4),epsi06)! divide by ice extt351 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , -vinfor(4) ) ) ! 352 vinfor(64) = zindb * vinfor(64) / MAX( vinfor(4) , epsi06 ) ! divide by ice extt 349 353 !! 2.2) Diagnostics dependent on age 350 354 !!------------------------------------ … … 368 372 ENDIF 369 373 END DO ! jl 370 IF ( (at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN374 IF ( at_i(ji,jj) > 0.15 .AND. zafy > zamy ) THEN 371 375 vinfor(22) = vinfor(22) + aire(ji,jj) * 1.e-12_wp ! Seasonal ice extent 372 376 ENDIF 373 IF ( (at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN377 IF ( at_i(ji,jj) > 0.15 .AND. zafy <= zamy ) THEN 374 378 vinfor(24) = vinfor(24) + aire(ji,jj) * 1.e-12_wp ! Perennial ice extent 375 379 ENDIF … … 377 381 END DO ! jj 378 382 END DO ! ji 379 zindb = 1.0 - MAX( 0.0,SIGN(1.0,-vinfor(26)))!=0 if no multiyear ice 1 if yes380 vinfor(50) = zindb *vinfor(50) / MAX(vinfor(26),epsi06)381 zindb = 1.0 - MAX( 0.0,SIGN(1.0,-vinfor(28))) !=0 if no multiyear ice 1 if yes382 vinfor(52) = zindb *vinfor(52) / MAX(vinfor(28),epsi06)383 zindb = 1.0 - MAX( 0.0,SIGN( 1._wp , -vinfor(26) ) ) !=0 if no multiyear ice 1 if yes 384 vinfor(50) = zindb * vinfor(50) / MAX( vinfor(26) , epsi06 ) 385 zindb = 1.0 - MAX( 0._wp , SIGN( 1._wp , -vinfor(28) ) ) !=0 if no multiyear ice 1 if yes 386 vinfor(52) = zindb * vinfor(52) / MAX( vinfor(28) , epsi06 ) 383 387 384 388 ! Accumulation before averaging -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90
r3294 r3625 15 15 !! lim_dyn_init : initialization and namelist read 16 16 !!---------------------------------------------------------------------- 17 USE phycst ! physical constants 18 USE dom_oce ! ocean space and time domain 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE sbc_ice ! Surface boundary condition: ice fields 21 USE ice ! LIM-3 variables 22 USE par_ice ! LIM-3 parameters 23 USE dom_ice ! LIM-3 domain 24 USE limrhg ! LIM-3 rheology 25 USE lbclnk ! lateral boundary conditions - MPP exchanges 26 USE lib_mpp ! MPP library 27 USE wrk_nemo ! work arrays 28 USE in_out_manager ! I/O manager 29 USE prtctl ! Print control 17 USE phycst ! physical constants 18 USE dom_oce ! ocean space and time domain 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE sbc_ice ! Surface boundary condition: ice fields 21 USE ice ! LIM-3 variables 22 USE par_ice ! LIM-3 parameters 23 USE dom_ice ! LIM-3 domain 24 USE limrhg ! LIM-3 rheology 25 USE lbclnk ! lateral boundary conditions - MPP exchanges 26 USE lib_mpp ! MPP library 27 USE wrk_nemo ! work arrays 28 USE in_out_manager ! I/O manager 29 USE prtctl ! Print control 30 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 30 31 31 32 IMPLICIT NONE … … 37 38 # include "vectopt_loop_substitute.h90" 38 39 !!---------------------------------------------------------------------- 39 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)40 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 40 41 !! $Id$ 41 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limhdf.F90
r3294 r3625 12 12 !! 'key_lim3' LIM3 sea-ice model 13 13 !!---------------------------------------------------------------------- 14 !! lim_hdf : diffusion trend on sea-ice variable14 !! lim_hdf : diffusion trend on sea-ice variable 15 15 !!---------------------------------------------------------------------- 16 USE dom_oce ! ocean domain 17 USE ice ! LIM-3: ice variables 18 USE lbclnk ! lateral boundary condition - MPP exchanges 19 USE lib_mpp ! MPP library 20 USE wrk_nemo ! work arrays 21 USE prtctl ! Print control 22 USE in_out_manager ! I/O manager 16 USE dom_oce ! ocean domain 17 USE ice ! LIM-3: ice variables 18 USE lbclnk ! lateral boundary condition - MPP exchanges 19 USE lib_mpp ! MPP library 20 USE wrk_nemo ! work arrays 21 USE prtctl ! Print control 22 USE in_out_manager ! I/O manager 23 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 23 24 24 25 IMPLICIT NONE … … 34 35 # include "vectopt_loop_substitute.h90" 35 36 !!---------------------------------------------------------------------- 36 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2010)37 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2010) 37 38 !! $Id$ 38 39 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90
r3610 r3625 26 26 USE lib_mpp ! MPP library 27 27 USE wrk_nemo ! work arrays 28 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 28 29 29 30 IMPLICIT NONE … … 48 49 49 50 !!---------------------------------------------------------------------- 50 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)51 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 51 52 !! $Id$ 52 53 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r3294 r3625 10 10 #if defined key_lim3 11 11 !!---------------------------------------------------------------------- 12 !! 'key_lim3' : LIM3 sea-ice model12 !! 'key_lim3' LIM-3 sea-ice model 13 13 !!---------------------------------------------------------------------- 14 USE par_oce ! ocean parameters 15 USE dom_oce ! ocean domain 16 USE phycst ! physical constants (ocean directory) 17 USE sbc_oce ! surface boundary condition: ocean fields 18 USE thd_ice ! LIM thermodynamics 19 USE ice ! LIM variables 20 USE par_ice ! LIM parameters 21 USE dom_ice ! LIM domain 22 USE limthd_lac ! LIM 23 USE limvar ! LIM 24 USE limcons ! LIM 25 USE in_out_manager ! I/O manager 26 USE lbclnk ! lateral boundary condition - MPP exchanges 27 USE lib_mpp ! MPP library 28 USE wrk_nemo ! work arrays 29 USE prtctl ! Print control 14 USE par_oce ! ocean parameters 15 USE dom_oce ! ocean domain 16 USE phycst ! physical constants (ocean directory) 17 USE sbc_oce ! surface boundary condition: ocean fields 18 USE thd_ice ! LIM thermodynamics 19 USE ice ! LIM variables 20 USE par_ice ! LIM parameters 21 USE dom_ice ! LIM domain 22 USE limthd_lac ! LIM 23 USE limvar ! LIM 24 USE limcons ! LIM 25 USE in_out_manager ! I/O manager 26 USE lbclnk ! lateral boundary condition - MPP exchanges 27 USE lib_mpp ! MPP library 28 USE wrk_nemo ! work arrays 29 USE prtctl ! Print control 30 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 30 31 31 32 IMPLICIT NONE … … 38 39 PUBLIC lim_itd_me_alloc ! called by iceini.F90 39 40 40 REAL(wp) 41 REAL(wp) 42 REAL(wp) 41 REAL(wp) :: epsi11 = 1.e-11_wp ! constant values 42 REAL(wp) :: epsi10 = 1.e-10_wp ! constant values 43 REAL(wp) :: epsi06 = 1.e-06_wp ! constant values 43 44 44 45 !----------------------------------------------------------------------- … … 47 48 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: asum ! sum of total ice and open water area 48 49 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: aksum ! ratio of area removed to area ridged 49 50 ! 50 51 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: athorn ! participation function; fraction of ridging/ 51 52 ! ! closing associated w/ category n 52 53 ! 53 54 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hrmin ! minimum ridge thickness 54 55 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hrmax ! maximum ridge thickness … … 70 71 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: dvirdgdt ! rate of ice volume ridged (m/s) 71 72 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: opening ! rate of opening due to divergence/shear (1/s) 73 ! 72 74 !!---------------------------------------------------------------------- 73 75 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) … … 126 128 INTEGER :: ji, jj, jk, jl ! dummy loop index 127 129 INTEGER :: niter, nitermax = 20 ! local integer 128 LOGICAL :: asum_error 129 INTEGER :: iterate_ridging 130 REAL(wp) :: w1, tmpfac , dti! local scalar130 LOGICAL :: asum_error ! flag for asum .ne. 1 131 INTEGER :: iterate_ridging ! if true, repeat the ridging 132 REAL(wp) :: w1, tmpfac ! local scalar 131 133 CHARACTER (len = 15) :: fieldid 132 134 REAL(wp), POINTER, DIMENSION(:,:) :: closing_net ! net rate at which area is removed (1/s) … … 152 154 ! 1) Thickness categories boundaries, ice / o.w. concentrations, init_ons 153 155 !-----------------------------------------------------------------------------! 154 ! Set hi_max(ncat) to a big value to ensure that all ridged ice 155 ! is thinner than hi_max(ncat). 156 ! Set hi_max(ncat) to a big value to ensure that all ridged ice is thinner than hi_max(ncat). 156 157 157 158 hi_max(jpl) = 999.99 158 159 159 Cp = 0.5 * grav * (rau0-rhoic) * rhoic / rau0 ! proport const for PE 160 CALL lim_itd_me_ridgeprep ! prepare ridging 161 160 Cp = 0.5 * grav * (rau0-rhoic) * rhoic / rau0 ! proport const for PE 161 ! 162 CALL lim_itd_me_ridgeprep ! prepare ridging 163 ! 162 164 IF( con_i) CALL lim_column_sum( jpl, v_i, vt_i_init ) ! conservation check 163 165 … … 166 168 msnow_mlt(ji,jj) = 0._wp 167 169 esnow_mlt(ji,jj) = 0._wp 168 dardg1dt (ji,jj) 169 dardg2dt (ji,jj) 170 dvirdgdt (ji,jj) 171 opening (ji,jj) 170 dardg1dt (ji,jj) = 0._wp 171 dardg2dt (ji,jj) = 0._wp 172 dvirdgdt (ji,jj) = 0._wp 173 opening (ji,jj) = 0._wp 172 174 173 175 !-----------------------------------------------------------------------------! … … 201 203 ! to give asum = 1.0 after ridging. 202 204 203 divu_adv(ji,jj) = ( 1._wp - asum(ji,jj) ) / rdt_ice ! asum found in ridgeprep205 divu_adv(ji,jj) = ( 1._wp - asum(ji,jj) ) * r1_rdtice ! asum found in ridgeprep 204 206 205 207 IF( divu_adv(ji,jj) < 0._wp ) closing_net(ji,jj) = MAX( closing_net(ji,jj), -divu_adv(ji,jj) ) … … 207 209 ! 2.3 opning 208 210 !------------ 209 ! Compute the (non-negative) opening rate that will give 210 ! asum = 1.0 after ridging. 211 ! Compute the (non-negative) opening rate that will give asum = 1.0 after ridging. 211 212 opning(ji,jj) = closing_net(ji,jj) + divu_adv(ji,jj) 212 213 END DO … … 257 258 IF ( a_i(ji,jj,jl) > epsi11 .AND. athorn(ji,jj,jl) > 0._wp )THEN 258 259 w1 = athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice 259 IF ( w1 >a_i(ji,jj,jl) ) THEN260 IF ( w1 > a_i(ji,jj,jl) ) THEN 260 261 tmpfac = a_i(ji,jj,jl) / w1 261 262 closing_gross(ji,jj) = closing_gross(ji,jj) * tmpfac … … 291 292 ELSE 292 293 iterate_ridging = 1 293 divu_adv (ji,jj) = (1._wp - asum(ji,jj)) / rdt_ice294 divu_adv (ji,jj) = (1._wp - asum(ji,jj)) * r1_rdtice 294 295 closing_net(ji,jj) = MAX( 0._wp, -divu_adv(ji,jj) ) 295 296 opning (ji,jj) = MAX( 0._wp, divu_adv(ji,jj) ) … … 308 309 309 310 IF( iterate_ridging == 1 ) THEN 310 IF( niter .GT.nitermax ) THEN311 IF( niter > nitermax ) THEN 311 312 WRITE(numout,*) ' ALERTE : non-converging ridging scheme ' 312 313 WRITE(numout,*) ' niter, iterate_ridging ', niter, iterate_ridging … … 323 324 ! Update fresh water and heat fluxes due to snow melt. 324 325 325 dti = 1._wp / rdt_ice326 327 326 asum_error = .false. 328 327 … … 330 329 DO ji = 1, jpi 331 330 332 IF (ABS(asum(ji,jj) - 1.0) .GT. epsi11)asum_error = .true.333 334 dardg1dt(ji,jj) = dardg1dt(ji,jj) * dti335 dardg2dt(ji,jj) = dardg2dt(ji,jj) * dti336 dvirdgdt(ji,jj) = dvirdgdt(ji,jj) * dti337 opening (ji,jj) = opening (ji,jj) * dti331 IF( ABS( asum(ji,jj) - 1.0 ) > epsi11 ) asum_error = .true. 332 333 dardg1dt(ji,jj) = dardg1dt(ji,jj) * r1_rdtice 334 dardg2dt(ji,jj) = dardg2dt(ji,jj) * r1_rdtice 335 dvirdgdt(ji,jj) = dvirdgdt(ji,jj) * r1_rdtice 336 opening (ji,jj) = opening (ji,jj) * r1_rdtice 338 337 339 338 !-----------------------------------------------------------------------------! 340 339 ! 5) Heat, salt and freshwater fluxes 341 340 !-----------------------------------------------------------------------------! 342 fmmec(ji,jj) = fmmec(ji,jj) + msnow_mlt(ji,jj) * dti! fresh water source for ocean343 fhmec(ji,jj) = fhmec(ji,jj) + esnow_mlt(ji,jj) * dti! heat sink for ocean341 fmmec(ji,jj) = fmmec(ji,jj) + msnow_mlt(ji,jj) * r1_rdtice ! fresh water source for ocean 342 fhmec(ji,jj) = fhmec(ji,jj) + esnow_mlt(ji,jj) * r1_rdtice ! heat sink for ocean 344 343 345 344 END DO … … 349 348 DO jj = 1, jpj 350 349 DO ji = 1, jpi 351 IF (ABS(asum(ji,jj) - 1.0) .GT. epsi11) THEN! there is a bug350 IF( ABS( asum(ji,jj) - 1._wp ) > epsi11 ) THEN ! there is a bug 352 351 WRITE(numout,*) ' ' 353 352 WRITE(numout,*) ' ALERTE : Ridging error: total area = ', asum(ji,jj) … … 391 390 d_oa_i_trp (:,:,:) = oa_i (:,:,:) - old_oa_i (:,:,:) 392 391 d_smv_i_trp(:,:,:) = 0._wp 393 IF( num_sal == 2 .OR. num_sal == 4) d_smv_i_trp(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:)392 IF( num_sal == 2 ) d_smv_i_trp(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:) 394 393 395 394 IF(ln_ctl) THEN ! Control print … … 430 429 431 430 ! update of fields will be made later in lim update 432 u_ice(:,:) 433 v_ice(:,:) 434 a_i(:,:,:) 435 v_s(:,:,:) 436 v_i(:,:,:) 437 e_s(:,:,:,:) 438 e_i(:,:,:,:) 439 oa_i(:,:,:) 440 IF( num_sal == 2 .OR. num_sal == 4 ) smv_i(:,:,:)= old_smv_i(:,:,:)431 u_ice(:,:) = old_u_ice(:,:) 432 v_ice(:,:) = old_v_ice(:,:) 433 a_i(:,:,:) = old_a_i(:,:,:) 434 v_s(:,:,:) = old_v_s(:,:,:) 435 v_i(:,:,:) = old_v_i(:,:,:) 436 e_s(:,:,:,:) = old_e_s(:,:,:,:) 437 e_i(:,:,:,:) = old_e_i(:,:,:,:) 438 oa_i(:,:,:) = old_oa_i(:,:,:) 439 IF( num_sal == 2 ) smv_i(:,:,:) = old_smv_i(:,:,:) 441 440 442 441 !----------------------------------------------------! … … 465 464 DO jj = 1, jpj 466 465 DO ji = 1, jpi 467 IF ( ( old_v_i(ji,jj,jl) < epsi06 ).AND. &468 ( d_v_i_trp(ji,jj,jl) > epsi06 )) THEN469 old_v_i (ji,jj,jl)= d_v_i_trp(ji,jj,jl)470 d_v_i_trp (ji,jj,jl) = 0._wp471 old_a_i (ji,jj,jl)= d_a_i_trp(ji,jj,jl)472 d_a_i_trp (ji,jj,jl) = 0._wp473 old_v_s (ji,jj,jl)= d_v_s_trp(ji,jj,jl)474 d_v_s_trp (ji,jj,jl) = 0._wp475 old_e_s (ji,jj,1,jl)= d_e_s_trp(ji,jj,1,jl)476 d_e_s_trp (ji,jj,1,jl) = 0._wp477 old_oa_i (ji,jj,jl)= d_oa_i_trp(ji,jj,jl)478 d_oa_i_trp(ji,jj,jl) = 0._wp479 IF( num_sal == 2 .OR. num_sal == 4 ) old_smv_i(ji,jj,jl)= d_smv_i_trp(ji,jj,jl)480 d_smv_i_trp(ji,jj,jl) = 0._wp466 IF( old_v_i (ji,jj,jl) < epsi06 .AND. & 467 d_v_i_trp(ji,jj,jl) > epsi06 ) THEN 468 old_v_i (ji,jj,jl) = d_v_i_trp(ji,jj,jl) 469 d_v_i_trp (ji,jj,jl) = 0._wp 470 old_a_i (ji,jj,jl) = d_a_i_trp(ji,jj,jl) 471 d_a_i_trp (ji,jj,jl) = 0._wp 472 old_v_s (ji,jj,jl) = d_v_s_trp(ji,jj,jl) 473 d_v_s_trp (ji,jj,jl) = 0._wp 474 old_e_s (ji,jj,1,jl) = d_e_s_trp(ji,jj,1,jl) 475 d_e_s_trp (ji,jj,1,jl) = 0._wp 476 old_oa_i (ji,jj,jl) = d_oa_i_trp(ji,jj,jl) 477 d_oa_i_trp(ji,jj,jl) = 0._wp 478 IF( num_sal == 2 ) old_smv_i(ji,jj,jl) = d_smv_i_trp(ji,jj,jl) 479 d_smv_i_trp(ji,jj,jl) = 0._wp 481 480 ENDIF 482 481 END DO 483 482 END DO 484 483 END DO 485 484 ! 486 485 CALL wrk_dealloc( jpi, jpj, closing_net, divu_adv, opning, closing_gross, msnow_mlt, esnow_mlt, vt_i_init, vt_i_final ) 487 486 ! … … 612 611 ! present 613 612 zworka(ji,jj) = 4.0 * strength(ji,jj) & 614 + strength(ji-1,jj) * tms(ji-1,jj) &615 + strength(ji+1,jj) * tms(ji+1,jj) &616 + strength(ji,jj-1) * tms(ji,jj-1) &617 + strength(ji,jj+1) * tms(ji,jj+1)613 & + strength(ji-1,jj) * tms(ji-1,jj) & 614 & + strength(ji+1,jj) * tms(ji+1,jj) & 615 & + strength(ji,jj-1) * tms(ji,jj-1) & 616 & + strength(ji,jj+1) * tms(ji,jj+1) 618 617 619 618 zw1 = 4.0 + tms(ji-1,jj) + tms(ji+1,jj) + tms(ji,jj-1) + tms(ji,jj+1) 620 619 zworka(ji,jj) = zworka(ji,jj) / zw1 621 620 ELSE 622 zworka(ji,jj) = 0. 0621 zworka(ji,jj) = 0._wp 623 622 ENDIF 624 623 END DO … … 1048 1047 DO jj = 1, jpj 1049 1048 DO ji = 1, jpi 1050 IF (aicen_init(ji,jj,jl1) .GT. epsi11 .AND. athorn(ji,jj,jl1) .GT. 0.0&1051 .AND. closing_gross(ji,jj) > 0. 0) THEN1049 IF( aicen_init(ji,jj,jl1) > epsi11 .AND. athorn(ji,jj,jl1) > 0._wp & 1050 .AND. closing_gross(ji,jj) > 0._wp ) THEN 1052 1051 icells = icells + 1 1053 1052 indxi(icells) = ji … … 1130 1129 ! Salinity 1131 1130 !------------- 1132 smsw(ji,jj) = sss_m(ji,jj) * vsw(ji,jj) * rhoic / rau0 ! salt content of water frozen in voids1131 smsw(ji,jj) = sss_m(ji,jj) * vsw(ji,jj) * rhoic / rau0 ! salt content of seawater frozen in voids 1133 1132 1134 1133 zsrdg2 = srdg1(ji,jj) + smsw(ji,jj) ! salt content of new ridge … … 1137 1136 1138 1137 ! ! excess of salt is flushed into the ocean 1139 fsalt_rpo(ji,jj) = fsalt_rpo(ji,jj) + ( zsrdg2 - srdg2(ji,jj) ) * rhoic / rdt_ice 1140 1138 sfx_mec(ji,jj) = sfx_mec(ji,jj) + ( zsrdg2 - srdg2(ji,jj) ) * rhoic * r1_rdtice 1139 1140 rdm_ice(ji,jj) = rdm_ice(ji,jj) + vsw(ji,jj) * rhoic / rau0 ! increase in ice volume du to seawater frozen in voids 1141 1141 1142 !------------------------------------ 1142 1143 ! 3.6 Increment ridging diagnostics … … 1148 1149 dardg1dt (ji,jj) = dardg1dt(ji,jj) + ardg1(ji,jj) + arft1(ji,jj) 1149 1150 dardg2dt (ji,jj) = dardg2dt(ji,jj) + ardg2(ji,jj) + arft2(ji,jj) 1150 diag_dyn_gr(ji,jj) = diag_dyn_gr(ji,jj) + ( vrdg2(ji,jj) + virft(ji,jj) ) / rdt_ice1151 opening (ji,jj) = opening (ji,jj) + opning(ji,jj) *rdt_ice1151 diag_dyn_gr(ji,jj) = diag_dyn_gr(ji,jj) + ( vrdg2(ji,jj) + virft(ji,jj) ) * r1_rdtice 1152 opening (ji,jj) = opening (ji,jj) + opning(ji,jj) * rdt_ice 1152 1153 1153 1154 IF( con_i ) vice_init(ji,jj) = vice_init(ji,jj) + vrdg2(ji,jj) - vrdg1(ji,jj) … … 1156 1157 ! 3.7 Put the snow somewhere in the ocean 1157 1158 !------------------------------------------ 1158 1159 1159 ! Place part of the snow lost by ridging into the ocean. 1160 1160 ! Note that esnow_mlt < 0; the ocean must cool to melt snow. … … 1179 1179 ! ij looping 1-icells 1180 1180 1181 dhr (ji,jj)= hrmax(ji,jj,jl1) - hrmin(ji,jj,jl1)1181 dhr (ji,jj) = hrmax(ji,jj,jl1) - hrmin(ji,jj,jl1) 1182 1182 dhr2(ji,jj) = hrmax(ji,jj,jl1) * hrmax(ji,jj,jl1) - hrmin(ji,jj,jl1) * hrmin(ji,jj,jl1) 1183 1184 1183 1185 1184 END DO ! ij … … 1211 1210 1212 1211 ! heat flux 1213 fheat_ rpo(ji,jj) = fheat_rpo(ji,jj) + ( ersw(ji,jj,jk) - zdummy0 ) / rdt_ice1212 fheat_mec(ji,jj) = fheat_mec(ji,jj) + ( ersw(ji,jj,jk) - zdummy0 ) * r1_rdtice 1214 1213 1215 1214 ! Correct dimensions to avoid big values … … 1275 1274 ! Transfer area, volume, and energy accordingly. 1276 1275 1277 IF (hrmin(ji,jj,jl1) .GE. hi_max(jl2).OR. &1278 hrmax(ji,jj,jl1) .LE. hi_max(jl2-1)) THEN1279 hL = 0. 01280 hR = 0. 01276 IF( hrmin(ji,jj,jl1) >= hi_max(jl2) .OR. & 1277 hrmax(ji,jj,jl1) <= hi_max(jl2-1) ) THEN 1278 hL = 0._wp 1279 hR = 0._wp 1281 1280 ELSE 1282 hL = MAX (hrmin(ji,jj,jl1), hi_max(jl2-1))1283 hR = MIN (hrmax(ji,jj,jl1), hi_max(jl2))1281 hL = MAX( hrmin(ji,jj,jl1), hi_max(jl2-1) ) 1282 hR = MIN( hrmax(ji,jj,jl1), hi_max(jl2) ) 1284 1283 ENDIF 1285 1284 1286 1285 ! fraction of ridged ice area and volume going to n2 1287 farea = ( hR-hL) / dhr(ji,jj)1288 fvol(ji,jj) = ( hR*hR - hL*hL) / dhr2(ji,jj)1289 1290 a_i (ji,jj ,jl2) = a_i (ji,jj,jl2)+ ardg2 (ji,jj) * farea1291 v_i (ji,jj ,jl2) = v_i (ji,jj,jl2)+ vrdg2 (ji,jj) * fvol(ji,jj)1292 v_s (ji,jj ,jl2) = v_s (ji,jj,jl2)+ vsrdg (ji,jj) * fvol(ji,jj) * fsnowrdg1286 farea = ( hR - hL ) / dhr(ji,jj) 1287 fvol(ji,jj) = ( hR*hR - hL*hL ) / dhr2(ji,jj) 1288 1289 a_i (ji,jj ,jl2) = a_i (ji,jj ,jl2) + ardg2 (ji,jj) * farea 1290 v_i (ji,jj ,jl2) = v_i (ji,jj ,jl2) + vrdg2 (ji,jj) * fvol(ji,jj) 1291 v_s (ji,jj ,jl2) = v_s (ji,jj ,jl2) + vsrdg (ji,jj) * fvol(ji,jj) * fsnowrdg 1293 1292 e_s (ji,jj,1,jl2) = e_s (ji,jj,1,jl2) + esrdg (ji,jj) * fvol(ji,jj) * fsnowrdg 1294 smv_i(ji,jj ,jl2) = smv_i(ji,jj,jl2)+ srdg2 (ji,jj) * fvol(ji,jj)1295 oa_i (ji,jj ,jl2) = oa_i (ji,jj,jl2)+ oirdg2(ji,jj) * farea1293 smv_i(ji,jj ,jl2) = smv_i(ji,jj ,jl2) + srdg2 (ji,jj) * fvol(ji,jj) 1294 oa_i (ji,jj ,jl2) = oa_i (ji,jj ,jl2) + oirdg2(ji,jj) * farea 1296 1295 1297 1296 END DO ! ij … … 1317 1316 ! Compute the fraction of rafted ice area and volume going to 1318 1317 ! thickness category jl2, transfer area, volume, and energy accordingly. 1319 1320 IF (hraft(ji,jj,jl1) .LE. hi_max(jl2).AND. &1321 hraft(ji,jj,jl1) .GT. hi_max(jl2-1)) THEN1322 a_i (ji,jj,jl2) = a_i(ji,jj,jl2) + arft2(ji,jj)1323 v_i (ji,jj,jl2) = v_i(ji,jj,jl2) + virft(ji,jj)1324 v_s (ji,jj,jl2) = v_s(ji,jj,jl2) + vsrft(ji,jj)*fsnowrft1325 e_s (ji,jj,1,jl2) = e_s(ji,jj,1,jl2) + esrft(ji,jj)*fsnowrft1326 smv_i(ji,jj ,jl2) = smv_i(ji,jj,jl2) + smrft(ji,jj)1327 oa_i (ji,jj,jl2) = oa_i(ji,jj,jl2)+ oirft2(ji,jj)1318 ! 1319 IF( hraft(ji,jj,jl1) <= hi_max(jl2) .AND. & 1320 hraft(ji,jj,jl1) > hi_max(jl2-1) ) THEN 1321 a_i (ji,jj ,jl2) = a_i (ji,jj ,jl2) + arft2 (ji,jj) 1322 v_i (ji,jj ,jl2) = v_i (ji,jj ,jl2) + virft (ji,jj) 1323 v_s (ji,jj ,jl2) = v_s (ji,jj ,jl2) + vsrft (ji,jj) * fsnowrft 1324 e_s (ji,jj,1,jl2) = e_s (ji,jj,1,jl2) + esrft (ji,jj) * fsnowrft 1325 smv_i(ji,jj ,jl2) = smv_i(ji,jj ,jl2) + smrft (ji,jj) 1326 oa_i (ji,jj ,jl2) = oa_i (ji,jj ,jl2) + oirft2(ji,jj) 1328 1327 ENDIF ! hraft 1329 1328 ! 1330 1329 END DO ! ij 1331 1330 … … 1336 1335 ji = indxi(ij) 1337 1336 jj = indxj(ij) 1338 IF (hraft(ji,jj,jl1) .LE. hi_max(jl2).AND. &1339 hraft(ji,jj,jl1) .GT. hi_max(jl2-1)) THEN1337 IF( hraft(ji,jj,jl1) <= hi_max(jl2) .AND. & 1338 hraft(ji,jj,jl1) > hi_max(jl2-1) ) THEN 1340 1339 e_i(ji,jj,jk,jl2) = e_i(ji,jj,jk,jl2) + eirft(ji,jj,jk) 1341 1340 ENDIF … … 1504 1503 DO jj = 1 , jpj 1505 1504 DO ji = 1 , jpi 1506 !!gm xtmp = e_i(ji,jj,jk,jl) / area(ji,jj) / rdt_ice1505 !!gm xtmp = e_i(ji,jj,jk,jl) / area(ji,jj) * r1_rdtice 1507 1506 !!gm xtmp = xtmp * unit_fac 1508 1507 ! fheat_res(ji,jj) = fheat_res(ji,jj) - xtmp … … 1524 1523 ! fluxes are positive to the ocean 1525 1524 ! here the flux has to be negative for the ocean 1526 !!gm xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) / rdt_ice1525 !!gm xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) * r1_rdtice 1527 1526 ! fheat_res(ji,jj) = fheat_res(ji,jj) - xtmp 1528 1527 1529 !!gm xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) / rdt_ice !RB ???????1528 !!gm xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) * r1_rdtice !RB ??????? 1530 1529 1531 1530 t_s(ji,jj,1,jl) = rtt * zmask(ji,jj) + t_s(ji,jj,1,jl) * ( 1 - zmask(ji,jj) ) … … 1536 1535 1537 1536 ! xtmp = (rhoi*vicen(i,j,n) + rhos*vsnon(i,j,n)) / dt 1538 ! fresh(i,j) = fresh(i,j) + xtmp 1539 ! fresh_hist(i,j) = fresh_hist(i,j) + xtmp 1540 1541 ! fsalt_res(ji,jj) = fsalt_res(ji,jj) + ( sss_m(ji,jj) ) * & 1542 ! rhosn * v_s(ji,jj,jl) / rdt_ice 1543 1544 ! fsalt_res(ji,jj) = fsalt_res(ji,jj) + ( sss_m(ji,jj) - sm_i(ji,jj,jl) ) * & 1545 ! rhoic * v_i(ji,jj,jl) / rdt_ice 1546 1547 ! emps(i,j) = emps(i,j) + xtmp 1548 ! fsalt_hist(i,j) = fsalt_hist(i,j) + xtmp 1537 ! sfx_res(ji,jj) = sfx_res(ji,jj) + ( sss_m(ji,jj) ) & 1538 ! * rhosn * v_s(ji,jj,jl) * r1_rdtice 1539 ! sfx_res(ji,jj) = sfx_res(ji,jj) + ( sss_m(ji,jj) - sm_i(ji,jj,jl) ) & 1540 ! * rhoic * v_i(ji,jj,jl) * r1_rdtice 1541 ! sfx (i,j) = sfx (i,j) + xtmp 1549 1542 1550 1543 ato_i(ji,jj) = a_i (ji,jj,jl) * zmask(ji,jj) + ato_i(ji,jj) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90
r3294 r3625 2 2 !!====================================================================== 3 3 !! *** MODULE limitd_th *** 4 !! Thermodynamics of ice thickness distribution 5 !! computation of changes in g(h) 4 !! LIM3 ice model : ice thickness distribution: Thermodynamics 6 5 !!====================================================================== 7 6 !! History : - ! (W. H. Lipscomb and E.C. Hunke) CICE (c) original code … … 20 19 !! lim_itd_shiftice : 21 20 !!---------------------------------------------------------------------- 22 USE dom_ice ! LIM-3 domain 23 USE par_oce ! ocean parameters 24 USE dom_oce ! ocean domain 25 USE phycst ! physical constants (ocean directory) 26 USE thd_ice ! LIM-3 thermodynamic variables 27 USE ice ! LIM-3 variables 28 USE par_ice ! LIM-3 parameters 29 USE limthd_lac ! LIM-3 lateral accretion 30 USE limvar ! LIM-3 variables 31 USE limcons ! LIM-3 conservation 32 USE prtctl ! Print control 33 USE in_out_manager ! I/O manager 34 USE lib_mpp ! MPP library 35 USE wrk_nemo ! work arrays 21 USE par_oce ! ocean parameters 22 USE dom_oce ! ocean domain 23 USE phycst ! physical constants (ocean directory) 24 USE ice ! LIM-3 variables 25 USE par_ice ! LIM-3 parameters 26 USE dom_ice ! LIM-3 domain 27 USE thd_ice ! LIM-3 thermodynamic variables 28 USE limthd_lac ! LIM-3 lateral accretion 29 USE limvar ! LIM-3 variables 30 USE limcons ! LIM-3 conservation 31 USE prtctl ! Print control 32 USE in_out_manager ! I/O manager 33 USE lib_mpp ! MPP library 34 USE wrk_nemo ! work arrays 35 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 36 36 37 37 IMPLICIT NONE … … 49 49 50 50 !!---------------------------------------------------------------------- 51 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2010)51 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2010) 52 52 !! $Id$ 53 53 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 101 101 102 102 !- Trend terms 103 d_a_i_thd (:,:,:)= a_i(:,:,:) - old_a_i(:,:,:)104 d_v_s_thd (:,:,:)= v_s(:,:,:) - old_v_s(:,:,:)105 d_v_i_thd (:,:,:)= v_i(:,:,:) - old_v_i(:,:,:)103 d_a_i_thd(:,:,:) = a_i(:,:,:) - old_a_i(:,:,:) 104 d_v_s_thd(:,:,:) = v_s(:,:,:) - old_v_s(:,:,:) 105 d_v_i_thd(:,:,:) = v_i(:,:,:) - old_v_i(:,:,:) 106 106 d_e_s_thd(:,:,:,:) = e_s(:,:,:,:) - old_e_s(:,:,:,:) 107 107 d_e_i_thd(:,:,:,:) = e_i(:,:,:,:) - old_e_i(:,:,:,:) 108 108 109 109 d_smv_i_thd(:,:,:) = 0._wp 110 IF( num_sal == 2 .OR. num_sal == 4) d_smv_i_thd(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:)110 IF( num_sal == 2 ) d_smv_i_thd(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:) 111 111 112 112 IF(ln_ctl) THEN ! Control print … … 143 143 144 144 !- Recover Old values 145 a_i(:,:,:) = old_a_i 146 v_s(:,:,:) = old_v_s 147 v_i(:,:,:) = old_v_i 148 e_s(:,:,:,:) = old_e_s 149 e_i(:,:,:,:) = old_e_i 150 ! 151 IF( num_sal == 2 .OR. num_sal == 4 ) smv_i(:,:,:) = old_smv_i(:,:,:)145 a_i(:,:,:) = old_a_i(:,:,:) 146 v_s(:,:,:) = old_v_s(:,:,:) 147 v_i(:,:,:) = old_v_i(:,:,:) 148 e_s(:,:,:,:) = old_e_s(:,:,:,:) 149 e_i(:,:,:,:) = old_e_i(:,:,:,:) 150 ! 151 IF( num_sal == 2 ) smv_i(:,:,:) = old_smv_i(:,:,:) 152 152 ! 153 153 END SUBROUTINE lim_itd_th -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limmsh.F90
r2715 r3625 10 10 !! 'key_lim3' LIM3 sea-ice model 11 11 !!---------------------------------------------------------------------- 12 !! lim_msh : definition of the ice mesh12 !! lim_msh : definition of the ice mesh 13 13 !!---------------------------------------------------------------------- 14 14 USE phycst ! physical constants … … 18 18 USE lbclnk ! lateral boundary condition - MPP exchanges 19 19 USE lib_mpp ! MPP library 20 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 20 21 21 22 IMPLICIT NONE … … 25 26 26 27 !!---------------------------------------------------------------------- 27 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)28 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 28 29 !! $Id$ 29 30 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90
r3294 r3625 15 15 !! 'key_lim2' AND NOT 'key_lim2_vp' EVP LIM-2 sea-ice model 16 16 !!---------------------------------------------------------------------- 17 !! lim_rhg : computes ice velocities17 !! lim_rhg : computes ice velocities 18 18 !!---------------------------------------------------------------------- 19 USE phycst ! Physical constant 20 USE par_oce ! Ocean parameters 21 USE dom_oce ! Ocean domain 22 USE sbc_oce ! Surface boundary condition: ocean fields 23 USE sbc_ice ! Surface boundary condition: ice fields 24 USE lbclnk ! Lateral Boundary Condition / MPP link 25 USE lib_mpp ! MPP library 26 USE wrk_nemo ! work arrays 27 USE in_out_manager ! I/O manager 28 USE prtctl ! Print control 19 USE phycst ! Physical constant 20 USE oce , ONLY : snwice_mass, snwice_mass_b 21 USE par_oce ! Ocean parameters 22 USE dom_oce ! Ocean domain 23 USE sbc_oce ! Surface boundary condition: ocean fields 24 USE sbc_ice ! Surface boundary condition: ice fields 29 25 #if defined key_lim3 30 USE ice 31 USE dom_ice 32 USE limitd_me 26 USE ice ! LIM-3: ice variables 27 USE dom_ice ! LIM-3: ice domain 28 USE limitd_me ! LIM-3: 33 29 #else 34 USE ice_2 ! LIM2: ice variables35 USE dom_ice_2 ! LIM2: ice domain30 USE ice_2 ! LIM-2: ice variables 31 USE dom_ice_2 ! LIM-2: ice domain 36 32 #endif 33 USE lbclnk ! Lateral Boundary Condition / MPP link 34 USE lib_mpp ! MPP library 35 USE wrk_nemo ! work arrays 36 USE in_out_manager ! I/O manager 37 USE prtctl ! Print control 38 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 37 39 38 40 IMPLICIT NONE … … 47 49 # include "vectopt_loop_substitute.h90" 48 50 !!---------------------------------------------------------------------- 49 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)51 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 50 52 !! $Id$ 51 53 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 124 126 REAL(wp) :: zindb ! ice (1) or not (0) 125 127 REAL(wp) :: zdummy ! dummy argument 128 REAL(wp) :: zintb, zintn ! dummy argument 126 129 127 130 REAL(wp), POINTER, DIMENSION(:,:) :: zpresh ! temporary array for ice strength … … 144 147 REAL(wp), POINTER, DIMENSION(:,:) :: zs12 ! Non-diagonal stress tensor component zs12 145 148 REAL(wp), POINTER, DIMENSION(:,:) :: zu_ice, zv_ice, zresr ! Local error on velocity 149 REAL(wp), POINTER, DIMENSION(:,:) :: zpice ! array used for the calculation of ice surface slope: 150 ! ocean surface (ssh_m) if ice is not embedded 151 ! ice top surface if ice is embedded 146 152 147 153 !!------------------------------------------------------------------- … … 150 156 CALL wrk_alloc( jpi,jpj, zc1 , u_oce1, u_oce2, u_ice2, zusw , v_oce1 , v_oce2, v_ice1 ) 151 157 CALL wrk_alloc( jpi,jpj, zf1 , deltat, zu_ice, zf2 , deltac, zv_ice , zdd , zdt , zds ) 152 CALL wrk_alloc( jpi,jpj, zdd , zdt , zds , zs1 , zs2 , zs12 , zresr 158 CALL wrk_alloc( jpi,jpj, zdd , zdt , zds , zs1 , zs2 , zs12 , zresr , zpice ) 153 159 154 160 #if defined key_lim2 && ! defined key_lim2_vp … … 231 237 ! v_oce2: ocean v component on v points 232 238 239 IF( nn_ice_embd == 2 ) THEN !== embedded sea ice: compute representative ice top surface ==! 240 ! 241 ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[n/nn_fsbc], n=0,nn_fsbc-1} 242 ! = (1/nn_fsbc)^2 * {SUM[n], n=0,nn_fsbc-1} 243 zintn = REAL( nn_fsbc - 1 ) / REAL( nn_fsbc ) * 0.5_wp 244 ! 245 ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[1-n/nn_fsbc], n=0,nn_fsbc-1} 246 ! = (1/nn_fsbc)^2 * (nn_fsbc^2 - {SUM[n], n=0,nn_fsbc-1}) 247 zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp 248 ! 249 zpice(:,:) = ssh_m(:,:) + ( zintn * snwice_mass(:,:) + zintb * snwice_mass_b(:,:) ) * r1_rau0 250 ! 251 ELSE !== non-embedded sea ice: use ocean surface for slope calculation ==! 252 zpice(:,:) = ssh_m(:,:) 253 ENDIF 254 233 255 DO jj = k_j1+1, k_jpj-1 234 256 DO ji = fs_2, fs_jpim1 … … 273 295 ! include it later 274 296 275 zdsshx = ( ssh_m(ji+1,jj) - ssh_m(ji,jj) ) / e1u(ji,jj)276 zdsshy = ( ssh_m(ji,jj+1) - ssh_m(ji,jj) ) / e2v(ji,jj)297 zdsshx = ( zpice(ji+1,jj) - zpice(ji,jj) ) / e1u(ji,jj) 298 zdsshy = ( zpice(ji,jj+1) - zpice(ji,jj) ) / e2v(ji,jj) 277 299 278 300 za1ct(ji,jj) = ztagnx - zmass1(ji,jj) * grav * zdsshx … … 746 768 CALL wrk_dealloc( jpi,jpj, zc1 , u_oce1, u_oce2, u_ice2, zusw , v_oce1 , v_oce2, v_ice1 ) 747 769 CALL wrk_dealloc( jpi,jpj, zf1 , deltat, zu_ice, zf2 , deltac, zv_ice , zdd , zdt , zds ) 748 CALL wrk_dealloc( jpi,jpj, zdd , zdt , zds , zs1 , zs2 , zs12 , zresr 770 CALL wrk_dealloc( jpi,jpj, zdd , zdt , zds , zs1 , zs2 , zs12 , zresr , zpice ) 749 771 750 772 END SUBROUTINE lim_rhg -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limrst.F90
r3294 r3625 12 12 !! 'key_lim3' : LIM sea-ice model 13 13 !!---------------------------------------------------------------------- 14 !! lim_rst_opn : open ice restart file 15 !! lim_rst_write : write of the restart file 16 !! lim_rst_read : read the restart file 17 !!---------------------------------------------------------------------- 18 USE ice ! sea-ice variables 19 USE par_ice ! sea-ice parameters 20 USE dom_oce ! ocean domain 21 USE sbc_oce ! Surface boundary condition: ocean fields 22 USE sbc_ice ! Surface boundary condition: ice fields 23 USE in_out_manager ! I/O manager 24 USE iom ! I/O library 25 USE lib_mpp ! MPP library 26 USE wrk_nemo ! work arrays 14 !! lim_rst_opn : open ice restart file 15 !! lim_rst_write : write of the restart file 16 !! lim_rst_read : read the restart file 17 !!---------------------------------------------------------------------- 18 USE ice ! sea-ice variables 19 USE par_ice ! sea-ice parameters 20 USE dom_oce ! ocean domain 21 USE sbc_oce ! Surface boundary condition: ocean fields 22 USE sbc_ice ! Surface boundary condition: ice fields 23 USE in_out_manager ! I/O manager 24 USE iom ! I/O library 25 USE lib_mpp ! MPP library 26 USE wrk_nemo ! work arrays 27 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 27 28 28 29 IMPLICIT NONE … … 37 38 38 39 !!---------------------------------------------------------------------- 39 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)40 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 40 41 !! $Id$ 41 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 402 403 zsmax = 4.5_wp 403 404 zsmin = 3.5_wp 404 IF( sm_i(ji,jj,jl) .LT.zsmin ) THEN405 IF( sm_i(ji,jj,jl) < zsmin ) THEN 405 406 zalpha = 1._wp 406 ELSEIF( sm_i(ji,jj,jl) .LT.zsmax ) THEN407 ELSEIF( sm_i(ji,jj,jl) < zsmax ) THEN 407 408 zalpha = sm_i(ji,jj,jl) / ( zsmin - zsmax ) + zsmax / ( zsmax - zsmin ) 408 409 ELSE -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90
r3294 r3625 9 9 !! 3.3 ! 2010-05 (G. Madec) decrease ocean & ice reference salinities in the Baltic sea 10 10 !! ! + simplification of the ice-ocean stress calculation 11 !! 4.0 ! 2011-02 (G. Madec) dynamical allocation 11 !! 3.4 ! 2011-02 (G. Madec) dynamical allocation 12 !! 3.5 ! 2012-10 (A. Coward, G. Madec) salt fluxes ; ice+snow mass 12 13 !!---------------------------------------------------------------------- 13 14 #if defined key_lim3 … … 34 35 USE prtctl ! Print control 35 36 USE cpl_oasis3, ONLY : lk_cpl 37 USE oce, ONLY : sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass, sshu_b, sshv_b, sshu_n, sshv_n, sshf_n 38 USE dom_ice, ONLY : tms 39 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 36 40 37 41 IMPLICIT NONE … … 42 46 PUBLIC lim_sbc_tau ! called by sbc_ice_lim 43 47 44 REAL(wp) :: r1_rdtice ! = 1. / rdt_ice45 48 REAL(wp) :: epsi16 = 1.e-16_wp ! constant values 46 49 REAL(wp) :: rzero = 0._wp … … 54 57 # include "vectopt_loop_substitute.h90" 55 58 !!---------------------------------------------------------------------- 56 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)59 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 57 60 !! $Id$ 58 61 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 86 89 !! - qns : sea heat flux: non solar 87 90 !! - emp : freshwater budget: volume flux 88 !! - emps : freshwater budget: concentration/dillution91 !! - sfx : salt flux 89 92 !! - fr_i : ice fraction 90 93 !! - tn_ice : sea-ice surface temperature … … 97 100 ! 98 101 INTEGER :: ji, jj ! dummy loop indices 99 INTEGER :: ierr 100 INTEGER :: if vt, i1mfr, idfr ! some switches101 INTEGER :: iflt, ial, iadv, ifral, ifrdv102 REAL(wp) :: z inda, zfons, zpme ! local scalars103 REAL(wp) , POINTER, DIMENSION(:,:) :: zfcm1 , zfcm2 ! solar/non solar heat fluxes104 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace102 INTEGER :: ierr, ifvt, i1mfr, idfr ! local integer 103 INTEGER :: iflt, ial , iadv , ifral, ifrdv ! - - 104 REAL(wp) :: zinda, zemp, zemp_snow, zfmm ! local scalars 105 REAL(wp) :: zemp_snw ! - - 106 REAL(wp) :: zfcm1 , zfcm2 ! - - 107 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace 105 108 !!--------------------------------------------------------------------- 106 109 107 CALL wrk_alloc( jpi, jpj, zfcm1 , zfcm2 )108 110 IF( lk_cpl ) CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp ) 109 111 … … 139 141 140 142 ! computation the solar flux at ocean surface 141 zfcm1 (ji,jj) = pfrld(ji,jj) * qsr(ji,jj) + ( 1.- pfrld(ji,jj) ) * fstric(ji,jj)143 zfcm1 = pfrld(ji,jj) * qsr(ji,jj) + ( 1._wp - pfrld(ji,jj) ) * fstric(ji,jj) 142 144 ! fstric Solar flux transmitted trough the ice 143 145 ! qsr Net short wave heat flux on free ocean … … 146 148 147 149 ! computation the non solar heat flux at ocean surface 148 zfcm2(ji,jj) = - zfcm1(ji,jj) & 149 & + iflt * ( fscmbq(ji,jj) ) & ! total abl -> fscmbq is given to the ocean 150 ! fscmbq and ffltbif are obsolete 151 ! & + iflt * ffltbif(ji,jj) !!! only if one category is used 152 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice & 153 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice & 154 & + fhmec(ji,jj) & ! new contribution due to snow melt in ridging!! 155 & + fheat_rpo(ji,jj) & ! contribution from ridge formation 156 & + fheat_res(ji,jj) 157 ! fscmbq Part of the solar radiation transmitted through the ice and going to the ocean 158 ! computed in limthd_zdf.F90 159 ! ffltbif Total heat content of the ice (brine pockets+ice) / delta_t 150 zfcm2 = - zfcm1 & ! ??? 151 & + iflt * fscmbq(ji,jj) & ! total ablation: heat given to the ocean 152 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice & 153 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice & 154 & + fhmec(ji,jj) & ! snow melt when ridging 155 & + fheat_mec(ji,jj) & ! ridge formation 156 & + fheat_res(ji,jj) ! residual heat flux 160 157 ! qcmif Energy needed to bring the ocean surface layer until its freezing (ok) 161 158 ! qldif heat balance of the lead (or of the open ocean) 162 ! qfvbq i think this is wrong! 163 ! ---> Array used to store energy in case of total lateral ablation 164 ! qfvbq latent heat uptake/release after accretion/ablation 165 ! qdtcn Energy from the turbulent oceanic heat flux heat flux coming in the lead 166 167 IF ( num_sal == 2 ) zfcm2(ji,jj) = zfcm2(ji,jj) + & 168 fhbri(ji,jj) ! new contribution due to brine drainage 169 170 ! bottom radiative component is sent to the computation of the 171 ! oceanic heat flux 172 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) 159 ! qfvbq latent heat uptake/release after accretion/ablation 160 ! qdtcn Energy from the turbulent oceanic heat flux heat flux coming in the lead 161 162 IF( num_sal == 2 ) zfcm2 = zfcm2 + fhbri(ji,jj) ! add contribution due to brine drainage 163 164 ! bottom radiative component is sent to the computation of the oceanic heat flux 165 fsbbq(ji,jj) = ( 1._wp - ( ifvt + iflt ) ) * fscmbq(ji,jj) 173 166 174 167 ! used to compute the oceanic heat flux at the next time step 175 qsr(ji,jj) = zfcm1 (ji,jj)! solar heat flux176 qns(ji,jj) = zfcm2 (ji,jj)- fdtcn(ji,jj) ! non solar heat flux168 qsr(ji,jj) = zfcm1 ! solar heat flux 169 qns(ji,jj) = zfcm2 - fdtcn(ji,jj) ! non solar heat flux 177 170 ! ! fdtcn : turbulent oceanic heat flux 178 171 179 172 !!gm this IF prevents the vertorisation of the whole loop 180 173 IF ( ( ji == jiindx ) .AND. ( jj == jjindx) ) THEN 181 174 WRITE(numout,*) ' lim_sbc : heat fluxes ' 182 175 WRITE(numout,*) ' qsr : ', qsr(jiindx,jjindx) 183 WRITE(numout,*) ' zfcm1 : ', zfcm1(jiindx,jjindx)184 176 WRITE(numout,*) ' pfrld : ', pfrld(jiindx,jjindx) 185 177 WRITE(numout,*) ' fstric : ', fstric (jiindx,jjindx) 186 178 WRITE(numout,*) 187 179 WRITE(numout,*) ' qns : ', qns(jiindx,jjindx) 188 WRITE(numout,*) ' zfcm2 : ', zfcm2(jiindx,jjindx) 189 WRITE(numout,*) ' zfcm1 : ', zfcm1(jiindx,jjindx) 180 WRITE(numout,*) ' fdtcn : ', fdtcn(jiindx,jjindx) 190 181 WRITE(numout,*) ' ifral : ', ifral 191 182 WRITE(numout,*) ' ial : ', ial … … 202 193 WRITE(numout,*) ' fdtcn : ', fdtcn(jiindx,jjindx) 203 194 WRITE(numout,*) ' fhmec : ', fhmec(jiindx,jjindx) 204 WRITE(numout,*) ' fheat_ rpo : ', fheat_rpo(jiindx,jjindx)195 WRITE(numout,*) ' fheat_mec : ', fheat_mec(jiindx,jjindx) 205 196 WRITE(numout,*) ' fhbri : ', fhbri(jiindx,jjindx) 206 197 WRITE(numout,*) ' fheat_res : ', fheat_res(jiindx,jjindx) 207 198 ENDIF 208 199 !!gm end 209 200 END DO 210 201 END DO … … 227 218 228 219 ! computing freshwater exchanges at the ice/ocean interface 229 zpme = - emp(ji,jj) * ( 1.0 - at_i(ji,jj) ) & ! evaporation over oceanic fraction 230 & + tprecip(ji,jj) * at_i(ji,jj) & ! all precipitation reach the ocean 231 & - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) ) & ! except solid precip intercepted by sea-ice 232 & - rdmsnif(ji,jj) * r1_rdtice & ! freshwaterflux due to snow melting 233 & + fmmec(ji,jj) ! snow falling when ridging 234 235 236 ! computing salt exchanges at the ice/ocean interface 237 ! sice should be the same as computed with the ice model 238 zfons = ( soce_0(ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 239 ! SOCE 240 zfons = ( sss_m (ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 241 242 !CT useless ! salt flux for constant salinity 243 !CT useless fsalt(ji,jj) = zfons / ( sss_m(ji,jj) + epsi16 ) + fsalt_res(ji,jj) 244 ! salt flux for variable salinity 245 zinda = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) ) 246 ! correcting brine and salt fluxes 247 fsbri(ji,jj) = zinda*fsbri(ji,jj) 248 ! converting the salt fluxes from ice to a freshwater flux from ocean 249 fsalt_res(ji,jj) = fsalt_res(ji,jj) / ( sss_m(ji,jj) + epsi16 ) 250 fseqv(ji,jj) = fseqv(ji,jj) / ( sss_m(ji,jj) + epsi16 ) 251 fsbri(ji,jj) = fsbri(ji,jj) / ( sss_m(ji,jj) + epsi16 ) 252 fsalt_rpo(ji,jj) = fsalt_rpo(ji,jj) / ( sss_m(ji,jj) + epsi16 ) 253 254 ! freshwater mass exchange (positive to the ice, negative for the ocean ?) 255 ! actually it's a salt flux (so it's minus freshwater flux) 256 ! if sea ice grows, zfons is positive, fsalt also 257 ! POSITIVE SALT FLUX FROM THE ICE TO THE OCEAN 258 ! POSITIVE FRESHWATER FLUX FROM THE OCEAN TO THE ICE [kg.m-2.s-1] 259 260 emp(ji,jj) = - zpme 220 zemp = emp(ji,jj) * ( 1.0 - at_i(ji,jj) ) & ! evaporation over oceanic fraction 221 & - tprecip(ji,jj) * at_i(ji,jj) & ! all precipitation reach the ocean 222 & + sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) ) & ! except solid precip intercepted by sea-ice 223 & - fmmec(ji,jj) ! snow falling when ridging 224 225 ! mass flux at the ocean/ice interface (sea ice fraction) 226 zemp_snw = rdm_snw(ji,jj) * r1_rdtice ! snow melting = pure water that enters the ocean 227 zfmm = rdm_ice(ji,jj) * r1_rdtice ! Freezing minus mesting 228 229 emp(ji,jj) = zemp + zemp_snw + zfmm ! mass flux + F/M mass flux (always ice/ocean mass exchange) 230 231 ! correcting brine salt fluxes (zinda = 1 if pfrld=1 , =0 otherwise) 232 zinda = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) ) 233 sfx_bri(ji,jj) = zinda * sfx_bri(ji,jj) 261 234 END DO 262 235 END DO 263 236 237 !------------------------------------------! 238 ! salt flux at the ocean surface ! 239 !------------------------------------------! 240 264 241 IF( num_sal == 2 ) THEN ! variable ice salinity: brine drainage included in the salt flux 265 emps(:,:) = fsbri(:,:) + fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:) + emp(:,:)242 sfx(:,:) = sfx_thd(:,:) + sfx_res(:,:) + sfx_mec(:,:) + sfx_bri(:,:) 266 243 ELSE ! constant ice salinity: 267 emps(:,:) = fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:) + emp(:,:) 244 sfx(:,:) = sfx_thd(:,:) + sfx_res(:,:) + sfx_mec(:,:) 245 ENDIF 246 !-----------------------------------------------! 247 ! mass of snow and ice per unit area ! 248 !-----------------------------------------------! 249 IF( nn_ice_embd /= 0 ) THEN ! embedded sea-ice (mass required) 250 snwice_mass_b(:,:) = snwice_mass(:,:) ! save mass from the previous ice time step 251 ! ! new mass per unit area 252 snwice_mass (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:) ) 253 ! ! time evolution of snow+ice mass 254 snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_rdtice 268 255 ENDIF 269 256 … … 285 272 IF(ln_ctl) THEN 286 273 CALL prt_ctl( tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ' ) 287 CALL prt_ctl( tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2= emps, clinfo2=' emps: ' )274 CALL prt_ctl( tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=sfx , clinfo2=' sfx : ' ) 288 275 CALL prt_ctl( tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ' ) 289 276 CALL prt_ctl( tab3d_1=tn_ice, clinfo1=' lim_sbc: tn_ice : ', kdim=jpl ) 290 277 ENDIF 291 278 ! 292 CALL wrk_dealloc( jpi, jpj, zfcm1 , zfcm2 )293 279 IF( lk_cpl ) CALL wrk_dealloc( jpi, jpj, jpl, zalb, zalbp ) 294 280 ! … … 383 369 !!------------------------------------------------------------------- 384 370 ! 371 INTEGER :: ji, jj ! dummy loop indices 372 REAL(wp) :: zcoefu, zcoefv, zcoeff ! local scalar 385 373 IF(lwp) WRITE(numout,*) 386 374 IF(lwp) WRITE(numout,*) 'lim_sbc_init : LIM-3 sea-ice - surface boundary condition' … … 389 377 ! ! allocate lim_sbc array 390 378 IF( lim_sbc_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'lim_sbc_init : unable to allocate standard arrays' ) 391 !392 r1_rdtice = 1. / rdt_ice393 379 ! 394 380 soce_0(:,:) = soce ! constant SSS and ice salinity used in levitating sea-ice case … … 402 388 END WHERE 403 389 ENDIF 390 ! ! embedded sea ice 391 IF( nn_ice_embd /= 0 ) THEN ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass 392 snwice_mass (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:) ) 393 snwice_mass_b(:,:) = snwice_mass(:,:) 394 ELSE 395 snwice_mass (:,:) = 0.0_wp ! no mass exchanges 396 snwice_mass_b(:,:) = 0.0_wp ! no mass exchanges 397 ENDIF 398 IF( nn_ice_embd == 2 .AND. & ! full embedment (case 2) & no restart 399 & .NOT. ln_rstart ) THEN ! deplete the initial ssh below sea-ice area 400 sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0 401 sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0 402 ! 403 ! Note: Changed the initial values of sshb and sshn=> need to recompute ssh[u,v,f]_[b,n] 404 ! which were previously set in domvvl 405 IF ( lk_vvl ) THEN ! Is this necessary? embd 2 should be restricted to vvl only??? 406 DO jj = 1, jpjm1 407 DO ji = 1, jpim1 ! caution: use of Vector Opt. not possible 408 zcoefu = 0.5 * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) ) 409 zcoefv = 0.5 * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) ) 410 zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1) 411 sshu_b(ji,jj) = zcoefu * ( e1t(ji ,jj) * e2t(ji ,jj) * sshb(ji ,jj) & 412 & + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) ) 413 sshv_b(ji,jj) = zcoefv * ( e1t(ji,jj ) * e2t(ji,jj ) * sshb(ji,jj ) & 414 & + e1t(ji,jj+1) * e2t(ji,jj+1) * sshb(ji,jj+1) ) 415 sshu_n(ji,jj) = zcoefu * ( e1t(ji ,jj) * e2t(ji ,jj) * sshn(ji ,jj) & 416 & + e1t(ji+1,jj) * e2t(ji+1,jj) * sshn(ji+1,jj) ) 417 sshv_n(ji,jj) = zcoefv * ( e1t(ji,jj ) * e2t(ji,jj ) * sshn(ji,jj ) & 418 & + e1t(ji,jj+1) * e2t(ji,jj+1) * sshn(ji,jj+1) ) 419 END DO 420 END DO 421 CALL lbc_lnk( sshu_b, 'U', 1. ) ; CALL lbc_lnk( sshu_n, 'U', 1. ) 422 CALL lbc_lnk( sshv_b, 'V', 1. ) ; CALL lbc_lnk( sshv_n, 'V', 1. ) 423 DO jj = 1, jpjm1 424 DO ji = 1, jpim1 ! NO Vector Opt. 425 sshf_n(ji,jj) = 0.5 * umask(ji,jj,1) * umask(ji,jj+1,1) & 426 & / ( e1f(ji,jj ) * e2f(ji,jj ) ) & 427 & * ( e1u(ji,jj ) * e2u(ji,jj ) * sshu_n(ji,jj ) & 428 & + e1u(ji,jj+1) * e2u(ji,jj+1) * sshu_n(ji,jj+1) ) 429 END DO 430 END DO 431 CALL lbc_lnk( sshf_n, 'F', 1. ) 432 ENDIF 433 ENDIF 404 434 ! 405 435 END SUBROUTINE lim_sbc_init -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limtab.F90
r2715 r3625 2 2 !!====================================================================== 3 3 !! *** MODULE limtab *** 4 !! LIM : transform 1D (2D) array to a 2D (1D) table4 !! LIM ice model : transform 1D (2D) array to a 2D (1D) table 5 5 !!====================================================================== 6 6 #if defined key_lim3 … … 20 20 21 21 !!---------------------------------------------------------------------- 22 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2010)22 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2010) 23 23 !! $Id$ 24 24 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r3294 r3625 8 8 !! 3.0 ! 2005-11 (M. Vancoppenolle) LIM-3 : Multi-layer thermodynamics + salinity variations 9 9 !! - ! 2007-04 (M. Vancoppenolle) add lim_thd_glohec, lim_thd_con_dh and lim_thd_con_dif 10 !! 3.2 ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in rdm snif10 !! 3.2 ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in rdm_snw 11 11 !! 3.3 ! 2010-11 (G. Madec) corrected snow melting heat (due to factor betas) 12 12 !! 4.0 ! 2011-02 (G. Madec) dynamical allocation … … 16 16 !! 'key_lim3' LIM3 sea-ice model 17 17 !!---------------------------------------------------------------------- 18 !! lim_thd 19 !! lim_thd_init 18 !! lim_thd : thermodynamic of sea ice 19 !! lim_thd_init : initialisation of sea-ice thermodynamic 20 20 !!---------------------------------------------------------------------- 21 USE phycst ! physical constants 22 USE dom_oce ! ocean space and time domain variables 23 USE ice ! LIM: sea-ice variables 24 USE par_ice ! LIM: sea-ice parameters 25 USE sbc_oce ! Surface boundary condition: ocean fields 26 USE sbc_ice ! Surface boundary condition: ice fields 27 USE thd_ice ! LIM thermodynamic sea-ice variables 28 USE dom_ice ! LIM sea-ice domain 29 USE domvvl ! domain: variable volume level 30 USE limthd_dif ! LIM: thermodynamics, vertical diffusion 31 USE limthd_dh ! LIM: thermodynamics, ice and snow thickness variation 32 USE limthd_sal ! LIM: thermodynamics, ice salinity 33 USE limthd_ent ! LIM: thermodynamics, ice enthalpy redistribution 34 USE limtab ! LIM: 1D <==> 2D transformation 35 USE limvar ! LIM: sea-ice variables 36 USE lbclnk ! lateral boundary condition - MPP links 37 USE lib_mpp ! MPP library 38 USE wrk_nemo ! work arrays 39 USE in_out_manager ! I/O manager 40 USE prtctl ! Print control 21 USE phycst ! physical constants 22 USE dom_oce ! ocean space and time domain variables 23 USE ice ! LIM: sea-ice variables 24 USE par_ice ! LIM: sea-ice parameters 25 USE sbc_oce ! Surface boundary condition: ocean fields 26 USE sbc_ice ! Surface boundary condition: ice fields 27 USE thd_ice ! LIM thermodynamic sea-ice variables 28 USE dom_ice ! LIM sea-ice domain 29 USE domvvl ! domain: variable volume level 30 USE limthd_dif ! LIM: thermodynamics, vertical diffusion 31 USE limthd_dh ! LIM: thermodynamics, ice and snow thickness variation 32 USE limthd_sal ! LIM: thermodynamics, ice salinity 33 USE limthd_ent ! LIM: thermodynamics, ice enthalpy redistribution 34 USE limtab ! LIM: 1D <==> 2D transformation 35 USE limvar ! LIM: sea-ice variables 36 USE lbclnk ! lateral boundary condition - MPP links 37 USE lib_mpp ! MPP library 38 USE wrk_nemo ! work arrays 39 USE in_out_manager ! I/O manager 40 USE prtctl ! Print control 41 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 41 42 42 43 IMPLICIT NONE … … 110 111 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_i(ji,jj,jl) , epsi06 ) ) * nlay_i 111 112 !0 if no ice and 1 if yes 112 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_i(ji,jj,jl) ))113 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - ht_i(ji,jj,jl) ) ) 113 114 !convert units ! very important that this line is here 114 115 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * unit_fac * zindb … … 122 123 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_s(ji,jj,jl) , epsi06 ) ) * nlay_s 123 124 !0 if no ice and 1 if yes 124 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_s(ji,jj,jl) ))125 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - ht_s(ji,jj,jl) ) ) 125 126 !convert units ! very important that this line is here 126 127 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * unit_fac * zindb … … 140 141 ffltbif(:,:) = 0.e0 ! linked with fstric 141 142 qfvbq (:,:) = 0.e0 ! linked with fstric 142 rdm snif(:,:) = 0.e0 ! variation of snow mass per unit area143 rdm icif(:,:) = 0.e0 ! variation of ice mass per unit area143 rdm_snw(:,:) = 0.e0 ! variation of snow mass per unit area 144 rdm_ice(:,:) = 0.e0 ! variation of ice mass per unit area 144 145 hicifp (:,:) = 0.e0 ! daily thermodynamic ice production. 145 fsbri(:,:) = 0.e0 ! brine flux contribution to salt flux to the ocean146 sfx_bri(:,:) = 0.e0 ! brine flux contribution to salt flux to the ocean 146 147 fhbri (:,:) = 0.e0 ! brine flux contribution to heat flux to the ocean 147 fseqv(:,:) = 0.e0 ! equivalent salt flux to the ocean due to ice/growth decay148 sfx_thd(:,:) = 0.e0 ! equivalent salt flux to the ocean due to ice/growth decay 148 149 149 150 !----------------------------------- … … 273 274 CALL tab_2d_1d( nbpb, fr2_i0_1d (1:nbpb), fr2_i0 , jpi, jpj, npb(1:nbpb) ) 274 275 CALL tab_2d_1d( nbpb, qnsr_ice_1d(1:nbpb), qns_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 275 276 276 #if ! defined key_coupled 277 CALL tab_2d_1d( nbpb, qla_ice_1d (1:nbpb), qla_ice(:,:,jl) 278 CALL tab_2d_1d( nbpb, dqla_ice_1d(1:nbpb), dqla_ice(:,:,jl) 277 CALL tab_2d_1d( nbpb, qla_ice_1d (1:nbpb), qla_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 278 CALL tab_2d_1d( nbpb, dqla_ice_1d(1:nbpb), dqla_ice(:,:,jl), jpi, jpj, npb(1:nbpb) ) 279 279 #endif 280 281 CALL tab_2d_1d( nbpb, dqns_ice_1d(1:nbpb), dqns_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 282 CALL tab_2d_1d( nbpb, t_bo_b (1:nbpb), t_bo , jpi, jpj, npb(1:nbpb) ) 283 CALL tab_2d_1d( nbpb, sprecip_1d (1:nbpb), sprecip , jpi, jpj, npb(1:nbpb) ) 284 CALL tab_2d_1d( nbpb, fbif_1d (1:nbpb), fbif , jpi, jpj, npb(1:nbpb) ) 285 CALL tab_2d_1d( nbpb, qldif_1d (1:nbpb), qldif , jpi, jpj, npb(1:nbpb) ) 286 CALL tab_2d_1d( nbpb, rdmicif_1d (1:nbpb), rdmicif , jpi, jpj, npb(1:nbpb) ) 287 CALL tab_2d_1d( nbpb, rdmsnif_1d (1:nbpb), rdmsnif , jpi, jpj, npb(1:nbpb) ) 288 CALL tab_2d_1d( nbpb, dmgwi_1d (1:nbpb), dmgwi , jpi, jpj, npb(1:nbpb) ) 289 CALL tab_2d_1d( nbpb, qlbbq_1d (1:nbpb), zqlbsbq , jpi, jpj, npb(1:nbpb) ) 290 291 CALL tab_2d_1d( nbpb, fseqv_1d (1:nbpb), fseqv , jpi, jpj, npb(1:nbpb) ) 292 CALL tab_2d_1d( nbpb, fsbri_1d (1:nbpb), fsbri , jpi, jpj, npb(1:nbpb) ) 293 CALL tab_2d_1d( nbpb, fhbri_1d (1:nbpb), fhbri , jpi, jpj, npb(1:nbpb) ) 294 CALL tab_2d_1d( nbpb, fstbif_1d (1:nbpb), fstric , jpi, jpj, npb(1:nbpb) ) 295 CALL tab_2d_1d( nbpb, qfvbq_1d (1:nbpb), qfvbq , jpi, jpj, npb(1:nbpb) ) 280 CALL tab_2d_1d( nbpb, dqns_ice_1d(1:nbpb), dqns_ice(:,:,jl), jpi, jpj, npb(1:nbpb) ) 281 CALL tab_2d_1d( nbpb, t_bo_b (1:nbpb), t_bo , jpi, jpj, npb(1:nbpb) ) 282 CALL tab_2d_1d( nbpb, sprecip_1d (1:nbpb), sprecip , jpi, jpj, npb(1:nbpb) ) 283 CALL tab_2d_1d( nbpb, fbif_1d (1:nbpb), fbif , jpi, jpj, npb(1:nbpb) ) 284 CALL tab_2d_1d( nbpb, qldif_1d (1:nbpb), qldif , jpi, jpj, npb(1:nbpb) ) 285 CALL tab_2d_1d( nbpb, rdm_ice_1d (1:nbpb), rdm_ice , jpi, jpj, npb(1:nbpb) ) 286 CALL tab_2d_1d( nbpb, rdm_snw_1d (1:nbpb), rdm_snw , jpi, jpj, npb(1:nbpb) ) 287 CALL tab_2d_1d( nbpb, dmgwi_1d (1:nbpb), dmgwi , jpi, jpj, npb(1:nbpb) ) 288 CALL tab_2d_1d( nbpb, qlbbq_1d (1:nbpb), zqlbsbq , jpi, jpj, npb(1:nbpb) ) 289 290 CALL tab_2d_1d( nbpb, sfx_thd_1d (1:nbpb), sfx_thd , jpi, jpj, npb(1:nbpb) ) 291 CALL tab_2d_1d( nbpb, sfx_bri_1d (1:nbpb), sfx_bri , jpi, jpj, npb(1:nbpb) ) 292 CALL tab_2d_1d( nbpb, fhbri_1d (1:nbpb), fhbri , jpi, jpj, npb(1:nbpb) ) 293 CALL tab_2d_1d( nbpb, fstbif_1d (1:nbpb), fstric , jpi, jpj, npb(1:nbpb) ) 294 CALL tab_2d_1d( nbpb, qfvbq_1d (1:nbpb), qfvbq , jpi, jpj, npb(1:nbpb) ) 296 295 297 296 !-------------------------------- … … 331 330 !-------------------------------- 332 331 333 CALL tab_1d_2d( nbpb, at_i , npb, at_i_b(1:nbpb), jpi, jpj ) 334 CALL tab_1d_2d( nbpb, ht_i(:,:,jl), npb, ht_i_b(1:nbpb), jpi, jpj ) 335 CALL tab_1d_2d( nbpb, ht_s(:,:,jl), npb, ht_s_b(1:nbpb), jpi, jpj ) 336 CALL tab_1d_2d( nbpb, a_i (:,:,jl), npb, a_i_b(1:nbpb) , jpi, jpj ) 337 CALL tab_1d_2d( nbpb, t_su(:,:,jl), npb, t_su_b(1:nbpb), jpi, jpj ) 338 CALL tab_1d_2d( nbpb, sm_i(:,:,jl), npb, sm_i_b(1:nbpb), jpi, jpj ) 339 332 CALL tab_1d_2d( nbpb, at_i , npb, at_i_b (1:nbpb) , jpi, jpj ) 333 CALL tab_1d_2d( nbpb, ht_i(:,:,jl) , npb, ht_i_b (1:nbpb) , jpi, jpj ) 334 CALL tab_1d_2d( nbpb, ht_s(:,:,jl) , npb, ht_s_b (1:nbpb) , jpi, jpj ) 335 CALL tab_1d_2d( nbpb, a_i (:,:,jl) , npb, a_i_b (1:nbpb) , jpi, jpj ) 336 CALL tab_1d_2d( nbpb, t_su(:,:,jl) , npb, t_su_b (1:nbpb) , jpi, jpj ) 337 CALL tab_1d_2d( nbpb, sm_i(:,:,jl) , npb, sm_i_b (1:nbpb) , jpi, jpj ) 340 338 DO jk = 1, nlay_s 341 CALL tab_1d_2d( nbpb, t_s(:,:,jk,jl), npb, t_s_b (1:nbpb,jk), jpi, jpj)342 CALL tab_1d_2d( nbpb, e_s(:,:,jk,jl), npb, q_s_b (1:nbpb,jk), jpi, jpj)339 CALL tab_1d_2d( nbpb, t_s(:,:,jk,jl), npb, t_s_b (1:nbpb,jk), jpi, jpj) 340 CALL tab_1d_2d( nbpb, e_s(:,:,jk,jl), npb, q_s_b (1:nbpb,jk), jpi, jpj) 343 341 END DO 344 345 342 DO jk = 1, nlay_i 346 CALL tab_1d_2d( nbpb, t_i(:,:,jk,jl), npb, t_i_b (1:nbpb,jk), jpi, jpj)347 CALL tab_1d_2d( nbpb, e_i(:,:,jk,jl), npb, q_i_b (1:nbpb,jk), jpi, jpj)348 CALL tab_1d_2d( nbpb, s_i(:,:,jk,jl), npb, s_i_b (1:nbpb,jk), jpi, jpj)343 CALL tab_1d_2d( nbpb, t_i(:,:,jk,jl), npb, t_i_b (1:nbpb,jk), jpi, jpj) 344 CALL tab_1d_2d( nbpb, e_i(:,:,jk,jl), npb, q_i_b (1:nbpb,jk), jpi, jpj) 345 CALL tab_1d_2d( nbpb, s_i(:,:,jk,jl), npb, s_i_b (1:nbpb,jk), jpi, jpj) 349 346 END DO 350 351 CALL tab_1d_2d( nbpb, fstric , npb, fstbif_1d (1:nbpb), jpi, jpj ) 352 CALL tab_1d_2d( nbpb, qldif , npb, qldif_1d (1:nbpb), jpi, jpj ) 353 CALL tab_1d_2d( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb), jpi, jpj ) 354 CALL tab_1d_2d( nbpb, rdmicif, npb, rdmicif_1d(1:nbpb), jpi, jpj ) 355 CALL tab_1d_2d( nbpb, rdmsnif, npb, rdmsnif_1d(1:nbpb), jpi, jpj ) 356 CALL tab_1d_2d( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb), jpi, jpj ) 357 CALL tab_1d_2d( nbpb, rdvosif, npb, dvsbq_1d (1:nbpb), jpi, jpj ) 358 CALL tab_1d_2d( nbpb, rdvobif, npb, dvbbq_1d (1:nbpb), jpi, jpj ) 359 CALL tab_1d_2d( nbpb, fdvolif, npb, dvlbq_1d (1:nbpb), jpi, jpj ) 360 CALL tab_1d_2d( nbpb, rdvonif, npb, dvnbq_1d (1:nbpb), jpi, jpj ) 361 CALL tab_1d_2d( nbpb, fseqv , npb, fseqv_1d (1:nbpb), jpi, jpj ) 347 CALL tab_1d_2d( nbpb, fstric , npb, fstbif_1d (1:nbpb) , jpi, jpj ) 348 CALL tab_1d_2d( nbpb, qldif , npb, qldif_1d (1:nbpb) , jpi, jpj ) 349 CALL tab_1d_2d( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb) , jpi, jpj ) 350 CALL tab_1d_2d( nbpb, rdm_ice , npb, rdm_ice_1d(1:nbpb) , jpi, jpj ) 351 CALL tab_1d_2d( nbpb, rdm_snw , npb, rdm_snw_1d(1:nbpb) , jpi, jpj ) 352 CALL tab_1d_2d( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb) , jpi, jpj ) 353 CALL tab_1d_2d( nbpb, rdvosif , npb, dvsbq_1d (1:nbpb) , jpi, jpj ) 354 CALL tab_1d_2d( nbpb, rdvobif , npb, dvbbq_1d (1:nbpb) , jpi, jpj ) 355 CALL tab_1d_2d( nbpb, fdvolif , npb, dvlbq_1d (1:nbpb) , jpi, jpj ) 356 CALL tab_1d_2d( nbpb, rdvonif , npb, dvnbq_1d (1:nbpb) , jpi, jpj ) 357 CALL tab_1d_2d( nbpb, sfx_thd , npb, sfx_thd_1d(1:nbpb) , jpi, jpj ) 362 358 ! 363 359 IF( num_sal == 2 ) THEN 364 CALL tab_1d_2d( nbpb, fsbri, npb, fsbri_1d(1:nbpb), jpi, jpj )365 CALL tab_1d_2d( nbpb, fhbri , npb, fhbri_1d(1:nbpb), jpi, jpj )360 CALL tab_1d_2d( nbpb, sfx_bri , npb, sfx_bri_1d(1:nbpb) , jpi, jpj ) 361 CALL tab_1d_2d( nbpb, fhbri , npb, fhbri_1d (1:nbpb) , jpi, jpj ) 366 362 ENDIF 367 363 ! 368 !+++++ 369 !temporary stuff for a dummy version 364 !+++++ temporary stuff for a dummy version 370 365 CALL tab_1d_2d( nbpb, dh_i_surf2D, npb, dh_i_surf(1:nbpb) , jpi, jpj ) 371 366 CALL tab_1d_2d( nbpb, dh_i_bott2D, npb, dh_i_bott(1:nbpb) , jpi, jpj ) … … 389 384 ! 5.1) Ice heat content 390 385 !------------------------ 391 ! Enthalpies are global variables we have to readjust the units 386 ! Enthalpies are global variables we have to readjust the units (heat content in 10^9 Joules) 392 387 zcoef = 1._wp / ( unit_fac * REAL( nlay_i ) ) 393 388 DO jl = 1, jpl 394 389 DO jk = 1, nlay_i 395 ! Multiply by volume, divide by nlayers so that heat content in 10^9 Joules396 390 e_i(:,:,jk,jl) = e_i(:,:,jk,jl) * area(:,:) * a_i(:,:,jl) * ht_i(:,:,jl) * zcoef 397 391 END DO … … 401 395 ! 5.2) Snow heat content 402 396 !------------------------ 403 ! Enthalpies are global variables we have to readjust the units 397 ! Enthalpies are global variables we have to readjust the units (heat content in 10^9 Joules) 404 398 zcoef = 1._wp / ( unit_fac * REAL( nlay_s ) ) 405 399 DO jl = 1, jpl 406 400 DO jk = 1, nlay_s 407 ! Multiply by volume, so that heat content in 10^9 Joules408 401 e_s(:,:,jk,jl) = e_s(:,:,jk,jl) * area(:,:) * a_i(:,:,jl) * ht_s(:,:,jl) * zcoef 409 402 END DO … … 419 412 !-------------------------------------------- 420 413 d_v_i_thd(:,:,:) = v_i (:,:,:) - old_v_i(:,:,:) ! ice volumes 421 dv_dt_thd(:,:,:) = d_v_i_thd(:,:,:) / rdt_ice * 86400.0414 dv_dt_thd(:,:,:) = d_v_i_thd(:,:,:) * r1_rdtice * rday 422 415 423 416 IF( con_i ) fbif(:,:) = fbif(:,:) + zqlbsbq(:,:) … … 488 481 ! 489 482 IF(lwp) WRITE(numout,*) ' lim_thd_glohec ' 490 IF(lwp) WRITE(numout,*) ' qt_i_in : ', eti(jiindex_1d,jl) / rdt_ice491 IF(lwp) WRITE(numout,*) ' qt_s_in : ', ets(jiindex_1d,jl) / rdt_ice492 IF(lwp) WRITE(numout,*) ' qt_in : ', ( eti(jiindex_1d,jl) + ets(jiindex_1d,jl) ) / rdt_ice483 IF(lwp) WRITE(numout,*) ' qt_i_in : ', eti(jiindex_1d,jl) * r1_rdtice 484 IF(lwp) WRITE(numout,*) ' qt_s_in : ', ets(jiindex_1d,jl) * r1_rdtice 485 IF(lwp) WRITE(numout,*) ' qt_in : ', ( eti(jiindex_1d,jl) + ets(jiindex_1d,jl) ) * r1_rdtice 493 486 ! 494 487 END SUBROUTINE lim_thd_glohec … … 538 531 !-------------------- 539 532 DO ji = kideb, kiut 540 cons_error(ji,jl) = ABS( dq_i(ji,jl) / rdt_ice + sum_fluxq(ji,jl) )533 cons_error(ji,jl) = ABS( dq_i(ji,jl) * r1_rdtice + sum_fluxq(ji,jl) ) 541 534 END DO 542 535 … … 597 590 WRITE(numout,*) ' cons_error : ', cons_error(ji,jl) 598 591 WRITE(numout,*) ' surf_error : ', surf_error(ji,jl) 599 WRITE(numout,*) ' dq_i : ', - dq_i(ji,jl) / rdt_ice592 WRITE(numout,*) ' dq_i : ', - dq_i(ji,jl) * r1_rdtice 600 593 WRITE(numout,*) ' Fdt : ', sum_fluxq(ji,jl) 601 594 WRITE(numout,*) … … 631 624 WRITE(numout,*) 632 625 WRITE(numout,*) ' Layer by layer ... ' 633 WRITE(numout,*) ' dq_snow : ', ( qt_s_fin(ji,jl) - qt_s_in(ji,jl) ) / rdt_ice626 WRITE(numout,*) ' dq_snow : ', ( qt_s_fin(ji,jl) - qt_s_in(ji,jl) ) * r1_rdtice 634 627 WRITE(numout,*) ' dfc_snow : ', fc_s(ji,1) - fc_s(ji,0) 635 628 DO jk = 1, nlay_i 636 629 WRITE(numout,*) ' layer : ', jk 637 WRITE(numout,*) ' dq_ice : ', dq_i_layer(ji,jk) / rdt_ice630 WRITE(numout,*) ' dq_ice : ', dq_i_layer(ji,jk) * r1_rdtice 638 631 WRITE(numout,*) ' radab : ', radab(ji,jk) 639 632 WRITE(numout,*) ' dfc_i : ', fc_i(ji,jk) - fc_i(ji,jk-1) … … 681 674 fatm (ji,jl) = qnsr_ice_1d(ji) + qsr_ice_1d(ji) ! total heat flux 682 675 sum_fluxq (ji,jl) = fatm(ji,jl) + fbif_1d(ji) - ftotal_fin(ji) - fstroc(zji,zjj,jl) 683 cons_error(ji,jl) = ABS( dq_i(ji,jl) / rdt_ice + sum_fluxq(ji,jl) )676 cons_error(ji,jl) = ABS( dq_i(ji,jl) * r1_rdtice + sum_fluxq(ji,jl) ) 684 677 END DO 685 678 … … 688 681 !-------------------- 689 682 DO ji = kideb, kiut 690 cons_error(ji,jl) = ABS( dq_i(ji,jl) / rdt_ice + sum_fluxq(ji,jl) )683 cons_error(ji,jl) = ABS( dq_i(ji,jl) * r1_rdtice + sum_fluxq(ji,jl) ) 691 684 END DO 692 685 … … 722 715 WRITE(numout,*) ' * ' 723 716 WRITE(numout,*) ' Ftotal : ', sum_fluxq(ji,jl) 724 WRITE(numout,*) ' dq_t : ', - dq_i(ji,jl) / rdt_ice725 WRITE(numout,*) ' dq_i : ', - ( qt_i_fin(ji,jl) - qt_i_in(ji,jl) ) / rdt_ice726 WRITE(numout,*) ' dq_s : ', - ( qt_s_fin(ji,jl) - qt_s_in(ji,jl) ) / rdt_ice717 WRITE(numout,*) ' dq_t : ', - dq_i(ji,jl) * r1_rdtice 718 WRITE(numout,*) ' dq_i : ', - ( qt_i_fin(ji,jl) - qt_i_in(ji,jl) ) * r1_rdtice 719 WRITE(numout,*) ' dq_s : ', - ( qt_s_fin(ji,jl) - qt_s_in(ji,jl) ) * r1_rdtice 727 720 WRITE(numout,*) ' cons_error : ', cons_error(ji,jl) 728 721 WRITE(numout,*) ' * ' … … 734 727 WRITE(numout,*) ' * ' 735 728 WRITE(numout,*) ' Heat contents --- : ' 736 WRITE(numout,*) ' qt_s_in : ', qt_s_in(ji,jl) / rdt_ice737 WRITE(numout,*) ' qt_i_in : ', qt_i_in(ji,jl) / rdt_ice738 WRITE(numout,*) ' qt_in : ', ( qt_i_in(ji,jl) + qt_s_in(ji,jl) ) / rdt_ice739 WRITE(numout,*) ' qt_s_fin : ', qt_s_fin(ji,jl) / rdt_ice740 WRITE(numout,*) ' qt_i_fin : ', qt_i_fin(ji,jl) / rdt_ice741 WRITE(numout,*) ' qt_fin : ', ( qt_i_fin(ji,jl) + qt_s_fin(ji,jl) ) / rdt_ice729 WRITE(numout,*) ' qt_s_in : ', qt_s_in(ji,jl) * r1_rdtice 730 WRITE(numout,*) ' qt_i_in : ', qt_i_in(ji,jl) * r1_rdtice 731 WRITE(numout,*) ' qt_in : ', ( qt_i_in(ji,jl) + qt_s_in(ji,jl) ) * r1_rdtice 732 WRITE(numout,*) ' qt_s_fin : ', qt_s_fin(ji,jl) * r1_rdtice 733 WRITE(numout,*) ' qt_i_fin : ', qt_i_fin(ji,jl) * r1_rdtice 734 WRITE(numout,*) ' qt_fin : ', ( qt_i_fin(ji,jl) + qt_s_fin(ji,jl) ) * r1_rdtice 742 735 WRITE(numout,*) ' * ' 743 736 WRITE(numout,*) ' Ice variables --- : ' -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r3294 r3625 7 7 !! ! 2005-06 (M. Vancoppenolle) 3D version 8 8 !! 3.2 ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in rdmsnif & rdmicif 9 !! 4.0 ! 2011-02 (G. Madec) dynamical allocation 9 !! 3.4 ! 2011-02 (G. Madec) dynamical allocation 10 !! 3.5 ! 2012-10 (G. Madec & co) salt flux + bug fixes 10 11 !!---------------------------------------------------------------------- 11 12 #if defined key_lim3 … … 13 14 !! 'key_lim3' LIM3 sea-ice model 14 15 !!---------------------------------------------------------------------- 15 !! lim_thd_dh : vertical accr./abl. and lateral ablation of sea ice 16 !!---------------------------------------------------------------------- 17 USE par_oce ! ocean parameters 18 USE phycst ! physical constants (OCE directory) 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE ice ! LIM variables 21 USE par_ice ! LIM parameters 22 USE thd_ice ! LIM thermodynamics 23 USE in_out_manager ! I/O manager 24 USE lib_mpp ! MPP library 25 USE wrk_nemo ! work arrays 16 !! lim_thd_dh : vertical accr./abl. and lateral ablation of sea ice 17 !!---------------------------------------------------------------------- 18 USE par_oce ! ocean parameters 19 USE phycst ! physical constants (OCE directory) 20 USE sbc_oce ! Surface boundary condition: ocean fields 21 USE ice ! LIM variables 22 USE par_ice ! LIM parameters 23 USE thd_ice ! LIM thermodynamics 24 USE in_out_manager ! I/O manager 25 USE lib_mpp ! MPP library 26 USE wrk_nemo ! work arrays 27 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 26 28 27 29 IMPLICIT NONE … … 37 39 38 40 !!---------------------------------------------------------------------- 39 !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2010)41 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 40 42 !! $Id$ 41 43 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 71 73 !! 72 74 INTEGER :: ji , jk ! dummy loop indices 73 INTEGER :: zji, zjj ! 2D corresponding indices to ji75 INTEGER :: ii, ij ! 2D corresponding indices to ji 74 76 INTEGER :: isnow ! switch for presence (1) or absence (0) of snow 75 77 INTEGER :: isnowic ! snow ice formation not … … 102 104 REAL(wp), POINTER, DIMENSION(:) :: zfmass_i ! 103 105 104 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_mel 105 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_pre 106 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_sub 107 REAL(wp), POINTER, DIMENSION(:) :: z fsalt_melt ! salt flux due to ice melt106 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_mel ! snow melt 107 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_pre ! snow precipitation 108 REAL(wp), POINTER, DIMENSION(:) :: zdh_s_sub ! snow sublimation 109 REAL(wp), POINTER, DIMENSION(:) :: zsfx_melt ! salt flux due to ice melt 108 110 109 111 REAL(wp), POINTER, DIMENSION(:,:) :: zdeltah … … 126 128 127 129 CALL wrk_alloc( jpij, zh_i, zh_s, ztfs, zhsold, zqprec, zqfont_su, zqfont_bo, z_f_surf, zhgnew, zfmass_i ) 128 CALL wrk_alloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, z fsalt_melt, zfdt_init, zfdt_final, zqt_i, zqt_s, zqt_dummy )130 CALL wrk_alloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zsfx_melt, zfdt_init, zfdt_final, zqt_i, zqt_s, zqt_dummy ) 129 131 CALL wrk_alloc( jpij, zinnermelt, zfbase, zdq_i ) 130 132 CALL wrk_alloc( jpij, jkmax, zdeltah, zqt_i_lay ) 131 133 132 z fsalt_melt(:)= 0._wp133 ftotal_fin(:) 134 zfdt_init (:)= 0._wp135 zfdt_final(:) 134 zsfx_melt (:) = 0._wp 135 ftotal_fin(:) = 0._wp 136 zfdt_init (:) = 0._wp 137 zfdt_final(:) = 0._wp 136 138 137 139 DO ji = kideb, kiut … … 145 147 ! 146 148 DO ji = kideb, kiut 147 isnow = INT( 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_s_b(ji) ) ))148 ztfs (ji)= isnow * rtt + ( 1.0 - isnow ) * rtt149 z_f_surf (ji)= qnsr_ice_1d(ji) + ( 1.0 - i0(ji) ) * qsr_ice_1d(ji) - fc_su(ji)150 z_f_surf (ji) = MAX( zzero , z_f_surf(ji) ) * MAX( zzero , SIGN( zone , t_su_b(ji) - ztfs(ji) ))149 isnow = INT( 1.0 - MAX( 0.0 , SIGN( 1.0 , - ht_s_b(ji) ) ) ) 150 ztfs (ji) = isnow * rtt + ( 1.0 - isnow ) * rtt 151 z_f_surf (ji) = qnsr_ice_1d(ji) + ( 1.0 - i0(ji) ) * qsr_ice_1d(ji) - fc_su(ji) 152 z_f_surf (ji) = MAX( zzero , z_f_surf(ji) ) * MAX( zzero , SIGN( zone , t_su_b(ji) - ztfs(ji) ) ) 151 153 zfdt_init(ji) = ( z_f_surf(ji) + MAX( fbif_1d(ji) + qlbbq_1d(ji) + fc_bo_i(ji),0.0 ) ) * rdt_ice 152 154 END DO ! ji … … 240 242 zhsnew = ht_s_b(ji) + dh_s_tot(ji) 241 243 ! If snow is still present zhn = 1, else zhn = 0 242 zhn = 1.0 - MAX( zzero , SIGN( zone , - zhsnew ))244 zhn = 1.0 - MAX( zzero , SIGN( zone , - zhsnew ) ) 243 245 ht_s_b(ji) = MAX( zzero , zhsnew ) 244 246 ! Volume and mass variations of snow 245 dvsbq_1d (ji) = a_i_b(ji) * ( ht_s_b(ji) - zhsold(ji) - zdh_s_ mel(ji) )247 dvsbq_1d (ji) = a_i_b(ji) * ( ht_s_b(ji) - zhsold(ji) - zdh_s_pre(ji) ) 246 248 dvsbq_1d (ji) = MIN( zzero, dvsbq_1d(ji) ) 247 rdm snif_1d(ji) = rdmsnif_1d(ji) + rhosn * dvsbq_1d(ji)249 rdm_snw_1d(ji) = rdm_snw_1d(ji) + rhosn * dvsbq_1d(ji) 248 250 END DO ! ji 249 251 … … 253 255 DO ji = kideb, kiut 254 256 dh_i_surf(ji) = 0._wp 255 z_f_surf (ji) = zqfont_su(ji) / rdt_ice! heat conservation test257 z_f_surf (ji) = zqfont_su(ji) * r1_rdtice ! heat conservation test 256 258 zdq_i (ji) = 0._wp 257 259 END DO ! ji … … 262 264 zdeltah (ji,jk) = - zqfont_su(ji) / q_i_b(ji,jk) 263 265 ! ! recompute heat available 264 zqfont_su(ji )= MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk)266 zqfont_su(ji ) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk) 265 267 ! ! melt of layer jk cannot be higher than its thickness 266 268 zdeltah (ji,jk) = MAX( zdeltah(ji,jk) , - zh_i(ji) ) 267 269 ! ! update surface melt 268 dh_i_surf(ji )= dh_i_surf(ji) + zdeltah(ji,jk)270 dh_i_surf(ji ) = dh_i_surf(ji) + zdeltah(ji,jk) 269 271 ! ! for energy conservation 270 zdq_i (ji ) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) / rdt_ice272 zdq_i (ji ) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) * r1_rdtice 271 273 ! 272 ! contribution to ice-ocean salt flux 273 zji = MOD( npb(ji) - 1 , jpi ) + 1 274 zjj = ( npb(ji) - 1 ) / jpi + 1 275 zfsalt_melt(ji) = zfsalt_melt(ji) + ( sss_m(zji,zjj) - sm_i_b(ji) ) * a_i_b(ji) & 276 & * MIN( zdeltah(ji,jk) , 0.e0 ) * rhoic / rdt_ice 274 ! ! contribution to ice-ocean salt flux 275 zsfx_melt(ji) = zsfx_melt(ji) - sm_i_b(ji) * a_i_b(ji) * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic * r1_rdtice 277 276 END DO 278 277 END DO … … 290 289 IF( z_f_surf(ji) + zdq_i(ji) .GE. 1.0e-3 ) THEN! 291 290 WRITE(numout,*) ' ALERTE heat loss for surface melt ' 292 WRITE(numout,*) ' zji, zjj, jl :', zji, zjj, jl291 WRITE(numout,*) ' ii, ij, jl :', ii, ij, jl 293 292 WRITE(numout,*) ' ht_i_b : ', ht_i_b(ji) 294 293 WRITE(numout,*) ' z_f_surf : ', z_f_surf(ji) … … 299 298 WRITE(numout,*) ' qlbbq_1d : ', qlbbq_1d(ji) 300 299 WRITE(numout,*) ' s_i_new : ', s_i_new(ji) 301 WRITE(numout,*) ' sss_m : ', sss_m( zji,zjj)300 WRITE(numout,*) ' sss_m : ', sss_m(ii,ij) 302 301 ENDIF 303 302 END DO … … 338 337 DO ji = kideb, kiut 339 338 ! In case of disparition of the snow, we have to update the snow temperatures 340 zhisn = MAX( zzero , SIGN( zone, - ht_s_b(ji) ))339 zhisn = MAX( zzero , SIGN( zone, - ht_s_b(ji) ) ) 341 340 t_s_b(ji,jk) = ( 1.0 - zhisn ) * t_s_b(ji,jk) + zhisn * rtt 342 341 q_s_b(ji,jk) = ( 1.0 - zhisn ) * q_s_b(ji,jk) … … 358 357 ! 4.1 Basal growth - (a) salinity not varying in time 359 358 !----------------------------------------------------- 360 IF( num_sal /= 2 .AND. num_sal /= 4 ) THEN361 DO ji = kideb, kiut 362 IF( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) < 0. 0) THEN359 IF( num_sal /= 2 ) THEN ! ice salinity constant in time 360 DO ji = kideb, kiut 361 IF( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) < 0._wp ) THEN 363 362 s_i_new(ji) = sm_i_b(ji) 364 363 ! Melting point in K … … 371 370 & - rcp * ( ztmelts - rtt ) ) 372 371 ! Basal growth rate = - F*dt / q 373 dh_i_bott(ji) = - rdt_ice *( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1)372 dh_i_bott(ji) = - rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1) 374 373 ENDIF 375 374 END DO … … 379 378 ! 4.1 Basal growth - (b) salinity varying in time 380 379 !------------------------------------------------- 381 IF( num_sal == 2 .OR. num_sal == 4 ) THEN 382 ! the growth rate (dh_i_bott) is function of the new ice 383 ! heat content (q_i_b(nlay_i+1)). q_i_b depends on the new ice 384 ! salinity (snewice). snewice depends on dh_i_bott 385 ! it converges quickly, so, no problem 380 IF( num_sal == 2 ) THEN 381 ! the growth rate (dh_i_bott) is function of the new ice heat content (q_i_b(nlay_i+1)). 382 ! q_i_b depends on the new ice salinity (snewice). 383 ! snewice depends on dh_i_bott ; it converges quickly, so, no problem 386 384 ! See Vancoppenolle et al., OM08 for more info on this 387 385 … … 394 392 DO ji = kideb, kiut 395 393 IF( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) < 0.e0 ) THEN 396 zji = MOD( npb(ji) - 1, jpi ) + 1397 zjj = ( npb(ji) - 1 ) / jpi + 1394 ii = MOD( npb(ji) - 1, jpi ) + 1 395 ij = ( npb(ji) - 1 ) / jpi + 1 398 396 ! Melting point in K 399 397 ztmelts = - tmut * s_i_new(ji) + rtt … … 408 406 ! zswi12 (1) if dh_i_bott/rdt .LT. 3.6e-7 and .GT. 2.0e-8 409 407 ! zswi1 (1) if dh_i_bott/rdt .LT. 2.0e-8 410 zgrr = MIN( 1.0e-3, MAX ( dh_i_bott(ji) / rdt_ice , epsi13 ) )408 zgrr = MIN( 1.0e-3, MAX ( dh_i_bott(ji) * r1_rdtice , epsi13 ) ) 411 409 zswi2 = MAX( zzero , SIGN( zone , zgrr - 3.6e-7 ) ) 412 410 zswi12 = MAX( zzero , SIGN( zone , zgrr - 2.0e-8 ) ) * ( 1.0 - zswi2 ) … … 414 412 zfracs = zswi1 * 0.12 + zswi12 * ( 0.8925 + 0.0568 * LOG( 100.0 * zgrr ) ) & 415 413 & + zswi2 * 0.26 / ( 0.26 + 0.74 * EXP ( - 724300.0 * zgrr ) ) 416 zds = zfracs * sss_m( zji,zjj) - s_i_new(ji)417 s_i_new(ji) = zfracs * sss_m( zji,zjj)414 zds = zfracs * sss_m(ii,ij) - s_i_new(ji) 415 s_i_new(ji) = zfracs * sss_m(ii,ij) 418 416 ENDIF ! fc_bo_i 419 417 END DO ! ji … … 432 430 & - rcp * ( ztmelts - rtt ) ) 433 431 ! Basal growth rate = - F*dt / q 434 dh_i_bott(ji) = - rdt_ice *( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1)432 dh_i_bott(ji) = - rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1) 435 433 ! Salinity update 436 434 ! entrapment during bottom growth … … 453 451 s_i_new(ji) = s_i_b(ji,nlay_i) 454 452 zqfont_bo(ji) = rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) 455 zfbase(ji) = zqfont_bo(ji) / rdt_ice ! heat conservation test453 zfbase(ji) = zqfont_bo(ji) * r1_rdtice ! heat conservation test 456 454 zdq_i(ji) = 0._wp 457 455 dh_i_bott(ji) = 0._wp … … 461 459 DO jk = nlay_i, 1, -1 462 460 DO ji = kideb, kiut 463 IF ( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) .GE. 0.0 ) THEN 464 ztmelts = - tmut * s_i_b(ji,jk) + rtt 465 IF( t_i_b(ji,jk) >= ztmelts ) THEN 466 zdeltah(ji,jk) = - zh_i(ji) 467 dh_i_bott(ji) = dh_i_bott(ji) + zdeltah(ji,jk) 468 zinnermelt(ji) = 1._wp 469 ELSE ! normal ablation 470 zdeltah(ji,jk) = - zqfont_bo(ji) / q_i_b(ji,jk) 471 zqfont_bo(ji) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk) 472 zdeltah(ji,jk) = MAX(zdeltah(ji,jk), - zh_i(ji) ) 473 dh_i_bott(ji) = dh_i_bott(ji) + zdeltah(ji,jk) 474 zdq_i(ji) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) / rdt_ice 475 ! contribution to salt flux 476 zji = MOD( npb(ji) - 1, jpi ) + 1 477 zjj = ( npb(ji) - 1 ) / jpi + 1 478 zfsalt_melt(ji) = zfsalt_melt(ji) + ( sss_m(zji,zjj) - sm_i_b(ji) ) * a_i_b(ji) & 479 & * MIN( zdeltah(ji,jk) , 0.0 ) * rhoic / rdt_ice 461 IF( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) >= 0._wp ) THEN 462 ztmelts = - tmut * s_i_b(ji,jk) + rtt 463 IF( t_i_b(ji,jk) >= ztmelts ) THEN !!gm : a comment is needed 464 zdeltah (ji,jk) = - zh_i(ji) 465 dh_i_bott (ji ) = dh_i_bott(ji) + zdeltah(ji,jk) 466 zinnermelt(ji ) = 1._wp 467 ELSE ! normal ablation 468 zdeltah (ji,jk) = - zqfont_bo(ji) / q_i_b(ji,jk) 469 zqfont_bo(ji ) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk) 470 zdeltah (ji,jk) = MAX(zdeltah(ji,jk), - zh_i(ji) ) 471 dh_i_bott(ji ) = dh_i_bott(ji) + zdeltah(ji,jk) 472 zdq_i (ji ) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) * r1_rdtice 480 473 ENDIF 474 ! contribution to salt flux 475 zsfx_melt(ji) = zsfx_melt(ji) - sm_i_b(ji) * a_i_b(ji) * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic * r1_rdtice 481 476 ENDIF 482 477 END DO ! ji … … 493 488 ENDIF 494 489 IF ( zfbase(ji) + zdq_i(ji) .GE. 1.0e-3 ) THEN 495 WRITE(numout,*) ' ALERTE heat loss for basal melt : zji, zjj, jl :', zji, zjj, jl490 WRITE(numout,*) ' ALERTE heat loss for basal melt : ii, ij, jl :', ii, ij, jl 496 491 WRITE(numout,*) ' ht_i_b : ', ht_i_b(ji) 497 492 WRITE(numout,*) ' zfbase : ', zfbase(ji) … … 502 497 WRITE(numout,*) ' qlbbq_1d : ', qlbbq_1d(ji) 503 498 WRITE(numout,*) ' s_i_new : ', s_i_new(ji) 504 WRITE(numout,*) ' sss_m : ', sss_m( zji,zjj)499 WRITE(numout,*) ' sss_m : ', sss_m(ii,ij) 505 500 WRITE(numout,*) ' dh_i_bott : ', dh_i_bott(ji) 506 501 WRITE(numout,*) ' innermelt : ', INT( zinnermelt(ji) ) … … 531 526 ! ! excessive energy is sent to lateral ablation 532 527 fsup (ji) = rhoic * lfus * at_i_b(ji) / MAX( 1.0 - at_i_b(ji) , epsi13 ) & 533 & * ( zdhbf - dh_i_bott(ji) ) / rdt_ice528 & * ( zdhbf - dh_i_bott(ji) ) * r1_rdtice 534 529 dh_i_bott(ji) = zdhbf 535 530 ! !since ice volume is only used for outputs, we keep it global for all categories … … 538 533 zhgnew (ji) = ht_i_b(ji) + dh_i_surf(ji) + dh_i_bott(ji) 539 534 ! ! diagnostic ( bottom ice growth ) 540 zji = MOD( npb(ji) - 1, jpi ) + 1541 zjj = ( npb(ji) - 1 ) / jpi + 1542 diag_bot_gr( zji,zjj) = diag_bot_gr(zji,zjj) + MAX(dh_i_bott(ji),0.0)*a_i_b(ji) / rdt_ice543 diag_sur_me( zji,zjj) = diag_sur_me(zji,zjj) + MIN(dh_i_surf(ji),0.0)*a_i_b(ji) / rdt_ice544 diag_bot_me( zji,zjj) = diag_bot_me(zji,zjj) + MIN(dh_i_bott(ji),0.0)*a_i_b(ji) / rdt_ice535 ii = MOD( npb(ji) - 1, jpi ) + 1 536 ij = ( npb(ji) - 1 ) / jpi + 1 537 diag_bot_gr(ii,ij) = diag_bot_gr(ii,ij) + MAX(dh_i_bott(ji),0.0)*a_i_b(ji) * r1_rdtice 538 diag_sur_me(ii,ij) = diag_sur_me(ii,ij) + MIN(dh_i_surf(ji),0.0)*a_i_b(ji) * r1_rdtice 539 diag_bot_me(ii,ij) = diag_bot_me(ii,ij) + MIN(dh_i_bott(ji),0.0)*a_i_b(ji) * r1_rdtice 545 540 END DO 546 541 … … 548 543 ! 5.2 More than available ice melts 549 544 !----------------------------------- 550 ! then heat applied minus heat content at previous time step 551 ! should equal heat remaining 545 ! then heat applied minus heat content at previous time step should equal heat remaining 552 546 ! 553 547 DO ji = kideb, kiut 554 548 ! Adapt the remaining energy if too much ice melts 555 549 !-------------------------------------------------- 556 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) !1 if ice550 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) ! =1 if ice 557 551 ! 0 if no more ice 558 552 zhgnew (ji) = zihgnew * zhgnew(ji) ! ice thickness is put to 0 … … 562 556 ! If snow remains, energy is used to melt snow 563 557 zhni = ht_s_b(ji) ! snow depth at previous time step 564 zihg = MAX( zzero , SIGN ( zone , - ht_s_b(ji) ) ) !0 if snow558 zihg = MAX( zzero , SIGN ( zone , - ht_s_b(ji) ) ) ! =0 if snow 565 559 566 560 ! energy of melting of remaining snow 567 561 zqt_s(ji) = ( 1. - zihg ) * zqt_s(ji) / MAX( zhni, epsi13 ) 568 562 zdhnm = - ( 1. - zihg ) * ( 1. - zihgnew ) * zfdt_final(ji) / MAX( zqt_s(ji) , epsi13 ) 569 zhnfi 563 zhnfi = zhni + zdhnm 570 564 zfdt_final(ji) = MAX( zfdt_final(ji) + zqt_s(ji) * zdhnm , 0.0 ) 571 565 ht_s_b(ji) = MAX( zzero , zhnfi ) … … 581 575 ! 582 576 ! ! mass variation cumulated over category 583 rdm snif_1d(ji) = rdmsnif_1d(ji) + zzfmass_s ! snow584 rdm icif_1d(ji) = rdmicif_1d(ji) + zzfmass_i ! ice577 rdm_snw_1d(ji) = rdm_snw_1d(ji) + zzfmass_s ! snow 578 rdm_ice_1d(ji) = rdm_ice_1d(ji) + zzfmass_i ! ice 585 579 586 580 ! Remaining heat to the ocean 587 581 !--------------------------------- 588 focea(ji) = - zfdt_final(ji) / rdt_ice ! focea is in W.m-2 * dt589 590 END DO 591 592 ftotal_fin (:) = zfdt_final(:) / rdt_ice582 focea(ji) = - zfdt_final(ji) * r1_rdtice ! focea is in W.m-2 * dt 583 584 END DO 585 586 ftotal_fin (:) = zfdt_final(:) * r1_rdtice 593 587 594 588 !--------------------------- … … 596 590 !--------------------------- 597 591 DO ji = kideb, kiut 598 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) ! 1 if ice599 592 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) ! =1 if ice 593 ! 600 594 ! Salt flux 601 zji = MOD( npb(ji) - 1, jpi ) + 1 602 zjj = ( npb(ji) - 1 ) / jpi + 1 603 ! new lines 604 IF( num_sal == 4 ) THEN 605 fseqv_1d(ji) = fseqv_1d(ji) + zihgnew * zfsalt_melt(ji) & 606 & + (1.0 - zihgnew) * zfmass_i(ji) * ( sss_m(zji,zjj) - bulk_sal ) / rdt_ice 607 ELSE 608 fseqv_1d(ji) = fseqv_1d(ji) + zihgnew * zfsalt_melt(ji) & 609 & + (1.0 - zihgnew) * zfmass_i(ji) * ( sss_m(zji,zjj) - sm_i_b(ji) ) / rdt_ice 610 ENDIF 595 sfx_thd_1d(ji) = sfx_thd_1d(ji) + zihgnew * zsfx_melt(ji) & 596 & - (1.0 - zihgnew) * zfmass_i (ji) * sm_i_b(ji) * r1_rdtice 597 ! 611 598 ! Heat flux 612 599 ! excessive bottom ablation energy (fsup) - 0 except if jpl = 1 613 ! excessive total ablation energy (focea) sent to the ocean600 ! excessive total ablation energy (focea) sent to the ocean 614 601 qfvbq_1d(ji) = qfvbq_1d(ji) + fsup(ji) + ( 1.0 - zihgnew ) * focea(ji) * a_i_b(ji) * rdt_ice 615 602 616 zihic = 1.0 - MAX( zzero , SIGN( zone , -ht_i_b(ji) ) ) 617 ! equals 0 if ht_i = 0, 1 if ht_i gt 0 603 zihic = 1.0 - MAX( zzero , SIGN( zone , -ht_i_b(ji) ) ) ! equals 0 if ht_i = 0, 1 if ht_i gt 0 618 604 fscbq_1d(ji) = a_i_b(ji) * fstbif_1d(ji) 619 qldif_1d(ji) = qldif_1d(ji) + fsup(ji) + ( 1.0 - zihgnew ) * focea (ji)* a_i_b(ji) * rdt_ice &605 qldif_1d(ji) = qldif_1d(ji) + fsup(ji) + ( 1.0 - zihgnew ) * focea (ji) * a_i_b(ji) * rdt_ice & 620 606 & + ( 1.0 - zihic ) * fscbq_1d(ji) * rdt_ice 621 607 END DO ! ji … … 656 642 dmgwi_1d (ji) = dmgwi_1d(ji) + a_i_b(ji) * ( ht_s_b(ji) - zhnnew ) * rhosn 657 643 658 rdmicif_1d(ji) = rdmicif_1d(ji) + a_i_b(ji) * ( zhgnew(ji) - ht_i_b(ji) ) * rhoic 659 rdmsnif_1d(ji) = rdmsnif_1d(ji) + a_i_b(ji) * ( zhnnew - ht_s_b(ji) ) * rhosn 644 ! All snow is thrown in the ocean, and seawater is taken to replace the volume 645 rdm_ice_1d(ji) = rdm_ice_1d(ji) + a_i_b(ji) * ( zhgnew(ji) - ht_i_b(ji) ) * rhoic * ( 1. - rhosn / rhoic ) 646 rdm_snw_1d(ji) = rdm_snw_1d(ji) + a_i_b(ji) * ( zhnnew - ht_s_b(ji) ) * rhosn 660 647 661 648 ! Equivalent salt flux (1) Snow-ice formation component 662 649 ! ----------------------------------------------------- 663 zji = MOD( npb(ji) - 1, jpi ) + 1664 zjj = ( npb(ji) - 1 ) / jpi + 1665 666 IF( num_sal /= 2 ) THEN ; zsm_snowice = sm_i_b(ji)667 ELSE ; zsm_snowice = ( rhoic - rhosn ) / rhoic * sss_m(zji,zjj)650 ii = MOD( npb(ji) - 1, jpi ) + 1 651 ij = ( npb(ji) - 1 ) / jpi + 1 652 653 IF( num_sal == 2 ) THEN ; zsm_snowice = sss_m(ii,ij) * ( rhoic - rhosn ) / rhoic 654 ELSE ; zsm_snowice = sm_i_b(ji) 668 655 ENDIF 669 IF( num_sal == 4 ) THEN 670 fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - bulk_sal ) * a_i_b(ji) & 671 & * ( zhgnew(ji) - ht_i_b(ji) ) * rhoic / rdt_ice 672 ELSE 673 fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - zsm_snowice ) * a_i_b(ji) & 674 & * ( zhgnew(ji) - ht_i_b(ji) ) * rhoic / rdt_ice 675 ENDIF 656 sfx_thd_1d(ji) = sfx_thd_1d(ji) - zsm_snowice * a_i_b(ji) * dh_snowice(ji) * rhoic * r1_rdtice 657 ! 676 658 ! entrapment during snow ice formation 677 659 i_ice_switch = 1.0 - MAX( 0.e0 , SIGN( 1.0 , - ht_i_b(ji) + 1.0e-6 ) ) 678 660 isnowic = 1.0 - MAX( 0.e0 , SIGN( 1.0 , - dh_snowice(ji) ) ) * i_ice_switch 679 IF( num_sal == 2 .OR. num_sal == 4) &680 dsm_i_si_1d(ji) = ( zsm_snowice*dh_snowice(ji)&681 & + sm_i_b(ji) * ht_i_b(ji) / MAX( ht_i_b(ji) + dh_snowice(ji), epsi13) &682 & - sm_i_b(ji)) * isnowic661 IF( num_sal == 2 ) & 662 dsm_i_si_1d(ji) = ( zsm_snowice * dh_snowice(ji) & 663 & + sm_i_b(ji) * ht_i_b(ji) / MAX( ht_i_b(ji) + dh_snowice(ji), epsi13 ) & 664 & - sm_i_b(ji) ) * isnowic 683 665 684 666 ! Actualize new snow and ice thickness. … … 690 672 691 673 ! diagnostic ( snow ice growth ) 692 zji = MOD( npb(ji) - 1, jpi ) + 1693 zjj = ( npb(ji) - 1 ) / jpi + 1694 diag_sni_gr( zji,zjj) = diag_sni_gr(zji,zjj) + dh_snowice(ji)*a_i_b(ji) / rdt_ice674 ii = MOD( npb(ji) - 1, jpi ) + 1 675 ij = ( npb(ji) - 1 ) / jpi + 1 676 diag_sni_gr(ii,ij) = diag_sni_gr(ii,ij) + dh_snowice(ji)*a_i_b(ji) * r1_rdtice 695 677 ! 696 678 END DO !ji 697 679 ! 698 680 CALL wrk_dealloc( jpij, zh_i, zh_s, ztfs, zhsold, zqprec, zqfont_su, zqfont_bo, z_f_surf, zhgnew, zfmass_i ) 699 CALL wrk_dealloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, z fsalt_melt, zfdt_init, zfdt_final, zqt_i, zqt_s, zqt_dummy )681 CALL wrk_dealloc( jpij, zdh_s_mel, zdh_s_pre, zdh_s_sub, zsfx_melt, zfdt_init, zfdt_final, zqt_i, zqt_s, zqt_dummy ) 700 682 CALL wrk_dealloc( jpij, zinnermelt, zfbase, zdq_i ) 701 683 CALL wrk_dealloc( jpij, jkmax, zdeltah, zqt_i_lay ) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limthd_dif.F90
r3610 r3625 15 15 !! 'key_lim3' LIM3 sea-ice model 16 16 !!---------------------------------------------------------------------- 17 USE par_oce ! ocean parameters 18 USE phycst ! physical constants (ocean directory) 19 USE ice ! LIM-3 variables 20 USE par_ice ! LIM-3 parameters 21 USE thd_ice ! LIM-3: thermodynamics 22 USE in_out_manager ! I/O manager 23 USE lib_mpp ! MPP library 24 USE wrk_nemo ! work arrays 17 USE par_oce ! ocean parameters 18 USE phycst ! physical constants (ocean directory) 19 USE ice ! LIM-3 variables 20 USE par_ice ! LIM-3 parameters 21 USE thd_ice ! LIM-3: thermodynamics 22 USE in_out_manager ! I/O manager 23 USE lib_mpp ! MPP library 24 USE wrk_nemo ! work arrays 25 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 25 26 26 27 IMPLICIT NONE … … 33 34 34 35 !!---------------------------------------------------------------------- 35 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)36 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 36 37 !! $Id$ 37 38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 147 148 REAL(wp), DIMENSION(kiut,jkmax+2) :: zdiagbis 148 149 REAL(wp), DIMENSION(kiut,jkmax+2,3) :: ztrid ! tridiagonal system terms 149 !!------------------------------------------------------------------ 150 150 !!------------------------------------------------------------------ 151 151 ! 152 152 !------------------------------------------------------------------------------! … … 156 156 DO ji = kideb , kiut 157 157 ! is there snow or not 158 isnow(ji)= INT( 1._wp - MAX( 0._wp , SIGN(1._wp, - ht_s_b(ji) )) )158 isnow(ji)= INT( 1._wp - MAX( 0._wp , SIGN(1._wp, - ht_s_b(ji) ) ) ) 159 159 ! surface temperature of fusion 160 160 !!gm ??? ztfs(ji) = rtt !!!???? … … 201 201 DO ji = kideb , kiut 202 202 ! switches 203 isnow(ji) = INT( 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_s_b(ji) )) )203 isnow(ji) = INT( 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_s_b(ji) ) ) ) 204 204 ! hs > 0, isnow = 1 205 205 zhsu (ji) = hnzst ! threshold for the computation of i0 … … 262 262 ! just to check energy conservation 263 263 DO ji = kideb, kiut 264 ii = MOD( npb(ji) - 1, jpi ) + 1265 ij =( npb(ji) - 1 ) / jpi + 1264 ii = MOD( npb(ji) - 1 , jpi ) + 1 265 ij = ( npb(ji) - 1 ) / jpi + 1 266 266 fstroc(ii,ij,jl) = zradtr_i(ji,nlay_i) 267 267 END DO … … 273 273 END DO 274 274 END DO 275 276 275 277 276 ! … … 662 661 663 662 ! surface temperature 664 isnow(ji) = INT( 1.0-max(0.0,sign(1.0,-ht_s_b(ji))))663 isnow(ji) = INT( 1.0 - MAX( 0.0 , SIGN( 1.0 , -ht_s_b(ji) ) ) ) 665 664 ztsuoldit(ji) = t_su_b(ji) 666 IF (t_su_b(ji) .LT. ztfs(ji))&665 IF( t_su_b(ji) < ztfs(ji) ) & 667 666 t_su_b(ji) = ( zindtbis(ji,numeqmin(ji)) - ztrid(ji,numeqmin(ji),3)* ( isnow(ji)*t_s_b(ji,1) & 668 667 & + (1.0-isnow(ji))*t_i_b(ji,1) ) ) / zdiagbis(ji,numeqmin(ji)) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limthd_ent.F90
r3294 r3625 16 16 !! 'key_lim3' LIM3 sea-ice model 17 17 !!---------------------------------------------------------------------- 18 !! lim_thd_ent : ice redistribution of enthalpy18 !! lim_thd_ent : ice redistribution of enthalpy 19 19 !!---------------------------------------------------------------------- 20 USE par_oce ! ocean parameters 21 USE dom_oce ! domain variables 22 USE domain ! 23 USE phycst ! physical constants 24 USE ice ! LIM variables 25 USE par_ice ! LIM parameters 26 USE thd_ice ! LIM thermodynamics 27 USE limvar ! LIM variables 28 USE in_out_manager ! I/O manager 29 USE lib_mpp ! MPP library 30 USE wrk_nemo ! work arrays 20 USE par_oce ! ocean parameters 21 USE dom_oce ! domain variables 22 USE domain ! 23 USE phycst ! physical constants 24 USE ice ! LIM variables 25 USE par_ice ! LIM parameters 26 USE thd_ice ! LIM thermodynamics 27 USE limvar ! LIM variables 28 USE in_out_manager ! I/O manager 29 USE lib_mpp ! MPP library 30 USE wrk_nemo ! work arrays 31 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 31 32 32 33 IMPLICIT NONE … … 43 44 44 45 !!---------------------------------------------------------------------- 45 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)46 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 46 47 !! $Id$ 47 48 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 408 409 IF ( con_i ) THEN 409 410 DO ji = kideb, kiut 410 IF ( ABS ( zqts_in(ji) - zqts_fin(ji) ) / rdt_ice .GT.1.0e-6 ) THEN411 IF ( ABS ( zqts_in(ji) - zqts_fin(ji) ) * r1_rdtice > 1.0e-6 ) THEN 411 412 zji = MOD( npb(ji) - 1, jpi ) + 1 412 413 zjj = ( npb(ji) - 1 ) / jpi + 1 413 414 WRITE(numout,*) ' violation of heat conservation : ', & 414 ABS ( zqts_in(ji) - zqts_fin(ji) ) / rdt_ice415 ABS ( zqts_in(ji) - zqts_fin(ji) ) * r1_rdtice 415 416 WRITE(numout,*) ' ji, jj : ', zji, zjj 416 417 WRITE(numout,*) ' ht_s_b : ', ht_s_b(ji) 417 WRITE(numout,*) ' zqts_in : ', zqts_in (ji) / rdt_ice418 WRITE(numout,*) ' zqts_fin : ', zqts_fin(ji) / rdt_ice418 WRITE(numout,*) ' zqts_in : ', zqts_in (ji) * r1_rdtice 419 WRITE(numout,*) ' zqts_fin : ', zqts_fin(ji) * r1_rdtice 419 420 WRITE(numout,*) ' dh_snowice : ', dh_snowice(ji) 420 421 WRITE(numout,*) ' dh_s_tot : ', dh_s_tot(ji) … … 526 527 ! bottom formation temperature 527 528 ztform = t_i_b(ji,nlay_i) 528 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )ztform = t_bo_b(ji)529 IF( num_sal == 2 ) ztform = t_bo_b(ji) 529 530 qm0(ji,nbot0(ji)) = ( 1.0 - icboswi(ji) )*qm0(ji,nbot0(ji)) & ! case of melting ice 530 531 & + icboswi(ji) * rhoic * ( cpic*(ztmelts-ztform) & ! case of forming ice … … 622 623 ! 623 624 DO ji = kideb, kiut 624 IF ( ABS ( zqti_in(ji) - zqti_fin(ji) ) / rdt_ice .GT.1.0e-6 ) THEN625 IF ( ABS ( zqti_in(ji) - zqti_fin(ji) ) * r1_rdtice > 1.0e-6 ) THEN 625 626 zji = MOD( npb(ji) - 1, jpi ) + 1 626 627 zjj = ( npb(ji) - 1 ) / jpi + 1 627 WRITE(numout,*) ' violation of heat conservation : ', ABS ( zqti_in(ji) - zqti_fin(ji) ) / rdt_ice628 WRITE(numout,*) ' violation of heat conservation : ', ABS ( zqti_in(ji) - zqti_fin(ji) ) * r1_rdtice 628 629 WRITE(numout,*) ' ji, jj : ', zji, zjj 629 630 WRITE(numout,*) ' ht_i_b : ', ht_i_b(ji) 630 WRITE(numout,*) ' zqti_in : ', zqti_in (ji) / rdt_ice631 WRITE(numout,*) ' zqti_fin : ', zqti_fin(ji) / rdt_ice631 WRITE(numout,*) ' zqti_in : ', zqti_in (ji) * r1_rdtice 632 WRITE(numout,*) ' zqti_fin : ', zqti_fin(ji) * r1_rdtice 632 633 WRITE(numout,*) ' dh_i_bott: ', dh_i_bott(ji) 633 634 WRITE(numout,*) ' dh_i_surf: ', dh_i_surf(ji) -
branches/2012/dev_NOC_2012_rev3555/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90
r3294 r3625 13 13 !! 'key_lim3' LIM3 sea-ice model 14 14 !!---------------------------------------------------------------------- 15 !! lim_lat_acr : lateral accretion of ice 16 !!---------------------------------------------------------------------- 17 USE par_oce ! ocean parameters 18 USE dom_oce ! domain variables 19 USE phycst ! physical constants 20 USE sbc_oce ! Surface boundary condition: ocean fields 21 USE sbc_ice ! Surface boundary condition: ice fields 22 USE thd_ice ! LIM thermodynamics 23 USE dom_ice ! LIM domain 24 USE par_ice ! LIM parameters 25 USE ice ! LIM variables 26 USE limtab ! LIM 2D <==> 1D 27 USE limcons ! LIM conservation 28 USE in_out_manager ! I/O manager 29 USE lib_mpp ! MPP library 30 USE wrk_nemo ! work arrays 15 !! lim_lat_acr : lateral accretion of ice 16 !!---------------------------------------------------------------------- 17 USE par_oce ! ocean parameters 18 USE dom_oce ! domain variables 19 USE phycst ! physical constants 20 USE sbc_oce ! Surface boundary condition: ocean fields 21 USE sbc_ice ! Surface boundary condition: ice fields 22 USE thd_ice ! LIM thermodynamics 23 USE dom_ice ! LIM domain 24 USE par_ice ! LIM parameters 25 USE ice ! LIM variables 26 USE limtab ! LIM 2D <==> 1D 27 USE limcons ! LIM conservation 28 USE in_out_manager ! I/O manager 29 USE lib_mpp ! MPP library 30 USE wrk_nemo ! work arrays 31 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 31 32 32 33 IMPLICIT NONE … … 45 46 46 47 !!---------------------------------------------------------------------- 47 !! NEMO/LIM3 4.0, UCL - NEMO Consortium (2011)48 !! NEMO/LIM3 3.4 , UCL - NEMO Consortium (2011) 48 49 !! $Id$ 49 50 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 77 78 !! update ht_s_b, ht_i_b and tbif_1d(:,:) 78 79 !!------------------------------------------------------------------------ 79 INTEGER :: ji,jj,jk,jl,jm ! dummy loop indices80 INTEGER :: layer, nbpac ! local integers81 INTEGER :: zji, zjj, iter ! - -82 REAL(wp) :: ztmelts, zdv, zqold, zfrazb, zweight, zalphai, zindb, zde! local scalars80 INTEGER :: ji,jj,jk,jl,jm ! dummy loop indices 81 INTEGER :: layer, nbpac ! local integers 82 INTEGER :: zji, zjj, iter ! - - 83 REAL(wp) :: ztmelts, zdv, zqold, zfrazb, zweight, zalphai, zindb, zde ! local scalars 83 84 REAL(wp) :: zgamafr, zvfrx, zvgx, ztaux, ztwogp, zf , zhicol_new ! - - 84 85 REAL(wp) :: ztenagm, zvfry, zvgy, ztauy, zvrel2, zfp, zsqcd , zhicrit ! - - 86 REAL(wp) :: zcoef ! - - 85 87 LOGICAL :: iterate_frazil ! iterate frazil ice collection thickness 86 88 CHARACTER (len = 15) :: fieldid … … 143 145 ! 1) Conservation check and changes in each ice category 144 146 !------------------------------------------------------------------------------! 145 IF 146 CALL lim_column_sum (jpl, v_i, vt_i_init)147 CALL lim_column_sum (jpl, v_s, vt_s_init)148 CALL lim_column_sum_energy ( jpl, nlay_i, e_i, et_i_init)149 CALL lim_column_sum (jpl,e_s(:,:,1,:) , et_s_init)147 IF( con_i ) THEN 148 CALL lim_column_sum ( jpl, v_i , vt_i_init) 149 CALL lim_column_sum ( jpl, v_s , vt_s_init) 150 CALL lim_column_sum_energy ( jpl, nlay_i , e_i , et_i_init) 151 CALL lim_column_sum ( jpl, e_s(:,:,1,:) , et_s_init) 150 152 ENDIF 151 153 … … 158 160 DO ji = 1, jpi 159 161 !Energy of melting q(S,T) [J.m-3] 160 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / & 161 MAX( area(ji,jj) * v_i(ji,jj,jl) , epsi10 ) * & 162 nlay_i 163 zindb = 1.0-MAX(0.0,SIGN(1.0,-v_i(ji,jj,jl))) !0 if no ice and 1 if yes 164 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl)*unit_fac*zindb 162 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / MAX( area(ji,jj) * v_i(ji,jj,jl) , epsi10 ) * nlay_i 163 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , -v_i(ji,jj,jl) ) ) !0 if no ice and 1 if yes 164 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * unit_fac * zindb 165 165 END DO 166 166 END DO … … 182 182 ! 183 183 184 zvrel(:,:) = 0. 0184 zvrel(:,:) = 0._wp 185 185 186 186 ! Default new ice thickness 187 DO jj = 1, jpj 188 DO ji = 1, jpi 189 hicol(ji,jj) = hiccrit(1) 190 END DO 191 END DO 192 193 IF (fraz_swi.eq.1.0) THEN 187 hicol(:,:) = hiccrit(1) 188 189 IF( fraz_swi == 1._wp ) THEN 194 190 195 191 !-------------------- 196 192 ! Physical constants 197 193 !-------------------- 198 hicol(:,:) = 0. 0194 hicol(:,:) = 0._wp 199 195 200 196 zhicrit = 0.04 ! frazil ice thickness … … 211 207 !------------- 212 208 ! C-grid wind stress components 213 ztaux = ( utau_ice(ji-1,jj ) * tmu(ji-1,jj ) &214 & + utau_ice(ji ,jj ) * tmu(ji ,jj ) ) / 2.0215 ztauy = ( vtau_ice(ji ,jj-1) * tmv(ji ,jj-1) &216 & + vtau_ice(ji ,jj ) * tmv(ji ,jj ) ) / 2.0209 ztaux = ( utau_ice(ji-1,jj ) * tmu(ji-1,jj ) & 210 & + utau_ice(ji ,jj ) * tmu(ji ,jj ) ) * 0.5_wp 211 ztauy = ( vtau_ice(ji ,jj-1) * tmv(ji ,jj-1) & 212 & + vtau_ice(ji ,jj ) * tmv(ji ,jj ) ) * 0.5_wp 217 213 ! Square root of wind stress 218 214 ztenagm = SQRT( SQRT( ztaux * ztaux + ztauy * ztauy ) ) … … 228 224 !------------------- 229 225 ! C-grid ice velocity 230 zindb = MAX( 0.0, SIGN(1.0, at_i(ji,jj) ))231 zvgx = zindb * ( u_ice(ji-1,jj ) * tmu(ji-1,jj )&232 + u_ice(ji,jj ) * tmu(ji ,jj ) ) / 2.0233 zvgy = zindb * ( v_ice(ji ,jj-1) * tmv(ji ,jj-1)&234 + v_ice(ji,jj ) * tmv(ji ,jj ) ) / 2.0226 zindb = MAX( 0._wp, SIGN( 1._wp , at_i(ji,jj) ) ) 227 zvgx = zindb * ( u_ice(ji-1,jj ) * tmu(ji-1,jj ) & 228 & + u_ice(ji,jj ) * tmu(ji ,jj ) ) * 0.5_wp 229 zvgy = zindb * ( v_ice(ji ,jj-1) * tmv(ji ,jj-1) & 230 & + v_ice(ji,jj ) * tmv(ji ,jj ) ) * 0.5_wp 235 231 236 232 !----------------------------------- … … 238 234 !----------------------------------- 239 235 ! absolute relative velocity 240 zvrel2 = MAX( ( zvfrx - zvgx ) * ( zvfrx - zvgx ) + & 241 ( zvfry - zvgy ) * ( zvfry - zvgy ) & 242 , 0.15 * 0.15 ) 243 zvrel(ji,jj) = SQRT(zvrel2) 236 zvrel2 = MAX( ( zvfrx - zvgx ) * ( zvfrx - zvgx ) & 237 & + ( zvfry - zvgy ) * ( zvfry - zvgy ) , 0.15 * 0.15 ) 238 zvrel(ji,jj) = SQRT( zvrel2 ) 244 239 245 240 !--------------------- … … 247 242 !--------------------- 248 243 hicol(ji,jj) = zhicrit + 0.1 249 hicol(ji,jj) = zhicrit + hicol(ji,jj) / & 250 ( hicol(ji,jj) * hicol(ji,jj) - & 251 zhicrit * zhicrit ) * ztwogp * zvrel2 244 hicol(ji,jj) = zhicrit + hicol(ji,jj) & 245 & / ( hicol(ji,jj) * hicol(ji,jj) - zhicrit * zhicrit ) * ztwogp * zvrel2 246 247 !!gm better coding: above: hicol(ji,jj) * hicol(ji,jj) = (zhicrit + 0.1)*(zhicrit + 0.1) 248 !!gm = zhicrit**2 + 0.2*zhicrit +0.01 249 !!gm therefore the 2 lines with hicol can be replaced by 1 line: 250 !!gm hicol(ji,jj) = zhicrit + (zhicrit + 0.1) / ( 0.2 * zhicrit + 0.01 ) * ztwogp * zvrel2 251 !!gm further more (zhicrit + 0.1)/(0.2 * zhicrit + 0.01 )*ztwogp can be computed one for all outside the DO loop 252 252 253 253 iter = 1 … … 284 284 DO jj = 1, jpj 285 285 DO ji = 1, jpi 286 IF ( tms(ji,jj) * ( qcmif(ji,jj) - qldif(ji,jj) ) > 0.e0) THEN286 IF ( tms(ji,jj) * ( qcmif(ji,jj) - qldif(ji,jj) ) > 0._wp ) THEN 287 287 nbpac = nbpac + 1 288 288 npac( nbpac ) = (jj - 1) * jpi + ji 289 IF ( (ji.eq.jiindx).AND.(jj.eq.jjindx) ) THEN 290 jiindex_1d = nbpac 291 ENDIF 289 IF( ji == jiindx .AND. jj == jjindx ) jiindex_1d = nbpac 292 290 ENDIF 293 291 END DO 294 292 END DO 295 293 296 IF( ln_nicep ) THEN 297 WRITE(numout,*) 'lim_thd_lac : nbpac = ', nbpac 298 ENDIF 294 IF( ln_nicep ) WRITE(numout,*) 'lim_thd_lac : nbpac = ', nbpac 299 295 300 296 !------------------------------ … … 306 302 IF ( nbpac > 0 ) THEN 307 303 308 CALL tab_2d_1d( nbpac, zat_i_ac (1:nbpac) , at_i , & 309 jpi, jpj, npac(1:nbpac) ) 304 CALL tab_2d_1d( nbpac, zat_i_ac (1:nbpac) , at_i , jpi, jpj, npac(1:nbpac) ) 310 305 DO jl = 1, jpl 311 CALL tab_2d_1d( nbpac, za_i_ac(1:nbpac,jl) , a_i(:,:,jl) , & 312 jpi, jpj, npac(1:nbpac) ) 313 CALL tab_2d_1d( nbpac, zv_i_ac(1:nbpac,jl) , v_i(:,:,jl) , & 314 jpi, jpj, npac(1:nbpac) ) 315 CALL tab_2d_1d( nbpac, zoa_i_ac(1:nbpac,jl) , oa_i(:,:,jl) , & 316 jpi, jpj, npac(1:nbpac) ) 317 CALL tab_2d_1d( nbpac, zsmv_i_ac(1:nbpac,jl), smv_i(:,:,jl), & 318 jpi, jpj, npac(1:nbpac) ) 306 CALL tab_2d_1d( nbpac, za_i_ac (1:nbpac,jl), a_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 307 CALL tab_2d_1d( nbpac, zv_i_ac (1:nbpac,jl), v_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 308 CALL tab_2d_1d( nbpac, zoa_i_ac (1:nbpac,jl), oa_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 309 CALL tab_2d_1d( nbpac, zsmv_i_ac(1:nbpac,jl), smv_i(:,:,jl), jpi, jpj, npac(1:nbpac) ) 319 310 DO jk = 1, nlay_i 320 CALL tab_2d_1d( nbpac, ze_i_ac(1:nbpac,jk,jl), e_i(:,:,jk,jl) , & 321 jpi, jpj, npac(1:nbpac) ) 311 CALL tab_2d_1d( nbpac, ze_i_ac(1:nbpac,jk,jl), e_i(:,:,jk,jl) , jpi, jpj, npac(1:nbpac) ) 322 312 END DO ! jk 323 313 END DO ! jl 324 314 325 CALL tab_2d_1d( nbpac, qldif_1d (1:nbpac) , qldif , & 326 jpi, jpj, npac(1:nbpac) ) 327 CALL tab_2d_1d( nbpac, qcmif_1d (1:nbpac) , qcmif , & 328 jpi, jpj, npac(1:nbpac) ) 329 CALL tab_2d_1d( nbpac, t_bo_b (1:nbpac) , t_bo , & 330 jpi, jpj, npac(1:nbpac) ) 331 CALL tab_2d_1d( nbpac, fseqv_1d (1:nbpac) , fseqv , & 332 jpi, jpj, npac(1:nbpac) ) 333 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , & 334 jpi, jpj, npac(1:nbpac) ) 335 CALL tab_2d_1d( nbpac, zvrel_ac (1:nbpac) , zvrel , & 336 jpi, jpj, npac(1:nbpac) ) 315 CALL tab_2d_1d( nbpac, qldif_1d (1:nbpac) , qldif , jpi, jpj, npac(1:nbpac) ) 316 CALL tab_2d_1d( nbpac, qcmif_1d (1:nbpac) , qcmif , jpi, jpj, npac(1:nbpac) ) 317 CALL tab_2d_1d( nbpac, t_bo_b (1:nbpac) , t_bo , jpi, jpj, npac(1:nbpac) ) 318 CALL tab_2d_1d( nbpac, sfx_thd_1d(1:nbpac) , sfx_thd, jpi, jpj, npac(1:nbpac) ) 319 CALL tab_2d_1d( nbpac, rdm_ice_1d(1:nbpac) , rdm_ice, jpi, jpj, npac(1:nbpac) ) 320 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , jpi, jpj, npac(1:nbpac) ) 321 CALL tab_2d_1d( nbpac, zvrel_ac (1:nbpac) , zvrel , jpi, jpj, npac(1:nbpac) ) 337 322 338 323 !------------------------------------------------------------------------------! … … 344 329 !---------------------- 345 330 DO ji = 1, nbpac 346 zh_newice(ji) 347 END DO 348 IF ( fraz_swi .EQ. 1.0 )zh_newice(:) = hicol_b(:)331 zh_newice(ji) = hiccrit(1) 332 END DO 333 IF( fraz_swi == 1.0 ) zh_newice(:) = hicol_b(:) 349 334 350 335 !---------------------- … … 352 337 !---------------------- 353 338 354 IF ( num_sal .EQ. 1 ) THEN 355 zs_newice(:) = bulk_sal 356 ENDIF ! num_sal 357 358 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN 359 360 DO ji = 1, nbpac 361 zs_newice(ji) = MIN( 4.606 + 0.91 / zh_newice(ji) , s_i_max ) 362 zji = MOD( npac(ji) - 1, jpi ) + 1 363 zjj = ( npac(ji) - 1 ) / jpi + 1 364 zs_newice(ji) = MIN( 0.5*sss_m(zji,zjj) , zs_newice(ji) ) 365 END DO ! jl 366 367 ENDIF ! num_sal 368 369 IF ( num_sal .EQ. 3 ) THEN 370 zs_newice(:) = 2.3 371 ENDIF ! num_sal 339 SELECT CASE ( num_sal ) 340 CASE ( 1 ) ! Sice = constant 341 zs_newice(:) = bulk_sal 342 CASE ( 2 ) ! Sice = F(z,t) [Vancoppenolle et al (2005)] 343 DO ji = 1, nbpac 344 zji = MOD( npac(ji) - 1 , jpi ) + 1 345 zjj = ( npac(ji) - 1 ) / jpi + 1 346 zs_newice(ji) = MIN( 4.606 + 0.91 / zh_newice(ji) , s_i_max , 0.5 * sss_m(zji,zjj) ) 347 END DO 348 CASE ( 3 ) ! Sice = F(z) [multiyear ice] 349 zs_newice(:) = 2.3 350 END SELECT 351 372 352 373 353 !------------------------- … … 376 356 ! We assume that new ice is formed at the seawater freezing point 377 357 DO ji = 1, nbpac 378 ztmelts = - tmut * zs_newice(ji) + rtt ! Melting point (K) 379 ze_newice(ji) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 380 + lfus * ( 1.0 - ( ztmelts - rtt ) & 381 / ( t_bo_b(ji) - rtt ) ) & 382 - rcp * ( ztmelts-rtt ) ) 383 ze_newice(ji) = MAX( ze_newice(ji) , 0.0 ) + & 384 MAX( 0.0 , SIGN( 1.0 , - ze_newice(ji) ) ) & 385 * rhoic * lfus 358 ztmelts = - tmut * zs_newice(ji) + rtt ! Melting point (K) 359 ze_newice(ji) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 360 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( t_bo_b(ji) - rtt ) ) & 361 & - rcp * ( ztmelts - rtt ) ) 362 ze_newice(ji) = MAX( ze_newice(ji) , 0._wp ) & 363 & + MAX( 0.0 , SIGN( 1.0 , - ze_newice(ji) ) ) * rhoic * lfus 386 364 END DO ! ji 387 365 !---------------- … … 389 367 !---------------- 390 368 DO ji = 1, nbpac 391 zo_newice(ji) = 0.0369 zo_newice(ji) = 0._wp 392 370 END DO ! ji 393 371 … … 396 374 !-------------------------- 397 375 DO ji = 1, nbpac 398 zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji)!<0376 zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji) !<0 399 377 END DO ! ji 400 378 … … 403 381 !------------------- 404 382 DO ji = 1, nbpac 405 zv_newice(ji) 383 zv_newice(ji) = - zqbgow(ji) / ze_newice(ji) 406 384 407 385 ! A fraction zfrazb of frazil ice is accreted at the ice bottom 408 zfrazb = ( TANH ( Cfrazb * ( zvrel_ac(ji) - vfrazb ) ) & 409 + 1.0 ) / 2.0 * maxfrazb 410 zdh_frazb(ji) = zfrazb*zv_newice(ji) 386 zfrazb = ( TANH ( Cfrazb * ( zvrel_ac(ji) - vfrazb ) ) + 1.0 ) * 0.5 * maxfrazb 387 zdh_frazb(ji) = zfrazb * zv_newice(ji) 411 388 zv_newice(ji) = ( 1.0 - zfrazb ) * zv_newice(ji) 412 389 END DO … … 415 392 ! Salt flux due to new ice growth 416 393 !--------------------------------- 417 IF ( ( num_sal .EQ. 4 ) ) THEN 418 DO ji = 1, nbpac 419 zji = MOD( npac(ji) - 1, jpi ) + 1 420 zjj = ( npac(ji) - 1 ) / jpi + 1 421 fseqv_1d(ji) = fseqv_1d(ji) + & 422 ( sss_m(zji,zjj) - bulk_sal ) * rhoic * & 423 zv_newice(ji) / rdt_ice 424 END DO 425 ELSE 426 DO ji = 1, nbpac 427 zji = MOD( npac(ji) - 1, jpi ) + 1 428 zjj = ( npac(ji) - 1 ) / jpi + 1 429 fseqv_1d(ji) = fseqv_1d(ji) + & 430 ( sss_m(zji,zjj) - zs_newice(ji) ) * rhoic * & 431 zv_newice(ji) / rdt_ice 432 END DO ! ji 433 ENDIF 394 ! note that for constant salinity zs_newice() = bulk_sal (see top of the subroutine) 395 DO ji = 1, nbpac 396 sfx_thd_1d(ji) = sfx_thd_1d(ji) - zs_newice(ji) * rhoic * zv_newice(ji) * r1_rdtice 397 rdm_ice_1d(ji) = rdm_ice_1d(ji) + rhoic * zv_newice(ji) 398 END DO ! ji 434 399 435 400 !------------------------------------ … … 437 402 !------------------------------------ 438 403 DO ji = 1, nbpac 439 ! Volume440 zj i = MOD( npac(ji) - 1, jpi )+ 1441 zjj = ( npac(ji) - 1 ) / jpi + 1442 vt_i_init(zji,zjj) = vt_i_init(zji,zjj) +zv_newice(ji)443 ! Energy444 zde = ze_newice(ji) / unit_fac445 zde = zde * area(zji,zjj) * zv_newice(ji)446 et_i_init(zji,zjj) = et_i_init(zji,zjj) + zde 404 zji = MOD( npac(ji) - 1 , jpi ) + 1 405 zjj = ( npac(ji) - 1 ) / jpi + 1 406 ! 407 zde = ze_newice(ji) / unit_fac * area(zji,zjj) * zv_newice(ji) 408 ! 409 vt_i_init(zji,zjj) = vt_i_init(zji,zjj) + zv_newice(ji) ! volume 410 et_i_init(zji,zjj) = et_i_init(zji,zjj) + zde ! Energy 411 447 412 END DO 448 413 449 414 ! keep new ice volume in memory 450 CALL tab_1d_2d( nbpac, v_newice , npac(1:nbpac), zv_newice(1:nbpac) , & 451 jpi, jpj ) 415 CALL tab_1d_2d( nbpac, v_newice , npac(1:nbpac), zv_newice(1:nbpac) , jpi, jpj ) 452 416 453 417 !----------------- … … 455 419 !----------------- 456 420 DO ji = 1, nbpac 457 za_newice(ji) = zv_newice(ji) / zh_newice(ji) 458 ! diagnostic 459 zji = MOD( npac(ji) - 1, jpi ) + 1 460 zjj = ( npac(ji) - 1 ) / jpi + 1 461 &nb