Changeset 9935
- Timestamp:
- 2018-07-12T16:12:48+02:00 (5 years ago)
- Location:
- NEMO/trunk/src
- Files:
-
- 24 edited
-
ICE/icecor.F90 (modified) (3 diffs)
-
ICE/icectl.F90 (modified) (4 diffs)
-
ICE/icedia.F90 (modified) (2 diffs)
-
ICE/icedyn_adv.F90 (modified) (1 diff)
-
ICE/icedyn_rdgrft.F90 (modified) (4 diffs)
-
ICE/icedyn_rhg_evp.F90 (modified) (2 diffs)
-
ICE/iceistate.F90 (modified) (4 diffs)
-
ICE/icethd.F90 (modified) (1 diff)
-
ICE/icethd_da.F90 (modified) (2 diffs)
-
ICE/icethd_dh.F90 (modified) (24 diffs)
-
ICE/icethd_do.F90 (modified) (5 diffs)
-
ICE/icethd_pnd.F90 (modified) (4 diffs)
-
ICE/icethd_sal.F90 (modified) (2 diffs)
-
ICE/icethd_zdf_bl99.F90 (modified) (7 diffs)
-
ICE/iceupdate.F90 (modified) (2 diffs)
-
ICE/icevar.F90 (modified) (6 diffs)
-
ICE/icewri.F90 (modified) (5 diffs)
-
OCE/BDY/bdyice.F90 (modified) (2 diffs)
-
OCE/DOM/phycst.F90 (modified) (3 diffs)
-
OCE/SBC/sbcblk.F90 (modified) (10 diffs)
-
OCE/SBC/sbccpl.F90 (modified) (7 diffs)
-
OCE/SBC/sbcice_cice.F90 (modified) (5 diffs)
-
OCE/SBC/sbcisf.F90 (modified) (7 diffs)
-
OCE/SBC/sbcrnf.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icecor.F90
r9604 r9935 86 86 IF ( nn_icesal == 2 ) THEN ! salinity must stay in bounds [Simin,Simax] ! 87 87 ! !----------------------------------------------------- 88 zzc = rhoi c* r1_rdtice88 zzc = rhoi * r1_rdtice 89 89 DO jl = 1, jpl 90 90 DO jj = 1, jpj … … 139 139 & + SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) ) ) * r1_rdtice 140 140 ! ! salt, volume 141 diag_sice(ji,jj) = SUM( sv_i(ji,jj,:) - sv_i_b(ji,jj,:) ) * rhoi c* r1_rdtice142 diag_vice(ji,jj) = SUM( v_i (ji,jj,:) - v_i_b (ji,jj,:) ) * rhoi c* r1_rdtice143 diag_vsnw(ji,jj) = SUM( v_s (ji,jj,:) - v_s_b (ji,jj,:) ) * rhos n* r1_rdtice141 diag_sice(ji,jj) = SUM( sv_i(ji,jj,:) - sv_i_b(ji,jj,:) ) * rhoi * r1_rdtice 142 diag_vice(ji,jj) = SUM( v_i (ji,jj,:) - v_i_b (ji,jj,:) ) * rhoi * r1_rdtice 143 diag_vsnw(ji,jj) = SUM( v_s (ji,jj,:) - v_s_b (ji,jj,:) ) * rhos * r1_rdtice 144 144 END DO 145 145 END DO … … 160 160 & + SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) ) ) * r1_rdtice 161 161 ! ! salt, volume 162 diag_sice(ji,jj) = diag_sice(ji,jj) + SUM( sv_i(ji,jj,:) - sv_i_b(ji,jj,:) ) * rhoi c* r1_rdtice163 diag_vice(ji,jj) = diag_vice(ji,jj) + SUM( v_i (ji,jj,:) - v_i_b (ji,jj,:) ) * rhoi c* r1_rdtice164 diag_vsnw(ji,jj) = diag_vsnw(ji,jj) + SUM( v_s (ji,jj,:) - v_s_b (ji,jj,:) ) * rhos n* r1_rdtice162 diag_sice(ji,jj) = diag_sice(ji,jj) + SUM( sv_i(ji,jj,:) - sv_i_b(ji,jj,:) ) * rhoi * r1_rdtice 163 diag_vice(ji,jj) = diag_vice(ji,jj) + SUM( v_i (ji,jj,:) - v_i_b (ji,jj,:) ) * rhoi * r1_rdtice 164 diag_vsnw(ji,jj) = diag_vsnw(ji,jj) + SUM( v_s (ji,jj,:) - v_s_b (ji,jj,:) ) * rhos * r1_rdtice 165 165 END DO 166 166 END DO -
NEMO/trunk/src/ICE/icectl.F90
r9913 r9935 93 93 & ) * e1e2t(:,:) ) * zconv 94 94 95 pdiag_v = glob_sum( SUM( v_i * rhoi c + v_s * rhosn, dim=3 ) * e1e2t * zconv )96 97 pdiag_s = glob_sum( SUM( sv_i * rhoi c, dim=3 ) * e1e2t * zconv )95 pdiag_v = glob_sum( SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t * zconv ) 96 97 pdiag_s = glob_sum( SUM( sv_i * rhoi , dim=3 ) * e1e2t * zconv ) 98 98 99 99 pdiag_t = glob_sum( ( SUM( SUM( e_i(:,:,1:nlay_i,:), dim=4 ), dim=3 ) & … … 120 120 121 121 ! outputs 122 zv = ( ( glob_sum( SUM( v_i * rhoi c + v_s * rhosn, dim=3 ) * e1e2t ) * zconv &122 zv = ( ( glob_sum( SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t ) * zconv & 123 123 & - pdiag_v ) * r1_rdtice - zfv ) * rday 124 124 125 zs = ( ( glob_sum( SUM( sv_i * rhoi c, dim=3 ) * e1e2t ) * zconv &125 zs = ( ( glob_sum( SUM( sv_i * rhoi , dim=3 ) * e1e2t ) * zconv & 126 126 & - pdiag_s ) * r1_rdtice + zfs ) * rday 127 127 … … 131 131 132 132 ! zvtrp and zetrp must be close to 0 if the advection scheme is conservative 133 zvtrp = glob_sum( ( diag_trp_vi * rhoi c + diag_trp_vs * rhosn) * e1e2t ) * zconv * rday134 zetrp = glob_sum( ( diag_trp_ei + diag_trp_es) * e1e2t ) * zconv133 zvtrp = glob_sum( ( diag_trp_vi * rhoi + diag_trp_vs * rhos ) * e1e2t ) * zconv * rday 134 zetrp = glob_sum( ( diag_trp_ei + diag_trp_es ) * e1e2t ) * zconv 135 135 136 136 zvmin = glob_min( v_i ) … … 382 382 DO jj = 1, jpj 383 383 DO ji = 1, jpi 384 ztmelts = - tmut * sz_i(ji,jj,jk,jl) + rt0384 ztmelts = -rTmlt * sz_i(ji,jj,jk,jl) + rt0 385 385 IF( t_i(ji,jj,jk,jl) >= ztmelts .AND. v_i(ji,jj,jl) > 1.e-10 & 386 386 & .AND. a_i(ji,jj,jl) > 0._wp ) THEN -
NEMO/trunk/src/ICE/icedia.F90
r9912 r9935 110 110 ! 3 - Content variations ! 111 111 ! ----------------------- ! 112 zdiff_vol = r1_rau0 * glob_sum( ( rhoi c*vt_i(:,:) + rhosn*vt_s(:,:) - vol_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater trend (km3)113 zdiff_sal = r1_rau0 * glob_sum( ( rhoi c* SUM( sv_i(:,:,:), dim=3 )- sal_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! salt content trend (km3*pss)114 zdiff_tem = glob_sum( ( et_i(:,:) + et_s(:,:) 112 zdiff_vol = r1_rau0 * glob_sum( ( rhoi*vt_i(:,:) + rhos*vt_s(:,:) - vol_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater trend (km3) 113 zdiff_sal = r1_rau0 * glob_sum( ( rhoi* SUM( sv_i(:,:,:), dim=3 ) - sal_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! salt content trend (km3*pss) 114 zdiff_tem = glob_sum( ( et_i(:,:) + et_s(:,:) - tem_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-20 ! heat content trend (1.e20 J) 115 115 ! + SUM( qevap_ice * a_i_b, dim=3 ) !! clem: I think this term should not be there (but needs a check) 116 116 … … 246 246 frc_sal = 0._wp 247 247 ! record initial ice volume, salt and temp 248 vol_loc_ini(:,:) = rhoi c * vt_i(:,:) + rhosn* vt_s(:,:) ! ice/snow volume (kg/m2)249 tem_loc_ini(:,:) = et_i(:,:) + et_s(:,:) 250 sal_loc_ini(:,:) = rhoi c * SUM( sv_i(:,:,:), dim=3 )! ice salt content (pss*kg/m2)248 vol_loc_ini(:,:) = rhoi * vt_i(:,:) + rhos * vt_s(:,:) ! ice/snow volume (kg/m2) 249 tem_loc_ini(:,:) = et_i(:,:) + et_s(:,:) ! ice/snow heat content (J) 250 sal_loc_ini(:,:) = rhoi * SUM( sv_i(:,:,:), dim=3 ) ! ice salt content (pss*kg/m2) 251 251 ENDIF 252 252 ! -
NEMO/trunk/src/ICE/icedyn_adv.F90
r9888 r9935 119 119 diag_trp_vi(:,:) = SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_rdtice 120 120 diag_trp_vs(:,:) = SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_rdtice 121 IF( iom_use('icemtrp') ) CALL iom_put( "icemtrp" , diag_trp_vi * rhoi c) ! ice mass transport122 IF( iom_use('snwmtrp') ) CALL iom_put( "snwmtrp" , diag_trp_vs * rhos n) ! snw mass transport123 IF( iom_use('salmtrp') ) CALL iom_put( "salmtrp" , diag_trp_sv * rhoi c* 1.e-03 ) ! salt mass transport (kg/m2/s)124 IF( iom_use('dihctrp') ) CALL iom_put( "dihctrp" , -diag_trp_ei 125 IF( iom_use('dshctrp') ) CALL iom_put( "dshctrp" , -diag_trp_es 121 IF( iom_use('icemtrp') ) CALL iom_put( "icemtrp" , diag_trp_vi * rhoi ) ! ice mass transport 122 IF( iom_use('snwmtrp') ) CALL iom_put( "snwmtrp" , diag_trp_vs * rhos ) ! snw mass transport 123 IF( iom_use('salmtrp') ) CALL iom_put( "salmtrp" , diag_trp_sv * rhoi * 1.e-03 ) ! salt mass transport (kg/m2/s) 124 IF( iom_use('dihctrp') ) CALL iom_put( "dihctrp" , -diag_trp_ei ) ! advected ice heat content (W/m2) 125 IF( iom_use('dshctrp') ) CALL iom_put( "dshctrp" , -diag_trp_es ) ! advected snw heat content (W/m2) 126 126 127 127 ! controls -
NEMO/trunk/src/ICE/icedyn_rdgrft.F90
r9880 r9935 547 547 ! volume and enthalpy (J/m2, >0) of seawater trapped into ridges 548 548 vsw = v_i_2d(ji,jl1) * afrdg * rn_porordg 549 ersw(ji) = -rhoi c* vsw * rcp * sst_1d(ji) ! clem: if sst>0, then ersw <0 (is that possible?)549 ersw(ji) = -rhoi * vsw * rcp * sst_1d(ji) ! clem: if sst>0, then ersw <0 (is that possible?) 550 550 551 551 ! volume etc of ridging / rafting ice and new ridges (vi, vs, sm, oi, es, ei) … … 574 574 575 575 ! Ice-ocean exchanges associated with ice porosity 576 wfx_dyn_1d(ji) = wfx_dyn_1d(ji) - vsw * rhoi c* r1_rdtice ! increase in ice volume due to seawater frozen in voids577 sfx_dyn_1d(ji) = sfx_dyn_1d(ji) - vsw * sss_1d(ji) * rhoi c* r1_rdtice576 wfx_dyn_1d(ji) = wfx_dyn_1d(ji) - vsw * rhoi * r1_rdtice ! increase in ice volume due to seawater frozen in voids 577 sfx_dyn_1d(ji) = sfx_dyn_1d(ji) - vsw * sss_1d(ji) * rhoi * r1_rdtice 578 578 hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ersw(ji) * r1_rdtice ! > 0 [W.m-2] 579 579 580 580 ! Put the snow lost by ridging into the ocean 581 581 ! Note that esrdg > 0; the ocean must cool to melt snow. If the ocean temp = Tf already, new ice must grow. 582 wfx_snw_dyn_1d(ji) = wfx_snw_dyn_1d(ji) + ( rhos n* vsrdg(ji) * ( 1._wp - rn_fsnwrdg ) & ! fresh water source for ocean583 & + rhos n* vsrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_rdtice582 wfx_snw_dyn_1d(ji) = wfx_snw_dyn_1d(ji) + ( rhos * vsrdg(ji) * ( 1._wp - rn_fsnwrdg ) & ! fresh water source for ocean 583 & + rhos * vsrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_rdtice 584 584 585 585 ! Put the melt pond water into the ocean … … 587 587 ! is no net mass flux between melt ponds and the ocean (see icethd_pnd.F90 for ex.) 588 588 !IF ( ln_pnd_fwb ) THEN 589 ! wfx_pnd_1d(ji) = wfx_pnd_1d(ji) + ( rho fw * vprdg(ji) * ( 1._wp - rn_fpndrdg ) & ! fresh water source for ocean590 ! & + rho fw * vprft(ji) * ( 1._wp - rn_fpndrft ) ) * r1_rdtice589 ! wfx_pnd_1d(ji) = wfx_pnd_1d(ji) + ( rhow * vprdg(ji) * ( 1._wp - rn_fpndrdg ) & ! fresh water source for ocean 590 ! & + rhow * vprft(ji) * ( 1._wp - rn_fpndrft ) ) * r1_rdtice 591 591 !ENDIF 592 592 … … 594 594 IF( nn_icesal /= 2 ) THEN 595 595 sirdg2(ji) = sirdg2(ji) - vsw * ( sss_1d(ji) - s_i_1d(ji) ) ! ridge salinity = s_i 596 sfx_bri_1d(ji) = sfx_bri_1d(ji) + sss_1d(ji) * vsw * rhoi c* r1_rdtice & ! put back sss_m into the ocean597 & - s_i_1d(ji) * vsw * rhoi c* r1_rdtice ! and get s_i from the ocean596 sfx_bri_1d(ji) = sfx_bri_1d(ji) + sss_1d(ji) * vsw * rhoi * r1_rdtice & ! put back sss_m into the ocean 597 & - s_i_1d(ji) * vsw * rhoi * r1_rdtice ! and get s_i from the ocean 598 598 ENDIF 599 599 -
NEMO/trunk/src/ICE/icedyn_rhg_evp.F90
r9866 r9935 285 285 286 286 ! Ice/snow mass at U-V points 287 zm1 = ( rhos n * vt_s(ji ,jj ) + rhoic* vt_i(ji ,jj ) )288 zm2 = ( rhos n * vt_s(ji+1,jj ) + rhoic* vt_i(ji+1,jj ) )289 zm3 = ( rhos n * vt_s(ji ,jj+1) + rhoic* vt_i(ji ,jj+1) )287 zm1 = ( rhos * vt_s(ji ,jj ) + rhoi * vt_i(ji ,jj ) ) 288 zm2 = ( rhos * vt_s(ji+1,jj ) + rhoi * vt_i(ji+1,jj ) ) 289 zm3 = ( rhos * vt_s(ji ,jj+1) + rhoi * vt_i(ji ,jj+1) ) 290 290 zmassU = 0.5_wp * ( zm1 * e1e2t(ji,jj) + zm2 * e1e2t(ji+1,jj) ) * r1_e1e2u(ji,jj) * umask(ji,jj,1) 291 291 zmassV = 0.5_wp * ( zm1 * e1e2t(ji,jj) + zm3 * e1e2t(ji,jj+1) ) * r1_e1e2v(ji,jj) * vmask(ji,jj,1) … … 795 795 zfac_y = 0.5 * v_ice(ji,jj) * e1v(ji,jj) * rswitch 796 796 797 zdiag_xmtrp_ice(ji,jj) = rhoi c* zfac_x * ( vt_i(ji+1,jj) + vt_i(ji,jj) ) ! ice mass transport, X-component798 zdiag_ymtrp_ice(ji,jj) = rhoi c* zfac_y * ( vt_i(ji,jj+1) + vt_i(ji,jj) ) ! '' Y- ''799 800 zdiag_xmtrp_snw(ji,jj) = rhos n* zfac_x * ( vt_s(ji+1,jj) + vt_s(ji,jj) ) ! snow mass transport, X-component801 zdiag_ymtrp_snw(ji,jj) = rhos n* zfac_y * ( vt_s(ji,jj+1) + vt_s(ji,jj) ) ! '' Y- ''802 803 zdiag_xatrp(ji,jj) = zfac_x * ( at_i(ji+1,jj) + at_i(ji,jj) ) 804 zdiag_yatrp(ji,jj) = zfac_y * ( at_i(ji,jj+1) + at_i(ji,jj) ) 797 zdiag_xmtrp_ice(ji,jj) = rhoi * zfac_x * ( vt_i(ji+1,jj) + vt_i(ji,jj) ) ! ice mass transport, X-component 798 zdiag_ymtrp_ice(ji,jj) = rhoi * zfac_y * ( vt_i(ji,jj+1) + vt_i(ji,jj) ) ! '' Y- '' 799 800 zdiag_xmtrp_snw(ji,jj) = rhos * zfac_x * ( vt_s(ji+1,jj) + vt_s(ji,jj) ) ! snow mass transport, X-component 801 zdiag_ymtrp_snw(ji,jj) = rhos * zfac_y * ( vt_s(ji,jj+1) + vt_s(ji,jj) ) ! '' Y- '' 802 803 zdiag_xatrp(ji,jj) = zfac_x * ( at_i(ji+1,jj) + at_i(ji,jj) ) ! area transport, X-component 804 zdiag_yatrp(ji,jj) = zfac_y * ( at_i(ji,jj+1) + at_i(ji,jj) ) ! '' Y- '' 805 805 806 806 END DO -
NEMO/trunk/src/ICE/iceistate.F90
r9656 r9935 295 295 ! In case snow load is in excess that would lead to transformation from snow to ice 296 296 ! Then, transfer the snow excess into the ice (different from icethd_dh) 297 zdh = MAX( 0._wp, ( rhos n * h_s(ji,jj,jl) + ( rhoic- rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 )297 zdh = MAX( 0._wp, ( rhos * h_s(ji,jj,jl) + ( rhoi - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 ) 298 298 ! recompute h_i, h_s avoiding out of bounds values 299 299 h_i(ji,jj,jl) = MIN( hi_max(jl), h_i(ji,jj,jl) + zdh ) 300 h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoi c * r1_rhosn)300 h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoi * r1_rhos ) 301 301 ! 302 302 ! ice volume, salt content, age content … … 321 321 t_s(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0 322 322 ! Snow energy of melting 323 e_s(ji,jj,jk,jl) = zswitch(ji,jj) * rhos n * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus )323 e_s(ji,jj,jk,jl) = zswitch(ji,jj) * rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus ) 324 324 ! 325 325 ! Mutliply by volume, and divide by number of layers to get heat content in J/m2 … … 337 337 t_i (ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0 338 338 sz_i(ji,jj,jk,jl) = zswitch(ji,jj) * zsm_i_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rn_simin 339 ztmelts = - tmut * sz_i(ji,jj,jk,jl) + rt0 !Melting temperature in K339 ztmelts = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 !Melting temperature in K 340 340 ! 341 341 ! heat content per unit volume 342 e_i(ji,jj,jk,jl) = zswitch(ji,jj) * rhoi c * ( cpic * ( ztmelts - t_i(ji,jj,jk,jl) )&343 & + lfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0) , -epsi20 ) ) &342 e_i(ji,jj,jk,jl) = zswitch(ji,jj) * rhoi * ( rcpi * ( ztmelts - t_i(ji,jj,jk,jl) ) & 343 & + rLfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0) , -epsi20 ) ) & 344 344 & - rcp * ( ztmelts - rt0 ) ) 345 345 ! … … 410 410 ! 5) Snow-ice mass (case ice is fully embedded) 411 411 !---------------------------------------------- 412 snwice_mass (:,:) = tmask(:,:,1) * SUM( rhos n * v_s(:,:,:) + rhoic* v_i(:,:,:), dim=3 ) ! snow+ice mass412 snwice_mass (:,:) = tmask(:,:,1) * SUM( rhos * v_s(:,:,:) + rhoi * v_i(:,:,:), dim=3 ) ! snow+ice mass 413 413 snwice_mass_b(:,:) = snwice_mass(:,:) 414 414 ! -
NEMO/trunk/src/ICE/icethd.F90
r9922 r9935 288 288 DO jk = 1, nlay_i 289 289 DO ji = 1, npti 290 ztmelts = - tmut * sz_i_1d(ji,jk)290 ztmelts = -rTmlt * sz_i_1d(ji,jk) 291 291 ! Conversion q(S,T) -> T (second order equation) 292 zbbb = ( rcp - cpic ) * ztmelts + e_i_1d(ji,jk) * r1_rhoic - lfus293 zccc = SQRT( MAX( zbbb * zbbb - 4._wp * cpic * lfus * ztmelts, 0._wp ) )294 t_i_1d(ji,jk) = rt0 - ( zbbb + zccc ) * 0.5_wp * r1_ cpic292 zbbb = ( rcp - rcpi ) * ztmelts + e_i_1d(ji,jk) * r1_rhoi - rLfus 293 zccc = SQRT( MAX( zbbb * zbbb - 4._wp * rcpi * rLfus * ztmelts, 0._wp ) ) 294 t_i_1d(ji,jk) = rt0 - ( zbbb + zccc ) * 0.5_wp * r1_rcpi 295 295 296 296 ! mask temperature -
NEMO/trunk/src/ICE/icethd_da.F90
r9604 r9935 137 137 138 138 ! Contribution to salt flux 139 sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoi c* h_i_1d(ji) * zda * s_i_1d(ji) * r1_rdtice139 sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoi * h_i_1d(ji) * zda * s_i_1d(ji) * r1_rdtice 140 140 141 141 ! Contribution to heat flux into the ocean [W.m-2], (<0) … … 144 144 145 145 ! Contribution to mass flux 146 wfx_lam_1d(ji) = wfx_lam_1d(ji) + zda * r1_rdtice * ( rhoi c * h_i_1d(ji) + rhosn* h_s_1d(ji) )146 wfx_lam_1d(ji) = wfx_lam_1d(ji) + zda * r1_rdtice * ( rhoi * h_i_1d(ji) + rhos * h_s_1d(ji) ) 147 147 148 148 ! new concentration -
NEMO/trunk/src/ICE/icethd_dh.F90
r9922 r9935 76 76 REAL(wp) :: zgrr ! bottom growth rate 77 77 REAL(wp) :: zt_i_new ! bottom formation temperature 78 REAL(wp) :: z1_rho ! 1/(rhos n+rau0-rhoic)78 REAL(wp) :: z1_rho ! 1/(rhos+rau0-rhoi) 79 79 80 80 REAL(wp) :: zQm ! enthalpy exchanged with the ocean (J/m2), >0 towards the ocean … … 174 174 IF( t_s_1d(ji,jk) > rt0 ) THEN 175 175 hfx_res_1d (ji) = hfx_res_1d (ji) + e_s_1d(ji,jk) * zh_s(ji,jk) * a_i_1d(ji) * r1_rdtice ! heat flux to the ocean [W.m-2], < 0 176 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhos n* zh_s(ji,jk) * a_i_1d(ji) * r1_rdtice ! mass flux176 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhos * zh_s(ji,jk) * a_i_1d(ji) * r1_rdtice ! mass flux 177 177 ! updates 178 178 dh_s_mlt(ji) = dh_s_mlt(ji) - zh_s(ji,jk) … … 194 194 ! 195 195 ! --- precipitation --- 196 zdh_s_pre (ji) = zsnw(ji) * sprecip_1d(ji) * rdt_ice * r1_rhos n/ at_i_1d(ji) ! thickness change196 zdh_s_pre (ji) = zsnw(ji) * sprecip_1d(ji) * rdt_ice * r1_rhos / at_i_1d(ji) ! thickness change 197 197 zqprec (ji) = - qprec_ice_1d(ji) ! enthalpy of the precip (>0, J.m-3) 198 198 ! 199 199 hfx_spr_1d(ji) = hfx_spr_1d(ji) + zdh_s_pre(ji) * a_i_1d(ji) * zqprec(ji) * r1_rdtice ! heat flux from snow precip (>0, W.m-2) 200 wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhos n* a_i_1d(ji) * zdh_s_pre(ji) * r1_rdtice ! mass flux, <0200 wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhos * a_i_1d(ji) * zdh_s_pre(ji) * r1_rdtice ! mass flux, <0 201 201 202 202 ! --- melt of falling snow --- … … 205 205 zdeltah (ji,1) = MAX( - zdh_s_pre(ji), zdeltah(ji,1) ) ! bound melting 206 206 hfx_snw_1d (ji) = hfx_snw_1d (ji) - zdeltah(ji,1) * a_i_1d(ji) * zqprec(ji) * r1_rdtice ! heat used to melt snow (W.m-2, >0) 207 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos n* a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice ! snow melting only = water into the ocean (then without snow precip), >0207 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice ! snow melting only = water into the ocean (then without snow precip), >0 208 208 209 209 ! updates available heat + precipitations after melting … … 245 245 246 246 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_1d(ji) * e_s_1d (ji,jk) * r1_rdtice ! heat used to melt snow(W.m-2, >0) 247 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos n* a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! snow melting only = water into the ocean (then without snow precip)247 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! snow melting only = water into the ocean (then without snow precip) 248 248 249 249 ! updates available heat + thickness … … 265 265 IF( evap_ice_1d(ji) > 0._wp ) THEN 266 266 ! 267 zdh_s_sub (ji) = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhos n* rdt_ice )268 zevap_rema(ji) = evap_ice_1d(ji) * rdt_ice + zdh_s_sub(ji) * rhos n! remaining evap in kg.m-2 (used for ice melting later on)267 zdh_s_sub (ji) = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhos * rdt_ice ) 268 zevap_rema(ji) = evap_ice_1d(ji) * rdt_ice + zdh_s_sub(ji) * rhos ! remaining evap in kg.m-2 (used for ice melting later on) 269 269 zdeltah (ji,1) = MAX( zdh_s_sub(ji), - zdh_s_pre(ji) ) 270 270 … … 272 272 & ( zdeltah(ji,1) * zqprec(ji) + ( zdh_s_sub(ji) - zdeltah(ji,1) ) * e_s_1d(ji,1) ) & 273 273 & * a_i_1d(ji) * r1_rdtice 274 wfx_snw_sub_1d(ji) = wfx_snw_sub_1d(ji) - rhos n* a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice ! Mass flux by sublimation274 wfx_snw_sub_1d(ji) = wfx_snw_sub_1d(ji) - rhos * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice ! Mass flux by sublimation 275 275 276 276 ! new snow thickness … … 301 301 e_s_1d(ji,jk) = rswitch / MAX( h_s_1d(ji), epsi20 ) * & 302 302 & ( ( zdh_s_pre(ji) ) * zqprec(ji) + & 303 & ( h_s_1d(ji) - zdh_s_pre(ji) ) * rhos n * ( cpic * ( rt0 - t_s_1d(ji,jk) ) + lfus ) )303 & ( h_s_1d(ji) - zdh_s_pre(ji) ) * rhos * ( rcpi * ( rt0 - t_s_1d(ji,jk) ) + rLfus ) ) 304 304 END DO 305 305 END DO … … 314 314 DO jk = 1, nlay_i 315 315 DO ji = 1, npti 316 ztmelts = - tmut * sz_i_1d(ji,jk) ! Melting point of layer k [C]316 ztmelts = - rTmlt * sz_i_1d(ji,jk) ! Melting point of layer k [C] 317 317 318 318 IF( t_i_1d(ji,jk) >= (ztmelts+rt0) ) THEN !-- Internal melting 319 319 320 zEi = - e_i_1d(ji,jk) * r1_rhoi c! Specific enthalpy of layer k [J/kg, <0]320 zEi = - e_i_1d(ji,jk) * r1_rhoi ! Specific enthalpy of layer k [J/kg, <0] 321 321 zdE = 0._wp ! Specific enthalpy difference (J/kg, <0) 322 322 ! set up at 0 since no energy is needed to melt water...(it is already melted) 323 323 zdeltah(ji,jk) = MIN( 0._wp , - zh_i(ji,jk) ) ! internal melting occurs when the internal temperature is above freezing 324 324 ! this should normally not happen, but sometimes, heat diffusion leads to this 325 zfmdt = - zdeltah(ji,jk) * rhoi c! Mass flux x time step > 0325 zfmdt = - zdeltah(ji,jk) * rhoi ! Mass flux x time step > 0 326 326 327 327 dh_i_itm(ji) = dh_i_itm(ji) + zdeltah(ji,jk) ! Cumulate internal melting 328 328 329 zfmdt = - rhoi c * zdeltah(ji,jk)! Recompute mass flux [kg/m2, >0]329 zfmdt = - rhoi * zdeltah(ji,jk) ! Recompute mass flux [kg/m2, >0] 330 330 331 331 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_rdtice ! Heat flux to the ocean [W.m-2], <0 332 332 ! ice enthalpy zEi is "sent" to the ocean 333 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice! Salt flux333 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice ! Salt flux 334 334 ! using s_i_1d and not sz_i_1d(jk) is ok 335 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice! Mass flux335 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! Mass flux 336 336 337 337 ELSE !-- Surface melting 338 338 339 zEi = - e_i_1d(ji,jk) * r1_rhoi c! Specific enthalpy of layer k [J/kg, <0]339 zEi = - e_i_1d(ji,jk) * r1_rhoi ! Specific enthalpy of layer k [J/kg, <0] 340 340 zEw = rcp * ztmelts ! Specific enthalpy of resulting meltwater [J/kg, <0] 341 341 zdE = zEi - zEw ! Specific enthalpy difference < 0 … … 343 343 zfmdt = - zq_top(ji) / zdE ! Mass flux to the ocean [kg/m2, >0] 344 344 345 zdeltah(ji,jk) = - zfmdt * r1_rhoi c! Melt of layer jk [m, <0]345 zdeltah(ji,jk) = - zfmdt * r1_rhoi ! Melt of layer jk [m, <0] 346 346 347 347 zdeltah(ji,jk) = MIN( 0._wp , MAX( zdeltah(ji,jk) , - zh_i(ji,jk) ) ) ! Melt of layer jk cannot exceed the layer thickness [m, <0] 348 348 349 zq_top(ji) = MAX( 0._wp , zq_top(ji) - zdeltah(ji,jk) * rhoi c* zdE ) ! update available heat349 zq_top(ji) = MAX( 0._wp , zq_top(ji) - zdeltah(ji,jk) * rhoi * zdE ) ! update available heat 350 350 351 351 dh_i_sum(ji) = dh_i_sum(ji) + zdeltah(ji,jk) ! Cumulate surface melt 352 352 353 zfmdt = - rhoi c * zdeltah(ji,jk)! Recompute mass flux [kg/m2, >0]353 zfmdt = - rhoi * zdeltah(ji,jk) ! Recompute mass flux [kg/m2, >0] 354 354 355 355 zQm = zfmdt * zEw ! Energy of the melt water sent to the ocean [J/m2, <0] 356 356 357 sfx_sum_1d(ji) = sfx_sum_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice! Salt flux >0357 sfx_sum_1d(ji) = sfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice ! Salt flux >0 358 358 ! using s_i_1d and not sz_i_1d(jk) is ok) 359 359 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice ! Heat flux [W.m-2], < 0 360 360 hfx_sum_1d(ji) = hfx_sum_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice ! Heat flux used in this process [W.m-2], > 0 361 361 ! 362 wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice! Mass flux362 wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! Mass flux 363 363 364 364 END IF … … 366 366 ! Ice sublimation 367 367 ! --------------- 368 zdum = MAX( - ( zh_i(ji,jk) + zdeltah(ji,jk) ) , - zevap_rema(ji) * r1_rhoi c)368 zdum = MAX( - ( zh_i(ji,jk) + zdeltah(ji,jk) ) , - zevap_rema(ji) * r1_rhoi ) 369 369 zdeltah (ji,jk) = zdeltah (ji,jk) + zdum 370 370 dh_i_sub(ji) = dh_i_sub(ji) + zdum 371 371 372 sfx_sub_1d(ji) = sfx_sub_1d(ji) - rhoi c * a_i_1d(ji) * zdum * s_i_1d(ji) * r1_rdtice! Salt flux >0372 sfx_sub_1d(ji) = sfx_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * s_i_1d(ji) * r1_rdtice ! Salt flux >0 373 373 ! clem: flux is sent to the ocean for simplicity 374 374 ! but salt should remain in the ice except … … 376 376 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zdum * e_i_1d(ji,jk) * a_i_1d(ji) * r1_rdtice ! Heat flux [W.m-2], < 0 377 377 378 wfx_ice_sub_1d(ji) = wfx_ice_sub_1d(ji) - rhoi c * a_i_1d(ji) * zdum * r1_rdtice! Mass flux > 0378 wfx_ice_sub_1d(ji) = wfx_ice_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * r1_rdtice ! Mass flux > 0 379 379 380 380 ! update remaining mass flux 381 zevap_rema(ji) = zevap_rema(ji) + zdum * rhoi c381 zevap_rema(ji) = zevap_rema(ji) + zdum * rhoi 382 382 383 383 ! record which layers have disappeared (for bottom melting) … … 438 438 s_i_new(ji) = zswitch_sal * zfracs * sss_1d(ji) + ( 1. - zswitch_sal ) * s_i_1d(ji) ! New ice salinity 439 439 440 ztmelts = - tmut * s_i_new(ji)! New ice melting point (C)440 ztmelts = - rTmlt * s_i_new(ji) ! New ice melting point (C) 441 441 442 442 zt_i_new = zswitch_sal * t_bo_1d(ji) + ( 1. - zswitch_sal) * t_i_1d(ji, nlay_i) 443 443 444 zEi = cpic* ( zt_i_new - (ztmelts+rt0) ) & ! Specific enthalpy of forming ice (J/kg, <0)445 & - lfus * ( 1.0 - ztmelts / ( zt_i_new - rt0 ) ) + rcp * ztmelts444 zEi = rcpi * ( zt_i_new - (ztmelts+rt0) ) & ! Specific enthalpy of forming ice (J/kg, <0) 445 & - rLfus * ( 1.0 - ztmelts / ( zt_i_new - rt0 ) ) + rcp * ztmelts 446 446 447 447 zEw = rcp * ( t_bo_1d(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0) … … 449 449 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) 450 450 451 dh_i_bog(ji) = rdt_ice * MAX( 0._wp , zf_tt(ji) / ( zdE * rhoi c) )451 dh_i_bog(ji) = rdt_ice * MAX( 0._wp , zf_tt(ji) / ( zdE * rhoi ) ) 452 452 453 453 END DO 454 454 ! Contribution to Energy and Salt Fluxes 455 zfmdt = - rhoi c* dh_i_bog(ji) ! Mass flux x time step (kg/m2, < 0)455 zfmdt = - rhoi * dh_i_bog(ji) ! Mass flux x time step (kg/m2, < 0) 456 456 457 457 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_rdtice ! Heat flux to the ocean [W.m-2], >0 458 458 hfx_bog_1d(ji) = hfx_bog_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice ! Heat flux used in this process [W.m-2], <0 459 459 460 sfx_bog_1d(ji) = sfx_bog_1d(ji) - rhoi c * a_i_1d(ji) * dh_i_bog(ji) * s_i_new(ji) * r1_rdtice! Salt flux, <0461 462 wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoi c * a_i_1d(ji) * dh_i_bog(ji) * r1_rdtice! Mass flux, <0460 sfx_bog_1d(ji) = sfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * s_i_new(ji) * r1_rdtice ! Salt flux, <0 461 462 wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * r1_rdtice ! Mass flux, <0 463 463 464 464 ! update heat content (J.m-2) and layer thickness 465 eh_i_old(ji,nlay_i+1) = eh_i_old(ji,nlay_i+1) + dh_i_bog(ji) * (-zEi * rhoi c)465 eh_i_old(ji,nlay_i+1) = eh_i_old(ji,nlay_i+1) + dh_i_bog(ji) * (-zEi * rhoi) 466 466 h_i_old (ji,nlay_i+1) = h_i_old (ji,nlay_i+1) + dh_i_bog(ji) 467 467 … … 477 477 IF( zf_tt(ji) > 0._wp .AND. jk > icount(ji,jk) ) THEN ! do not calculate where layer has already disappeared by surface melting 478 478 479 ztmelts = - tmut * sz_i_1d(ji,jk) ! Melting point of layer jk (C)479 ztmelts = - rTmlt * sz_i_1d(ji,jk) ! Melting point of layer jk (C) 480 480 481 481 IF( t_i_1d(ji,jk) >= (ztmelts+rt0) ) THEN !-- Internal melting 482 482 483 zEi = - e_i_1d(ji,jk) * r1_rhoi c! Specific enthalpy of melting ice (J/kg, <0)483 zEi = - e_i_1d(ji,jk) * r1_rhoi ! Specific enthalpy of melting ice (J/kg, <0) 484 484 zdE = 0._wp ! Specific enthalpy difference (J/kg, <0) 485 485 ! set up at 0 since no energy is needed to melt water...(it is already melted) … … 489 489 dh_i_itm (ji) = dh_i_itm(ji) + zdeltah(ji,jk) 490 490 491 zfmdt = - zdeltah(ji,jk) * rhoi c! Mass flux x time step > 0491 zfmdt = - zdeltah(ji,jk) * rhoi ! Mass flux x time step > 0 492 492 493 493 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_rdtice ! Heat flux to the ocean [W.m-2], <0 494 494 ! ice enthalpy zEi is "sent" to the ocean 495 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice! Salt flux495 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice ! Salt flux 496 496 ! using s_i_1d and not sz_i_1d(jk) is ok 497 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice! Mass flux497 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! Mass flux 498 498 499 499 ! update heat content (J.m-2) and layer thickness … … 503 503 ELSE !-- Basal melting 504 504 505 zEi = - e_i_1d(ji,jk) * r1_rhoi c! Specific enthalpy of melting ice (J/kg, <0)505 zEi = - e_i_1d(ji,jk) * r1_rhoi ! Specific enthalpy of melting ice (J/kg, <0) 506 506 zEw = rcp * ztmelts ! Specific enthalpy of meltwater (J/kg, <0) 507 507 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) … … 509 509 zfmdt = - zq_bot(ji) / zdE ! Mass flux x time step (kg/m2, >0) 510 510 511 zdeltah(ji,jk) = - zfmdt * r1_rhoi c! Gross thickness change511 zdeltah(ji,jk) = - zfmdt * r1_rhoi ! Gross thickness change 512 512 513 513 zdeltah(ji,jk) = MIN( 0._wp , MAX( zdeltah(ji,jk), - zh_i(ji,jk) ) ) ! bound thickness change 514 514 515 zq_bot(ji) = MAX( 0._wp , zq_bot(ji) - zdeltah(ji,jk) * rhoi c * zdE )! update available heat. MAX is necessary for roundup errors515 zq_bot(ji) = MAX( 0._wp , zq_bot(ji) - zdeltah(ji,jk) * rhoi * zdE ) ! update available heat. MAX is necessary for roundup errors 516 516 517 517 dh_i_bom(ji) = dh_i_bom(ji) + zdeltah(ji,jk) ! Update basal melt 518 518 519 zfmdt = - zdeltah(ji,jk) * rhoi c! Mass flux x time step > 0519 zfmdt = - zdeltah(ji,jk) * rhoi ! Mass flux x time step > 0 520 520 521 521 zQm = zfmdt * zEw ! Heat exchanged with ocean … … 524 524 hfx_bom_1d(ji) = hfx_bom_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_rdtice ! Heat used in this process [W.m-2], >0 525 525 526 sfx_bom_1d(ji) = sfx_bom_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice! Salt flux526 sfx_bom_1d(ji) = sfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_rdtice ! Salt flux 527 527 ! using s_i_1d and not sz_i_1d(jk) is ok 528 528 529 wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoi c * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice! Mass flux529 wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_rdtice ! Mass flux 530 530 531 531 ! update heat content (J.m-2) and layer thickness … … 558 558 zq_rema(ji) = zq_rema(ji) + zdeltah(ji,1) * e_s_1d(ji,1) ! update available heat (J.m-2) 559 559 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_1d(ji) * e_s_1d(ji,1) * r1_rdtice ! Heat used to melt snow, W.m-2 (>0) 560 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos n * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice! Mass flux560 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice ! Mass flux 561 561 dh_s_mlt(ji) = dh_s_mlt(ji) + zdeltah(ji,1) 562 562 ! … … 572 572 ! When snow load excesses Archimede's limit, snow-ice interface goes down under sea-level, 573 573 ! flooding of seawater transforms snow into ice dh_snowice is positive for the ice 574 z1_rho = 1._wp / ( rhos n+rau0-rhoic)574 z1_rho = 1._wp / ( rhos+rau0-rhoi ) 575 575 DO ji = 1, npti 576 576 ! 577 dh_snowice(ji) = MAX( 0._wp , ( rhos n * h_s_1d(ji) + (rhoic-rau0) * h_i_1d(ji) ) * z1_rho )577 dh_snowice(ji) = MAX( 0._wp , ( rhos * h_s_1d(ji) + (rhoi-rau0) * h_i_1d(ji) ) * z1_rho ) 578 578 579 579 h_i_1d(ji) = h_i_1d(ji) + dh_snowice(ji) … … 581 581 582 582 ! Contribution to energy flux to the ocean [J/m2], >0 (if sst<0) 583 zfmdt = ( rhos n - rhoic) * dh_snowice(ji) ! <0583 zfmdt = ( rhos - rhoi ) * dh_snowice(ji) ! <0 584 584 zEw = rcp * sst_1d(ji) 585 585 zQm = zfmdt * zEw … … 592 592 IF( nn_icesal /= 2 ) THEN 593 593 sfx_bri_1d(ji) = sfx_bri_1d(ji) - sss_1d (ji) * a_i_1d(ji) * zfmdt * r1_rdtice & ! put back sss_m into the ocean 594 & - s_i_1d(ji) * a_i_1d(ji) * dh_snowice(ji) * rhoi c * r1_rdtice! and get rn_icesal from the ocean594 & - s_i_1d(ji) * a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_rdtice ! and get rn_icesal from the ocean 595 595 ENDIF 596 596 597 597 ! Mass flux: All snow is thrown in the ocean, and seawater is taken to replace the volume 598 wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_1d(ji) * dh_snowice(ji) * rhoi c* r1_rdtice599 wfx_snw_sni_1d(ji) = wfx_snw_sni_1d(ji) + a_i_1d(ji) * dh_snowice(ji) * rhos n* r1_rdtice598 wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_rdtice 599 wfx_snw_sni_1d(ji) = wfx_snw_sni_1d(ji) + a_i_1d(ji) * dh_snowice(ji) * rhos * r1_rdtice 600 600 601 601 ! update heat content (J.m-2) and layer thickness … … 619 619 e_s_1d(ji,jk) = rswitch * e_s_1d(ji,jk) 620 620 ! recalculate t_s_1d from e_s_1d 621 t_s_1d(ji,jk) = rt0 + rswitch * ( - e_s_1d(ji,jk) * r1_rhos n * r1_cpic + lfus * r1_cpic)621 t_s_1d(ji,jk) = rt0 + rswitch * ( - e_s_1d(ji,jk) * r1_rhos * r1_rcpi + rLfus * r1_rcpi ) 622 622 END DO 623 623 END DO -
NEMO/trunk/src/ICE/icethd_do.F90
r9604 r9935 140 140 ! Physical constants 141 141 zhicrit = 0.04 ! frazil ice thickness 142 ztwogp = 2. * rau0 / ( grav * 0.3 * ( rau0 - rhoi c ) )! reduced grav142 ztwogp = 2. * rau0 / ( grav * 0.3 * ( rau0 - rhoi ) ) ! reduced grav 143 143 zsqcd = 1.0 / SQRT( 1.3 * zcai ) ! 1/SQRT(airdensity*drag) 144 144 zgamafr = 0.03 … … 263 263 ! We assume that new ice is formed at the seawater freezing point 264 264 DO ji = 1, npti 265 ztmelts = - tmut * zs_newice(ji) ! Melting point (C)266 ze_newice(ji) = rhoi c * ( cpic* ( ztmelts - ( t_bo_1d(ji) - rt0 ) ) &267 & + lfus * ( 1.0 - ztmelts / MIN( t_bo_1d(ji) - rt0, -epsi10 ) ) &268 & - rcp* ztmelts )265 ztmelts = - rTmlt * zs_newice(ji) ! Melting point (C) 266 ze_newice(ji) = rhoi * ( rcpi * ( ztmelts - ( t_bo_1d(ji) - rt0 ) ) & 267 & + rLfus * ( 1.0 - ztmelts / MIN( t_bo_1d(ji) - rt0, -epsi10 ) ) & 268 & - rcp * ztmelts ) 269 269 END DO 270 270 … … 275 275 DO ji = 1, npti 276 276 277 zEi = - ze_newice(ji) * r1_rhoi c! specific enthalpy of forming ice [J/kg]277 zEi = - ze_newice(ji) * r1_rhoi ! specific enthalpy of forming ice [J/kg] 278 278 279 279 zEw = rcp * ( t_bo_1d(ji) - rt0 ) ! specific enthalpy of seawater at t_bo_1d [J/kg] … … 284 284 zfmdt = - qlead_1d(ji) / zdE ! Fm.dt [kg/m2] (<0) 285 285 ! clem: we use qlead instead of zqld (icethd) because we suppose we are at the freezing point 286 zv_newice(ji) = - zfmdt * r1_rhoi c286 zv_newice(ji) = - zfmdt * r1_rhoi 287 287 288 288 zQm = zfmdt * zEw ! heat to the ocean >0 associated with mass flux … … 293 293 hfx_opw_1d(ji) = hfx_opw_1d(ji) - zfmdt * zdE * r1_rdtice 294 294 ! mass flux 295 wfx_opw_1d(ji) = wfx_opw_1d(ji) - zv_newice(ji) * rhoi c* r1_rdtice295 wfx_opw_1d(ji) = wfx_opw_1d(ji) - zv_newice(ji) * rhoi * r1_rdtice 296 296 ! salt flux 297 sfx_opw_1d(ji) = sfx_opw_1d(ji) - zv_newice(ji) * rhoi c* zs_newice(ji) * r1_rdtice297 sfx_opw_1d(ji) = sfx_opw_1d(ji) - zv_newice(ji) * rhoi * zs_newice(ji) * r1_rdtice 298 298 END DO 299 299 -
NEMO/trunk/src/ICE/icethd_pnd.F90
r9750 r9935 133 133 REAL(wp) :: zdv_mlt ! available meltwater for melt ponding 134 134 REAL(wp) :: z1_Tp ! inverse reference temperature 135 REAL(wp) :: z1_rho fw! inverse freshwater density135 REAL(wp) :: z1_rhow ! inverse freshwater density 136 136 REAL(wp) :: z1_zpnd_aspect ! inverse pond aspect ratio 137 137 REAL(wp) :: zfac, zdum … … 139 139 INTEGER :: ji ! loop indices 140 140 !!------------------------------------------------------------------- 141 z1_rho fw = 1._wp / rhofw141 z1_rhow = 1._wp / rhow 142 142 z1_zpnd_aspect = 1._wp / zpnd_aspect 143 143 z1_Tp = 1._wp / zTp … … 157 157 ! 158 158 ! available meltwater for melt ponding [m, >0] and fraction 159 zdv_mlt = -( dh_i_sum(ji)*rhoi c + dh_s_mlt(ji)*rhosn ) * z1_rhofw * a_i_1d(ji)159 zdv_mlt = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow * a_i_1d(ji) 160 160 zfr_mlt = zrmin + ( zrmax - zrmin ) * a_i_1d(ji) ! from CICE doc 161 161 !zfr_mlt = zrmin + zrmax * a_i_1d(ji) ! from Holland paper … … 168 168 ! melt pond mass flux (<0) 169 169 IF( ln_pnd_fwb .AND. zdv_mlt > 0._wp ) THEN 170 zfac = zfr_mlt * zdv_mlt * rho fw * r1_rdtice170 zfac = zfr_mlt * zdv_mlt * rhow * r1_rdtice 171 171 wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac 172 172 ! -
NEMO/trunk/src/ICE/icethd_sal.F90
r9750 r9935 77 77 !--------------------------------------------------------- 78 78 IF( h_i_1d(ji) > 0._wp ) THEN 79 zs_sni = sss_1d(ji) * ( rhoic - rhosn ) * r1_rhoic! Salinity of snow ice79 zs_sni = sss_1d(ji) * ( rhoi - rhos ) * r1_rhoi ! Salinity of snow ice 80 80 zs_i_si = ( zs_sni - s_i_1d(ji) ) * dh_snowice(ji) / h_i_1d(ji) ! snow-ice 81 81 zs_i_bg = ( s_i_new(ji) - s_i_1d(ji) ) * dh_i_bog (ji) / h_i_1d(ji) ! bottom growth … … 98 98 99 99 ! Salt flux 100 sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoi c* a_i_1d(ji) * h_i_1d(ji) * ( zs_i_fl + zs_i_gd ) * r1_rdtice100 sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoi * a_i_1d(ji) * h_i_1d(ji) * ( zs_i_fl + zs_i_gd ) * r1_rdtice 101 101 ENDIF 102 102 END DO -
NEMO/trunk/src/ICE/icethd_zdf_bl99.F90
r9929 r9935 93 93 REAL(wp) :: zdti_bnd = 1.e-4_wp ! maximal authorized error on temperature 94 94 REAL(wp) :: zhs_min = 0.01_wp ! minimum snow thickness for conductivity calculation 95 REAL(wp) :: ztmelt _i! ice melting temperature95 REAL(wp) :: ztmelts ! ice melting temperature 96 96 REAL(wp) :: zdti_max ! current maximal error on temperature 97 97 REAL(wp) :: zcpi ! Ice specific heat … … 217 217 ! 218 218 DO ji = 1, npti 219 ztcond_i(ji,0) = rc dic+ zbeta * sz_i_1d(ji,1) / MIN( -epsi10, t_i_1d(ji,1) - rt0 )220 ztcond_i(ji,nlay_i) = rc dic+ zbeta * sz_i_1d(ji,nlay_i) / MIN( -epsi10, t_bo_1d(ji) - rt0 )219 ztcond_i(ji,0) = rcnd_i + zbeta * sz_i_1d(ji,1) / MIN( -epsi10, t_i_1d(ji,1) - rt0 ) 220 ztcond_i(ji,nlay_i) = rcnd_i + zbeta * sz_i_1d(ji,nlay_i) / MIN( -epsi10, t_bo_1d(ji) - rt0 ) 221 221 END DO 222 222 DO jk = 1, nlay_i-1 223 223 DO ji = 1, npti 224 ztcond_i(ji,jk) = rc dic+ zbeta * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) / &225 & MIN( -epsi10, 0.5_wp * (t_i_1d(ji,jk) + t_i_1d(ji,jk+1)) - rt0 )224 ztcond_i(ji,jk) = rcnd_i + zbeta * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) / & 225 & MIN( -epsi10, 0.5_wp * (t_i_1d(ji,jk) + t_i_1d(ji,jk+1)) - rt0 ) 226 226 END DO 227 227 END DO … … 230 230 ! 231 231 DO ji = 1, npti 232 ztcond_i(ji,0) = rc dic+ 0.09_wp * sz_i_1d(ji,1) / MIN( -epsi10, t_i_1d(ji,1) - rt0 ) &233 & - 0.011_wp * ( t_i_1d(ji,1) - rt0 )234 ztcond_i(ji,nlay_i) = rc dic+ 0.09_wp * sz_i_1d(ji,nlay_i) / MIN( -epsi10, t_bo_1d(ji) - rt0 ) &235 & - 0.011_wp * ( t_bo_1d(ji) - rt0 )232 ztcond_i(ji,0) = rcnd_i + 0.09_wp * sz_i_1d(ji,1) / MIN( -epsi10, t_i_1d(ji,1) - rt0 ) & 233 & - 0.011_wp * ( t_i_1d(ji,1) - rt0 ) 234 ztcond_i(ji,nlay_i) = rcnd_i + 0.09_wp * sz_i_1d(ji,nlay_i) / MIN( -epsi10, t_bo_1d(ji) - rt0 ) & 235 & - 0.011_wp * ( t_bo_1d(ji) - rt0 ) 236 236 END DO 237 237 DO jk = 1, nlay_i-1 238 238 DO ji = 1, npti 239 ztcond_i(ji,jk) = rc dic+ 0.09_wp * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) / &240 & MIN( -epsi10, 0.5_wp * ( t_i_1d (ji,jk) + t_i_1d (ji,jk+1) ) - rt0 ) &241 & - 0.011_wp * ( 0.5_wp * ( t_i_1d (ji,jk) + t_i_1d (ji,jk+1) ) - rt0 )239 ztcond_i(ji,jk) = rcnd_i + 0.09_wp * 0.5_wp * ( sz_i_1d(ji,jk) + sz_i_1d(ji,jk+1) ) / & 240 & MIN( -epsi10, 0.5_wp * ( t_i_1d (ji,jk) + t_i_1d (ji,jk+1) ) - rt0 ) & 241 & - 0.011_wp * ( 0.5_wp * ( t_i_1d (ji,jk) + t_i_1d (ji,jk+1) ) - rt0 ) 242 242 END DO 243 243 END DO … … 299 299 DO jk = 1, nlay_i 300 300 DO ji = 1, npti 301 zcpi = cpic+ zgamma * sz_i_1d(ji,jk) / MAX( ( t_i_1d(ji,jk) - rt0 ) * ( ztiold(ji,jk) - rt0 ), epsi10 )302 zeta_i(ji,jk) = rdt_ice * r1_rhoi c* z1_h_i(ji) / MAX( epsi10, zcpi )301 zcpi = rcpi + zgamma * sz_i_1d(ji,jk) / MAX( ( t_i_1d(ji,jk) - rt0 ) * ( ztiold(ji,jk) - rt0 ), epsi10 ) 302 zeta_i(ji,jk) = rdt_ice * r1_rhoi * z1_h_i(ji) / MAX( epsi10, zcpi ) 303 303 END DO 304 304 END DO … … 306 306 DO jk = 1, nlay_s 307 307 DO ji = 1, npti 308 zeta_s(ji,jk) = rdt_ice * r1_rhos n * r1_cpic* z1_h_s(ji)308 zeta_s(ji,jk) = rdt_ice * r1_rhos * r1_rcpi * z1_h_s(ji) 309 309 END DO 310 310 END DO … … 539 539 DO jk = 1, nlay_i 540 540 DO ji = 1, npti 541 ztmelt _i = -tmut * sz_i_1d(ji,jk) + rt0542 t_i_1d(ji,jk) = MAX( MIN( t_i_1d(ji,jk), ztmelt _i), rt0 - 100._wp )541 ztmelts = -rTmlt * sz_i_1d(ji,jk) + rt0 542 t_i_1d(ji,jk) = MAX( MIN( t_i_1d(ji,jk), ztmelts ), rt0 - 100._wp ) 543 543 zdti_max = MAX( zdti_max, ABS( t_i_1d(ji,jk) - ztib(ji,jk) ) ) 544 544 END DO … … 704 704 DO jk = 1, nlay_i 705 705 DO ji = 1, npti 706 ztmelt _i = -tmut * sz_i_1d(ji,jk) + rt0707 t_i_1d(ji,jk) = MAX( MIN( t_i_1d(ji,jk), ztmelt _i), rt0 - 100._wp )706 ztmelts = -rTmlt * sz_i_1d(ji,jk) + rt0 707 t_i_1d(ji,jk) = MAX( MIN( t_i_1d(ji,jk), ztmelts ), rt0 - 100._wp ) 708 708 zdti_max = MAX( zdti_max, ABS( t_i_1d(ji,jk) - ztib(ji,jk) ) ) 709 709 END DO -
NEMO/trunk/src/ICE/iceupdate.F90
r9922 r9935 171 171 snwice_mass_b(ji,jj) = snwice_mass(ji,jj) ! save mass from the previous ice time step 172 172 ! ! new mass per unit area 173 snwice_mass (ji,jj) = tmask(ji,jj,1) * ( rhos n * vt_s(ji,jj) + rhoic* vt_i(ji,jj) )173 snwice_mass (ji,jj) = tmask(ji,jj,1) * ( rhos * vt_s(ji,jj) + rhoi * vt_i(ji,jj) ) 174 174 ! ! time evolution of snow+ice mass 175 175 snwice_fmass (ji,jj) = ( snwice_mass(ji,jj) - snwice_mass_b(ji,jj) ) * r1_rdtice … … 431 431 ELSE ! start from rest 432 432 IF(lwp) WRITE(numout,*) ' ==>> previous run without snow-ice mass output then set it' 433 snwice_mass (:,:) = tmask(:,:,1) * ( rhos n * vt_s(:,:) + rhoic* vt_i(:,:) )433 snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) ) 434 434 snwice_mass_b(:,:) = snwice_mass(:,:) 435 435 ENDIF 436 436 ELSE !* Start from rest 437 437 IF(lwp) WRITE(numout,*) ' ==>> start from rest: set the snow-ice mass' 438 snwice_mass (:,:) = tmask(:,:,1) * ( rhos n * vt_s(:,:) + rhoic* vt_i(:,:) )438 snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) ) 439 439 snwice_mass_b(:,:) = snwice_mass(:,:) 440 440 ENDIF -
NEMO/trunk/src/ICE/icevar.F90
r9888 r9935 226 226 ! 227 227 ze_i = e_i (ji,jj,jk,jl) * z1_v_i(ji,jj,jl) * zlay_i ! Energy of melting e(S,T) [J.m-3] 228 ztmelts = - sz_i(ji,jj,jk,jl) * tmut ! Ice layer melt temperature [C]228 ztmelts = - sz_i(ji,jj,jk,jl) * rTmlt ! Ice layer melt temperature [C] 229 229 ! Conversion q(S,T) -> T (second order equation) 230 zbbb = ( rcp - cpic ) * ztmelts + ze_i * r1_rhoic - lfus231 zccc = SQRT( MAX( zbbb * zbbb - 4._wp * cpic * lfus * ztmelts , 0._wp) )232 t_i(ji,jj,jk,jl) = MAX( -100._wp , MIN( -( zbbb + zccc ) * 0.5_wp * r1_ cpic, ztmelts ) ) + rt0 ! [K] with bounds: -100 < t_i < ztmelts230 zbbb = ( rcp - rcpi ) * ztmelts + ze_i * r1_rhoi - rLfus 231 zccc = SQRT( MAX( zbbb * zbbb - 4._wp * rcpi * rLfus * ztmelts , 0._wp) ) 232 t_i(ji,jj,jk,jl) = MAX( -100._wp , MIN( -( zbbb + zccc ) * 0.5_wp * r1_rcpi , ztmelts ) ) + rt0 ! [K] with bounds: -100 < t_i < ztmelts 233 233 ! 234 234 ELSE !--- no ice … … 247 247 WHERE( v_s(:,:,:) > epsi20 ) !--- icy area 248 248 t_s(:,:,jk,:) = rt0 + MAX( -100._wp , & 249 & MIN( r1_ cpic * ( -r1_rhosn * ( e_s(:,:,jk,:) / v_s(:,:,:) * zlay_s ) + lfus ) , 0._wp ) )249 & MIN( r1_rcpi * ( -r1_rhos * ( e_s(:,:,jk,:) / v_s(:,:,:) * zlay_s ) + rLfus ) , 0._wp ) ) 250 250 ELSEWHERE !--- no ice 251 251 t_s(:,:,jk,:) = rt0 … … 498 498 DO ji = 1 , jpi 499 499 ! update exchanges with ocean 500 sfx_res(ji,jj) = sfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * sv_i(ji,jj,jl) * rhoi c* r1_rdtice501 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl) * rhoi c* r1_rdtice502 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl) * rhos n* r1_rdtice500 sfx_res(ji,jj) = sfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * sv_i(ji,jj,jl) * rhoi * r1_rdtice 501 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl) * rhoi * r1_rdtice 502 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl) * rhos * r1_rdtice 503 503 ! 504 504 !----------------------------------------------------------------- … … 669 669 ! In case snow load is in excess that would lead to transformation from snow to ice 670 670 ! Then, transfer the snow excess into the ice (different from icethd_dh) 671 zdh = MAX( 0._wp, ( rhos n * zh_s(ji,jl) + ( rhoic- rau0 ) * zh_i(ji,jl) ) * r1_rau0 )671 zdh = MAX( 0._wp, ( rhos * zh_s(ji,jl) + ( rhoi - rau0 ) * zh_i(ji,jl) ) * r1_rau0 ) 672 672 ! recompute h_i, h_s avoiding out of bounds values 673 673 zh_i(ji,jl) = MIN( hi_max(jl), zh_i(ji,jl) + zdh ) 674 zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoi c * r1_rhosn)674 zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoi * r1_rhos ) 675 675 ENDIF 676 676 END DO … … 827 827 DO jk = 1, nlay_i 828 828 WHERE( t_i(:,:,jk,jl) < rt0 - epsi10 ) 829 bv_i(:,:,jl) = bv_i(:,:,jl) - tmut * sz_i(:,:,jk,jl) * r1_nlay_i / ( t_i(:,:,jk,jl) - rt0 )829 bv_i(:,:,jl) = bv_i(:,:,jl) - rTmlt * sz_i(:,:,jk,jl) * r1_nlay_i / ( t_i(:,:,jk,jl) - rt0 ) 830 830 END WHERE 831 831 END DO … … 852 852 DO jk = 1, nlay_i ! Sea ice energy of melting 853 853 DO ji = 1, npti 854 ztmelts = - tmut * sz_i_1d(ji,jk)854 ztmelts = - rTmlt * sz_i_1d(ji,jk) 855 855 t_i_1d(ji,jk) = MIN( t_i_1d(ji,jk), ztmelts + rt0 ) ! Force t_i_1d to be lower than melting point => likely conservation issue 856 856 ! (sometimes zdf scheme produces abnormally high temperatures) 857 e_i_1d(ji,jk) = rhoi c * ( cpic* ( ztmelts - ( t_i_1d(ji,jk) - rt0 ) ) &858 & + lfus * ( 1._wp - ztmelts / ( t_i_1d(ji,jk) - rt0 ) ) &859 & - rcp *ztmelts )857 e_i_1d(ji,jk) = rhoi * ( rcpi * ( ztmelts - ( t_i_1d(ji,jk) - rt0 ) ) & 858 & + rLfus * ( 1._wp - ztmelts / ( t_i_1d(ji,jk) - rt0 ) ) & 859 & - rcp * ztmelts ) 860 860 END DO 861 861 END DO 862 862 DO jk = 1, nlay_s ! Snow energy of melting 863 863 DO ji = 1, npti 864 e_s_1d(ji,jk) = rhos n * ( cpic * ( rt0 - t_s_1d(ji,jk) ) + lfus )864 e_s_1d(ji,jk) = rhos * ( rcpi * ( rt0 - t_s_1d(ji,jk) ) + rLfus ) 865 865 END DO 866 866 END DO -
NEMO/trunk/src/ICE/icewri.F90
r9604 r9935 85 85 ! Standard outputs 86 86 !----------------- 87 zrho1 = ( rau0 - rhoi c ) * r1_rau0; zrho2 = rhosn* r1_rau087 zrho1 = ( rau0 - rhoi ) * r1_rau0; zrho2 = rhos * r1_rau0 88 88 ! masks 89 89 IF( iom_use('icemask' ) ) CALL iom_put( "icemask" , zmsk00 ) ! ice mask 0% … … 92 92 ! 93 93 ! general fields 94 IF( iom_use('icemass' ) ) CALL iom_put( "icemass", rhoi c * vt_i * zmsk00) ! Ice mass per cell area95 IF( iom_use('snwmass' ) ) CALL iom_put( "snwmass", rhos n * vt_s * zmsksn) ! Snow mass per cell area94 IF( iom_use('icemass' ) ) CALL iom_put( "icemass", rhoi * vt_i * zmsk00 ) ! Ice mass per cell area 95 IF( iom_use('snwmass' ) ) CALL iom_put( "snwmass", rhos * vt_s * zmsksn ) ! Snow mass per cell area 96 96 IF( iom_use('icepres' ) ) CALL iom_put( "icepres", zmsk00 ) ! Ice presence (1 or 0) 97 97 IF( iom_use('iceconc' ) ) CALL iom_put( "iceconc", at_i * zmsk00 ) ! ice concentration … … 115 115 ! salt 116 116 IF( iom_use('icesalt' ) ) CALL iom_put( "icesalt", sm_i * zmsk00 ) ! mean ice salinity 117 IF( iom_use('icesalm' ) ) CALL iom_put( "icesalm", SUM( sv_i, DIM = 3 ) * rhoi c* 1.0e-3 * zmsk00 ) ! Mass of salt in sea ice per cell area117 IF( iom_use('icesalm' ) ) CALL iom_put( "icesalm", SUM( sv_i, DIM = 3 ) * rhoi * 1.0e-3 * zmsk00 ) ! Mass of salt in sea ice per cell area 118 118 119 119 ! heat … … 164 164 ! trends 165 165 IF( iom_use('dmithd') ) CALL iom_put( "dmithd", - wfx_bog - wfx_bom - wfx_sum - wfx_sni - wfx_opw - wfx_lam - wfx_res ) ! Sea-ice mass change from thermodynamics 166 IF( iom_use('dmidyn') ) CALL iom_put( "dmidyn", - wfx_dyn + rhoi c * diag_trp_vi) ! Sea-ice mass change from dynamics(kg/m2/s)166 IF( iom_use('dmidyn') ) CALL iom_put( "dmidyn", - wfx_dyn + rhoi * diag_trp_vi ) ! Sea-ice mass change from dynamics(kg/m2/s) 167 167 IF( iom_use('dmiopw') ) CALL iom_put( "dmiopw", - wfx_opw ) ! Sea-ice mass change through growth in open water 168 168 IF( iom_use('dmibog') ) CALL iom_put( "dmibog", - wfx_bog ) ! Sea-ice mass change through basal growth … … 174 174 IF( iom_use('dmisub') ) CALL iom_put( "dmisub", - wfx_ice_sub ) ! Sea-ice mass change through sublimation 175 175 IF( iom_use('dmsspr') ) CALL iom_put( "dmsspr", - wfx_spr ) ! Snow mass change through snow fall 176 IF( iom_use('dmsssi') ) CALL iom_put( "dmsssi", wfx_sni*rhos n*r1_rhoic) ! Snow mass change through snow-to-ice conversion176 IF( iom_use('dmsssi') ) CALL iom_put( "dmsssi", wfx_sni*rhos*r1_rhoi ) ! Snow mass change through snow-to-ice conversion 177 177 IF( iom_use('dmsmel') ) CALL iom_put( "dmsmel", - wfx_snw_sum ) ! Snow mass change through melt 178 IF( iom_use('dmsdyn') ) CALL iom_put( "dmsdyn", - wfx_snw_dyn + rhos n * diag_trp_vs) ! Snow mass change through dynamics(kg/m2/s)178 IF( iom_use('dmsdyn') ) CALL iom_put( "dmsdyn", - wfx_snw_dyn + rhos * diag_trp_vs ) ! Snow mass change through dynamics(kg/m2/s) 179 179 180 180 ! Global ice diagnostics -
NEMO/trunk/src/OCE/BDY/bdyice.F90
r9927 r9935 125 125 126 126 ! Then, a) transfer the snow excess into the ice (different from icethd_dh) 127 zdh = MAX( 0._wp, ( rhos n * h_s(ji,jj,jl) + ( rhoic- rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 )127 zdh = MAX( 0._wp, ( rhos * h_s(ji,jj,jl) + ( rhoi - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 ) 128 128 ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment) 129 !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhos n / rhoic)129 !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhos / rhoi ) 130 130 131 131 ! recompute h_i, h_s 132 132 h_i(ji,jj,jl) = MIN( hi_max(jl), h_i(ji,jj,jl) + zdh ) 133 h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoi c / rhosn)133 h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoi / rhos ) 134 134 135 135 ENDDO … … 198 198 sv_i(ji,jj,jl) = MIN( s_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content 199 199 DO jk = 1, nlay_s 200 e_s(ji,jj,jk,jl) = rhos n * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) ! enthalpy in J/m3200 e_s(ji,jj,jk,jl) = rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus ) ! enthalpy in J/m3 201 201 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s ! enthalpy in J/m2 202 202 END DO 203 203 DO jk = 1, nlay_i 204 ztmelts = - tmut * sz_i(ji,jj,jk,jl)! Melting temperature in C204 ztmelts = - rTmlt * sz_i(ji,jj,jk,jl) ! Melting temperature in C 205 205 t_i(ji,jj,jk,jl) = MIN( t_i(ji,jj,jk,jl), ztmelts + rt0 ) ! Force t_i to be lower than melting point => likely conservation issue 206 206 ! 207 e_i(ji,jj,jk,jl) = rhoi c * ( cpic* ( ztmelts - ( t_i(ji,jj,jk,jl) - rt0 ) ) & ! enthalpy in J/m3208 & + lfus * ( 1._wp - ztmelts / ( t_i(ji,jj,jk,jl) - rt0 ) ) &209 & - rcp* ztmelts )207 e_i(ji,jj,jk,jl) = rhoi * ( rcpi * ( ztmelts - ( t_i(ji,jj,jk,jl) - rt0 ) ) & ! enthalpy in J/m3 208 & + rLfus * ( 1._wp - ztmelts / ( t_i(ji,jj,jk,jl) - rt0 ) ) & 209 & - rcp * ztmelts ) 210 210 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * v_i(ji,jj,jl) * r1_nlay_i ! enthalpy in J/m2 211 211 END DO -
NEMO/trunk/src/OCE/DOM/phycst.F90
r9929 r9935 46 46 REAL(wp), PUBLIC :: r1_rau0_rcp !: = 1. / ( rau0 * rcp ) 47 47 48 !clem: not sure these are needed for cice 49 #if defined key_cice 50 REAL(wp), PUBLIC :: rt0_snow = 273.15_wp !: melting point of snow [Kelvin] 51 REAL(wp), PUBLIC :: rt0_ice = 273.05_wp !: melting point of ice [Kelvin] 52 REAL(wp), PUBLIC :: emic = 0.97_wp !: emissivity of snow or ice 53 REAL(wp), PUBLIC :: xlsn !: = lfus*rhosn (volumetric latent heat fusion of snow) [J/m3] 54 #endif 48 REAL(wp), PUBLIC :: emic = 0.97_wp !: emissivity of snow or ice (not used?) 55 49 56 50 REAL(wp), PUBLIC :: sice = 6.0_wp !: salinity of ice (for pisces) [psu] 57 51 REAL(wp), PUBLIC :: soce = 34.7_wp !: salinity of sea (for pisces and isf) [psu] 58 REAL(wp), PUBLIC :: cevap = 2.5e+6_wp !: latent heat of evaporation (water) 59 REAL(wp), PUBLIC :: srgamma = 0.9_wp !: correction factor for solar radiation (Oberhuber, 1974) 52 REAL(wp), PUBLIC :: rLevap = 2.5e+6_wp !: latent heat of evaporation (water) 60 53 REAL(wp), PUBLIC :: vkarmn = 0.4_wp !: von Karman constant 61 54 REAL(wp), PUBLIC :: stefan = 5.67e-8_wp !: Stefan-Boltzmann constant 62 55 63 REAL(wp), PUBLIC :: rhosn = 330._wp !: volumic mass of snow [kg/m3] 64 REAL(wp), PUBLIC :: rhoic = 917._wp !: volumic mass of sea ice [kg/m3] 65 REAL(wp), PUBLIC :: rhofw = 1000._wp !: volumic mass of freshwater in melt ponds [kg/m3] 66 REAL(wp), PUBLIC :: rcdic = 2.034396_wp !: thermal conductivity of fresh ice [W/m/K] 67 #if defined key_cice 68 REAL(wp), PUBLIC :: rcdsn = 0.31_wp !: thermal conductivity of snow [W/m/K] 69 #endif 70 REAL(wp), PUBLIC :: cpic = 2067.0_wp !: specific heat of fresh ice [J/kg/K] 71 REAL(wp), PUBLIC :: lsub = 2.834e+6_wp !: pure ice latent heat of sublimation [J/kg] 72 REAL(wp), PUBLIC :: lfus = 0.334e+6_wp !: latent heat of fusion of fresh ice [J/kg] 73 REAL(wp), PUBLIC :: tmut = 0.054_wp !: decrease of seawater meltpoint with salinity 56 REAL(wp), PUBLIC :: rhos = 330._wp !: volumic mass of snow [kg/m3] 57 REAL(wp), PUBLIC :: rhoi = 917._wp !: volumic mass of sea ice [kg/m3] 58 REAL(wp), PUBLIC :: rhow = 1000._wp !: volumic mass of freshwater in melt ponds [kg/m3] 59 REAL(wp), PUBLIC :: rcnd_i = 2.034396_wp !: thermal conductivity of fresh ice [W/m/K] 60 REAL(wp), PUBLIC :: rcpi = 2067.0_wp !: specific heat of fresh ice [J/kg/K] 61 REAL(wp), PUBLIC :: rLsub = 2.834e+6_wp !: pure ice latent heat of sublimation [J/kg] 62 REAL(wp), PUBLIC :: rLfus = 0.334e+6_wp !: latent heat of fusion of fresh ice [J/kg] 63 REAL(wp), PUBLIC :: rTmlt = 0.054_wp !: decrease of seawater meltpoint with salinity 74 64 75 REAL(wp), PUBLIC :: r1_rhoi c !: 1 / rhoic76 REAL(wp), PUBLIC :: r1_rhos n !: 1 / rhosn77 REAL(wp), PUBLIC :: r1_ cpic !: 1 / cpic65 REAL(wp), PUBLIC :: r1_rhoi !: 1 / rhoi 66 REAL(wp), PUBLIC :: r1_rhos !: 1 / rhos 67 REAL(wp), PUBLIC :: r1_rcpi !: 1 / rcpi 78 68 !!---------------------------------------------------------------------- 79 69 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 99 89 #endif 100 90 101 #if defined key_cice 102 xlsn = lfus * rhosn 103 #endif 104 105 r1_rhoic = 1._wp / rhoic 106 r1_rhosn = 1._wp / rhosn 107 r1_cpic = 1._wp / cpic 91 r1_rhoi = 1._wp / rhoi 92 r1_rhos = 1._wp / rhos 93 r1_rcpi = 1._wp / rcpi 108 94 109 95 IF(lwp) THEN … … 129 115 WRITE(numout,*) ' reference density and heat capacity now defined in eosbn2.f90' 130 116 WRITE(numout,*) 131 #if defined key_cice 132 WRITE(numout,*) ' thermal conductivity of the snow = ', rcdsn , ' J/s/m/K' 133 #endif 134 WRITE(numout,*) ' thermal conductivity of pure ice = ', rcdic , ' J/s/m/K' 135 WRITE(numout,*) ' fresh ice specific heat = ', cpic , ' J/kg/K' 136 WRITE(numout,*) ' latent heat of fusion of fresh ice / snow = ', lfus , ' J/kg' 137 WRITE(numout,*) ' latent heat of subl. of fresh ice / snow = ', lsub , ' J/kg' 138 WRITE(numout,*) ' density of sea ice = ', rhoic , ' kg/m^3' 139 WRITE(numout,*) ' density of snow = ', rhosn , ' kg/m^3' 140 WRITE(numout,*) ' density of freshwater (in melt ponds) = ', rhofw , ' kg/m^3' 117 WRITE(numout,*) ' thermal conductivity of pure ice = ', rcnd_i , ' J/s/m/K' 118 WRITE(numout,*) ' thermal conductivity of snow is defined in a namelist ' 119 WRITE(numout,*) ' fresh ice specific heat = ', rcpi , ' J/kg/K' 120 WRITE(numout,*) ' latent heat of fusion of fresh ice / snow = ', rLfus , ' J/kg' 121 WRITE(numout,*) ' latent heat of subl. of fresh ice / snow = ', rLsub , ' J/kg' 122 WRITE(numout,*) ' density of sea ice = ', rhoi , ' kg/m^3' 123 WRITE(numout,*) ' density of snow = ', rhos , ' kg/m^3' 124 WRITE(numout,*) ' density of freshwater (in melt ponds) = ', rhow , ' kg/m^3' 141 125 WRITE(numout,*) ' salinity of ice (for pisces) = ', sice , ' psu' 142 126 WRITE(numout,*) ' salinity of sea (for pisces and isf) = ', soce , ' psu' 143 WRITE(numout,*) ' latent heat of evaporation (water) = ', cevap , ' J/m^3' 144 WRITE(numout,*) ' correction factor for solar radiation = ', srgamma 127 WRITE(numout,*) ' latent heat of evaporation (water) = ', rLevap , ' J/m^3' 145 128 WRITE(numout,*) ' von Karman constant = ', vkarmn 146 129 WRITE(numout,*) ' Stefan-Boltzmann constant = ', stefan , ' J/s/m^2/K^4' -
NEMO/trunk/src/OCE/SBC/sbcblk.F90
r9929 r9935 526 526 ! 527 527 qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) & ! Downward Non Solar 528 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * lfus& ! remove latent melting heat for solid precip528 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * rLfus & ! remove latent melting heat for solid precip 529 529 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST 530 530 & + ( sf(jp_prec)%fnow(:,:,1) - sf(jp_snow)%fnow(:,:,1) ) * rn_pfac & ! add liquid precip heat content at Tair 531 531 & * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & 532 532 & + sf(jp_snow)%fnow(:,:,1) * rn_pfac & ! add solid precip heat content at min(Tair,Tsnow) 533 & * ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * cpic533 & * ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * rcpi 534 534 qns(:,:) = qns(:,:) * tmask(:,:,1) 535 535 ! … … 659 659 zrv = pqa(ji,jj) / (1. - pqa(ji,jj)) 660 660 ziRT = 1. / (R_dry*ptak(ji,jj)) ! 1/RT 661 gamma_moist(ji,jj) = grav * ( 1. + cevap*zrv*ziRT ) / ( Cp_dry + cevap*cevap*zrv*reps0*ziRT/ptak(ji,jj) )661 gamma_moist(ji,jj) = grav * ( 1. + rLevap*zrv*ziRT ) / ( Cp_dry + rLevap*rLevap*zrv*reps0*ziRT/ptak(ji,jj) ) 662 662 END DO 663 663 END DO … … 792 792 REAL(wp) :: zst3 ! local variable 793 793 REAL(wp) :: zcoef_dqlw, zcoef_dqla ! - - 794 REAL(wp) :: zztmp, z1_ lsub ! - -794 REAL(wp) :: zztmp, z1_rLsub ! - - 795 795 REAL(wp) :: zfr1, zfr2 ! local variables 796 796 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_st ! inverse of surface temperature … … 868 868 869 869 ! --- evaporation --- ! 870 z1_ lsub = 1._wp /Lsub871 evap_ice (:,:,:) = rn_efac * qla_ice (:,:,:) * z1_ lsub ! sublimation872 devap_ice(:,:,:) = rn_efac * dqla_ice(:,:,:) * z1_ lsub ! d(sublimation)/dT873 zevap (:,:) = rn_efac * ( emp(:,:) + tprecip(:,:) ) ! evaporation over ocean870 z1_rLsub = 1._wp / rLsub 871 evap_ice (:,:,:) = rn_efac * qla_ice (:,:,:) * z1_rLsub ! sublimation 872 devap_ice(:,:,:) = rn_efac * dqla_ice(:,:,:) * z1_rLsub ! d(sublimation)/dT 873 zevap (:,:) = rn_efac * ( emp(:,:) + tprecip(:,:) ) ! evaporation over ocean 874 874 875 875 ! --- evaporation minus precipitation --- ! … … 884 884 & + ( tprecip(:,:) - sprecip(:,:) ) * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & ! liquid precip at Tair 885 885 & + sprecip(:,:) * ( 1._wp - zsnw ) * & ! solid precip at min(Tair,Tsnow) 886 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * cpic * tmask(:,:,1) - lfus )886 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * rcpi * tmask(:,:,1) - rLfus ) 887 887 qemp_ice(:,:) = sprecip(:,:) * zsnw * & ! solid precip (only) 888 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * cpic * tmask(:,:,1) - lfus )888 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * rcpi * tmask(:,:,1) - rLfus ) 889 889 890 890 ! --- total solar and non solar fluxes --- ! … … 894 894 895 895 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) --- ! 896 qprec_ice(:,:) = rhos n * ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * cpic * tmask(:,:,1) - lfus )896 qprec_ice(:,:) = rhos * ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0 ) - rt0 ) * rcpi * tmask(:,:,1) - rLfus ) 897 897 898 898 ! --- heat content of evap over ice in W/m2 (to be used in 1D-thermo) --- 899 899 DO jl = 1, jpl 900 qevap_ice(:,:,jl) = 0._wp ! should be -evap_ice(:,:,jl)*( ( Tice - rt0 ) * cpic* tmask(:,:,1) )900 qevap_ice(:,:,jl) = 0._wp ! should be -evap_ice(:,:,jl)*( ( Tice - rt0 ) * rcpi * tmask(:,:,1) ) 901 901 ! ! But we do not have Tice => consider it at 0degC => evap=0 902 902 END DO … … 971 971 CASE ( 1 , 2 ) 972 972 ! 973 zfac = 1._wp / ( rn_cnd_s + rc dic)973 zfac = 1._wp / ( rn_cnd_s + rcnd_i ) 974 974 zfac2 = EXP(1._wp) * 0.5_wp * zepsilon 975 975 zfac3 = 2._wp / zepsilon … … 978 978 DO jj = 1 , jpj 979 979 DO ji = 1, jpi 980 zhe = ( rn_cnd_s * phi(ji,jj,jl) + rc dic * phs(ji,jj,jl) ) * zfac! Effective thickness980 zhe = ( rn_cnd_s * phi(ji,jj,jl) + rcnd_i * phs(ji,jj,jl) ) * zfac ! Effective thickness 981 981 IF( zhe >= zfac2 ) zgfac(ji,jj,jl) = MIN( 2._wp, 0.5_wp * ( 1._wp + LOG( zhe * zfac3 ) ) ) ! Enhanced conduction factor 982 982 END DO … … 990 990 ! -------------------------------------------------------------! 991 991 ! 992 zfac = rc dic* rn_cnd_s992 zfac = rcnd_i * rn_cnd_s 993 993 ! 994 994 DO jl = 1, jpl … … 997 997 ! 998 998 zkeff_h = zfac * zgfac(ji,jj,jl) / & ! Effective conductivity of the snow-ice system divided by thickness 999 & ( rc dic* phs(ji,jj,jl) + rn_cnd_s * MAX( 0.01, phi(ji,jj,jl) ) )999 & ( rcnd_i * phs(ji,jj,jl) + rn_cnd_s * MAX( 0.01, phi(ji,jj,jl) ) ) 1000 1000 ztsu = ptsu(ji,jj,jl) ! Store current iteration temperature 1001 1001 ztsu0 = ptsu(ji,jj,jl) ! Store initial surface temperature -
NEMO/trunk/src/OCE/SBC/sbccpl.F90
r9921 r9935 1418 1418 zqns(:,:) = zqns(:,:) - zemp(:,:) * sst_m(:,:) * rcp ! remove heat content due to mass flux (assumed to be at SST) 1419 1419 IF( srcv(jpr_snow )%laction ) THEN 1420 zqns(:,:) = zqns(:,:) - frcv(jpr_snow)%z3(:,:,1) * lfus! energy for melting solid precipitation over the free ocean1420 zqns(:,:) = zqns(:,:) - frcv(jpr_snow)%z3(:,:,1) * rLfus ! energy for melting solid precipitation over the free ocean 1421 1421 ENDIF 1422 1422 ENDIF 1423 1423 ! 1424 IF( srcv(jpr_icb)%laction ) zqns(:,:) = zqns(:,:) - frcv(jpr_icb)%z3(:,:,1) * lfus ! remove heat content associated to iceberg melting1424 IF( srcv(jpr_icb)%laction ) zqns(:,:) = zqns(:,:) - frcv(jpr_icb)%z3(:,:,1) * rLfus ! remove heat content associated to iceberg melting 1425 1425 ! 1426 1426 IF( ln_mixcpl ) THEN ; qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:) … … 1811 1811 ! 1812 1812 ! --- calving (removed from qns_tot) --- ! 1813 IF( srcv(jpr_cal)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_cal)%z3(:,:,1) * lfus ! remove latent heat of calving1814 ! we suppose it melts at 0deg, though it should be temp. of surrounding ocean1813 IF( srcv(jpr_cal)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_cal)%z3(:,:,1) * rLfus ! remove latent heat of calving 1814 ! we suppose it melts at 0deg, though it should be temp. of surrounding ocean 1815 1815 ! --- iceberg (removed from qns_tot) --- ! 1816 IF( srcv(jpr_icb)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_icb)%z3(:,:,1) * lfus ! remove latent heat of iceberg melting1816 IF( srcv(jpr_icb)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_icb)%z3(:,:,1) * rLfus ! remove latent heat of iceberg melting 1817 1817 1818 1818 #if defined key_si3 … … 1823 1823 1824 1824 ! Heat content per unit mass of snow (J/kg) 1825 WHERE( SUM( a_i, dim=3 ) > 1.e-10 ) ; zcptsnw(:,:) = cpic* SUM( (tn_ice - rt0) * a_i, dim=3 ) / SUM( a_i, dim=3 )1825 WHERE( SUM( a_i, dim=3 ) > 1.e-10 ) ; zcptsnw(:,:) = rcpi * SUM( (tn_ice - rt0) * a_i, dim=3 ) / SUM( a_i, dim=3 ) 1826 1826 ELSEWHERE ; zcptsnw(:,:) = zcptn(:,:) 1827 1827 ENDWHERE … … 1830 1830 1831 1831 ! --- enthalpy of snow precip over ice in J/m3 (to be used in 1D-thermo) --- ! 1832 zqprec_ice(:,:) = rhos n * ( zcptsnw(:,:) - lfus )1832 zqprec_ice(:,:) = rhos * ( zcptsnw(:,:) - rLfus ) 1833 1833 1834 1834 ! --- heat content of evap over ice in W/m2 (to be used in 1D-thermo) --- ! 1835 1835 DO jl = 1, jpl 1836 zqevap_ice(:,:,jl) = 0._wp ! should be -evap * ( ( Tice - rt0 ) * cpic) but atm. does not take it into account1836 zqevap_ice(:,:,jl) = 0._wp ! should be -evap * ( ( Tice - rt0 ) * rcpi ) but atm. does not take it into account 1837 1837 END DO 1838 1838 … … 1840 1840 zqemp_oce(:,:) = - zevap_oce(:,:) * zcptn (:,:) & ! evap 1841 1841 & + ( ztprecip(:,:) - zsprecip(:,:) ) * zcptrain(:,:) & ! liquid precip 1842 & + zsprecip(:,:) * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - lfus )! solid precip over ocean + snow melting1843 zqemp_ice(:,:) = zsprecip(:,:) * zsnw * ( zcptsnw (:,:) - lfus )! solid precip over ice (qevap_ice=0 since atm. does not take it into account)1842 & + zsprecip(:,:) * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - rLfus ) ! solid precip over ocean + snow melting 1843 zqemp_ice(:,:) = zsprecip(:,:) * zsnw * ( zcptsnw (:,:) - rLfus ) ! solid precip over ice (qevap_ice=0 since atm. does not take it into account) 1844 1844 !! zqemp_ice(:,:) = - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * zcptsnw (:,:) & ! ice evap 1845 !! & + zsprecip(:,:) * zsnw * zqprec_ice(:,:) * r1_rhos n! solid precip over ice1845 !! & + zsprecip(:,:) * zsnw * zqprec_ice(:,:) * r1_rhos ! solid precip over ice 1846 1846 1847 1847 ! --- total non solar flux (including evap/precip) --- ! … … 1874 1874 1875 1875 ! clem: this formulation is certainly wrong... but better than it was... 1876 zqns_tot(:,:) = zqns_tot(:,:) & ! zqns_tot update over free ocean with:1877 & - ( ziceld(:,:) * zsprecip(:,:) * lfus ) & ! remove the latent heat flux of solid precip. melting1878 & - ( zemp_tot(:,:) & ! remove the heat content of mass flux (assumed to be at SST)1876 zqns_tot(:,:) = zqns_tot(:,:) & ! zqns_tot update over free ocean with: 1877 & - ( ziceld(:,:) * zsprecip(:,:) * rLfus ) & ! remove the latent heat flux of solid precip. melting 1878 & - ( zemp_tot(:,:) & ! remove the heat content of mass flux (assumed to be at SST) 1879 1879 & - zemp_ice(:,:) ) * zcptn(:,:) 1880 1880 … … 1892 1892 #endif 1893 1893 ! outputs 1894 IF ( srcv(jpr_cal)%laction ) CALL iom_put('hflx_cal_cea' , - frcv(jpr_cal)%z3(:,:,1) * lfus )! latent heat from calving1895 IF ( srcv(jpr_icb)%laction ) CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * lfus )! latent heat from icebergs melting1894 IF ( srcv(jpr_cal)%laction ) CALL iom_put('hflx_cal_cea' , - frcv(jpr_cal)%z3(:,:,1) * rLfus ) ! latent heat from calving 1895 IF ( srcv(jpr_icb)%laction ) CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * rLfus ) ! latent heat from icebergs melting 1896 1896 IF ( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea' , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) ) ! heat flux from rain (cell average) 1897 1897 IF ( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea' , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) & 1898 1898 & * picefr(:,:) ) * zcptn(:,:) * tmask(:,:,1) ) ! heat flux from evap (cell average) 1899 IF ( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea' , sprecip(:,:) * ( zcptsnw(:,:) - Lfus )) ! heat flux from snow (cell average)1900 IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - Lfus ) &1899 IF ( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea' , sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) ) ! heat flux from snow (cell average) 1900 IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & 1901 1901 & * ( 1._wp - zsnw(:,:) ) ) ! heat flux from snow (over ocean) 1902 IF ( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea', sprecip(:,:) * ( zcptsnw(:,:) - Lfus ) &1902 IF ( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & 1903 1903 & * zsnw(:,:) ) ! heat flux from snow (over ice) 1904 1904 ! note: hflx for runoff and iceshelf are done in sbcrnf and sbcisf resp. -
NEMO/trunk/src/OCE/SBC/sbcice_cice.F90
r9927 r9935 13 13 USE dom_oce ! ocean space and time domain 14 14 USE domvvl 15 USE phycst, only : rcp, rau0, r1_rau0, rhos n, rhoic15 USE phycst, only : rcp, rau0, r1_rau0, rhos, rhoi 16 16 USE in_out_manager ! I/O manager 17 17 USE iom, ONLY : iom_put,iom_use ! I/O manager library !!Joakim edit … … 222 222 CALL cice2nemo(vsno(:,:,:),ztmp1,'T', 1. ) 223 223 CALL cice2nemo(vice(:,:,:),ztmp2,'T', 1. ) 224 snwice_mass (:,:) = ( rhos n * ztmp1(:,:) + rhoic* ztmp2(:,:) )224 snwice_mass (:,:) = ( rhos * ztmp1(:,:) + rhoi * ztmp2(:,:) ) 225 225 snwice_mass_b(:,:) = snwice_mass(:,:) 226 226 … … 328 328 ELSE 329 329 ! emp_ice is set in sbc_cpl_ice_flx as sublimation-snow 330 qla_ice(:,:,1)= - ( emp_ice(:,:)+sprecip(:,:) ) * Lsub330 qla_ice(:,:,1)= - ( emp_ice(:,:)+sprecip(:,:) ) * rLsub 331 331 ! End of temporary code 332 332 DO jj=1,jpj … … 644 644 CALL cice2nemo(vsno(:,:,:),ztmp1,'T', 1. ) 645 645 CALL cice2nemo(vice(:,:,:),ztmp2,'T', 1. ) 646 snwice_mass (:,:) = ( rhos n * ztmp1(:,:) + rhoic* ztmp2(:,:) )646 snwice_mass (:,:) = ( rhos * ztmp1(:,:) + rhoi * ztmp2(:,:) ) 647 647 snwice_mass_b(:,:) = snwice_mass(:,:) 648 648 snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) / dt … … 801 801 tprecip(:,:) = sf(jp_snow)%fnow(:,:,1)+sf(jp_rain)%fnow(:,:,1) 802 802 ! May be better to do this conversion somewhere else 803 qla_ice(:,:,1) = - Lsub*sf(jp_sblm)%fnow(:,:,1)803 qla_ice(:,:,1) = -rLsub*sf(jp_sblm)%fnow(:,:,1) 804 804 topmelt(:,:,1) = sf(jp_top1)%fnow(:,:,1) 805 805 topmelt(:,:,2) = sf(jp_top2)%fnow(:,:,1) -
NEMO/trunk/src/OCE/SBC/sbcisf.F90
r9865 r9935 52 52 LOGICAL, PUBLIC :: l_isfcpl = .false. !: isf recieved from oasis 53 53 54 REAL(wp), PUBLIC, SAVE :: rcpi 54 REAL(wp), PUBLIC, SAVE :: rcpisf = 2000.0_wp !: specific heat of ice shelf [J/kg/K] 55 55 REAL(wp), PUBLIC, SAVE :: rkappa = 1.54e-6_wp !: heat diffusivity through the ice-shelf [m2/s] 56 56 REAL(wp), PUBLIC, SAVE :: rhoisf = 920.0_wp !: volumic mass of ice shelf [kg/m3] 57 57 REAL(wp), PUBLIC, SAVE :: tsurf = -20.0_wp !: air temperature on top of ice shelf [C] 58 REAL(wp), PUBLIC, SAVE :: r lfusisf = 0.334e6_wp !: latent heat of fusion of ice shelf [J/kg]58 REAL(wp), PUBLIC, SAVE :: rLfusisf = 0.334e6_wp !: latent heat of fusion of ice shelf [J/kg] 59 59 60 60 !: Variable used in fldread to read the forcing file (nn_isf == 4 .OR. nn_isf == 3) … … 114 114 ! compute fwf and heat flux 115 115 IF( .NOT.l_isfcpl ) THEN ; CALL sbc_isf_cav (kt) 116 ELSE ; qisf(:,:) = fwfisf(:,:) * r lfusisf ! heat flux116 ELSE ; qisf(:,:) = fwfisf(:,:) * rLfusisf ! heat flux 117 117 ENDIF 118 118 ! … … 127 127 fwfisf(:,:) = - sf_rnfisf(1)%fnow(:,:,1) ! fresh water flux from the isf (fwfisf <0 mean melting) 128 128 ENDIF 129 qisf(:,:) = fwfisf(:,:) * r lfusisf ! heat flux129 qisf(:,:) = fwfisf(:,:) * rLfusisf ! heat flux 130 130 stbl(:,:) = soce 131 131 ! … … 137 137 fwfisf(:,:) = -sf_fwfisf(1)%fnow(:,:,1) ! fwf 138 138 ENDIF 139 qisf(:,:) = fwfisf(:,:) * r lfusisf ! heat flux139 qisf(:,:) = fwfisf(:,:) * rLfusisf ! heat flux 140 140 stbl(:,:) = soce 141 141 ! … … 454 454 & * r1_e1e2t(ji,jj) * tmask(ji,jj,jk) 455 455 456 fwfisf(ji,jj) = qisf(ji,jj) / r lfusisf !fresh water flux kg/(m2s)456 fwfisf(ji,jj) = qisf(ji,jj) / rLfusisf !fresh water flux kg/(m2s) 457 457 fwfisf(ji,jj) = fwfisf(ji,jj) * ( soce / stbl(ji,jj) ) 458 458 !add to salinity trend … … 526 526 DO ji = 1, jpi 527 527 zhtflx(ji,jj) = zgammat(ji,jj)*rcp*rau0*(ttbl(ji,jj)-zfrz(ji,jj)) 528 zfwflx(ji,jj) = - zhtflx(ji,jj)/r lfusisf528 zfwflx(ji,jj) = - zhtflx(ji,jj)/rLfusisf 529 529 END DO 530 530 END DO … … 544 544 ! compute coeficient to solve the 2nd order equation 545 545 zeps1 = rcp*rau0*zgammat(ji,jj) 546 zeps2 = r lfusisf*rau0*zgammas(ji,jj)547 zeps3 = rhoisf*rcpi *rkappa/MAX(risfdep(ji,jj),zeps)546 zeps2 = rLfusisf*rau0*zgammas(ji,jj) 547 zeps3 = rhoisf*rcpisf*rkappa/MAX(risfdep(ji,jj),zeps) 548 548 zeps4 = zlamb2+zlamb3*risfdep(ji,jj) 549 549 zeps6 = zeps4-ttbl(ji,jj) -
NEMO/trunk/src/OCE/SBC/sbcrnf.F90
r9727 r9935 128 128 END WHERE 129 129 WHERE( sf_t_rnf(1)%fnow(:,:,1) == -222._wp ) ! where fwf comes from melting of ice shelves or iceberg 130 rnf_tsc(:,:,jp_tem) = ztfrz(:,:) * rnf(:,:) * r1_rau0 - rnf(:,:) * r lfusisf * r1_rau0_rcp130 rnf_tsc(:,:,jp_tem) = ztfrz(:,:) * rnf(:,:) * r1_rau0 - rnf(:,:) * rLfusisf * r1_rau0_rcp 131 131 END WHERE 132 132 ELSE ! use SST as runoffs temperature
Note: See TracChangeset
for help on using the changeset viewer.
