Changeset 4506 for branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM
- Timestamp:
- 2014-02-21T15:20:39+01:00 (10 years ago)
- Location:
- branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO
- Files:
-
- 12 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
r4332 r4506 266 266 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fbif !: Heat flux at the ice base 267 267 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_snw !: Variation of snow mass over 1 time step [Kg/m2] 268 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_snw !: Heat content associated with rdm_snw [J/m2]268 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_snw !: Heat per grid cell area released by snow melt [J/m2], >0 if to the ocean 269 269 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdm_ice !: Variation of ice mass over 1 time step [Kg/m2] 270 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_ice !: Heat content associated with rdm_ice [J/m2]270 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rdq_ice !: Heat per grid cell area taken or released by ice growth and melt [J/m2], >0 if to the ocean 271 271 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qldif !: heat balance of the lead (or of the open ocean) 272 272 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qcmif !: Energy needed to bring the ocean to freezing … … 284 284 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_res !: residual salt flux due to correction of ice thickness [PSU/m2/s] 285 285 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fhbri !: heat flux due to brine rejection 286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fheat_mec !: heat flux associated with porous ridged ice formation [???]287 286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fheat_res !: residual heat flux due to correction of ice thickness 288 287 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fmmec !: mass flux due to snow loss during compression [Kg/m2/s] … … 455 454 & fdtcn (jpi,jpj) , qdtcn (jpi,jpj) , fstric (jpi,jpj) , fscmbq (jpi,jpj) , & 456 455 & ffltbif (jpi,jpj) , fsbbq (jpi,jpj) , qfvbq (jpi,jpj) , dmgwi (jpi,jpj) , & 457 & sfx_res (jpi,jpj) , sfx_bri(jpi,jpj) , sfx_mec(jpi,jpj) , fh eat_mec(jpi,jpj) , &458 & f hbri (jpi,jpj) , fmmec (jpi,jpj) , sfx_thd(jpi,jpj) , fhmec (jpi,jpj) ,&456 & sfx_res (jpi,jpj) , sfx_bri(jpi,jpj) , sfx_mec(jpi,jpj) , fhbri (jpi,jpj) , & 457 & fmmec (jpi,jpj) , sfx_thd(jpi,jpj) , fhmec (jpi,jpj) , & 459 458 & fheat_res(jpi,jpj) , STAT=ierr(ii) ) 460 459 -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r4345 r4506 363 363 !-----------------------------------------------------------------------------! 364 364 fmmec(ji,jj) = fmmec(ji,jj) + msnow_mlt(ji,jj) * r1_rdtice ! fresh water source for ocean 365 ! there is a unit problem for e_snow_mlt that we should fix 365 366 fhmec(ji,jj) = fhmec(ji,jj) + esnow_mlt(ji,jj) * r1_rdtice ! heat sink for ocean 366 367 … … 908 909 INTEGER :: ij ! horizontal index, combines i and j loops 909 910 INTEGER :: icells ! number of cells with aicen > puny 910 REAL(wp) :: zindb, zsrdg2 ! local scalar911 REAL(wp) :: zindb, zsrdg2 ! local scalar 911 912 REAL(wp) :: hL, hR, farea, zdummy, zdummy0, ztmelts ! left and right limits of integration 913 REAL(wp) :: zsstK ! SST in Kelvin 912 914 913 915 INTEGER , POINTER, DIMENSION(:) :: indxi, indxj ! compressed indices … … 1098 1100 esrdg(ji,jj) = esnon_init(ji,jj,jl1) * afrac(ji,jj) 1099 1101 srdg1(ji,jj) = smv_i_init(ji,jj,jl1) * afrac(ji,jj) / ( 1._wp + ridge_por ) 1100 srdg2(ji,jj) = smv_i_init(ji,jj,jl1) * afrac(ji,jj) 1102 srdg2(ji,jj) = smv_i_init(ji,jj,jl1) * afrac(ji,jj) !! MV HC 2014 this line seems useless 1101 1103 1102 1104 ! rafting volumes, heat contents ... … … 1120 1122 ! Salinity 1121 1123 !------------- 1122 smsw(ji,jj) = sss_m(ji,jj) * vsw(ji,jj) * rhoic / rau0 ! salt content of seawater frozen in voids 1123 1124 smsw(ji,jj) = vsw(ji,jj) * sss_m(ji,jj) ! salt content of seawater frozen in voids !! MV HC2014 1124 1125 zsrdg2 = srdg1(ji,jj) + smsw(ji,jj) ! salt content of new ridge 1125 1126 … … 1128 1129 ! ! excess of salt is flushed into the ocean 1129 1130 !sfx_mec(ji,jj) = sfx_mec(ji,jj) + ( zsrdg2 - srdg2(ji,jj) ) * rhoic * r1_rdtice 1131 ! MV HC 2014 this previous line seems ok too... Maybe units are not that good 1130 1132 1131 1133 !rdm_ice(ji,jj) = rdm_ice(ji,jj) + vsw(ji,jj) * rhoic ! gurvan: increase in ice volume du to seawater frozen in voids 1134 ! MV HC 2014 this previous line seems ok, i'm not sure at this moment of the sign convention 1132 1135 1133 1136 !------------------------------------ … … 1158 1161 & + rhosn*vsrft(ji,jj)*(1.0-fsnowrft) 1159 1162 1163 ! MV HC 2014 - I think the units are wrong here... shit there is a 10-9 in the snow... 1164 ! we should therefore do something 1160 1165 esnow_mlt(ji,jj) = esnow_mlt(ji,jj) + esrdg(ji,jj)*(1.0-fsnowrdg) & !rafting included 1161 1166 & + esrft(ji,jj)*(1.0-fsnowrft) … … 1187 1192 eirft(ji,jj,jk) = eicen_init(ji,jj,jk,jl1) * afrft(ji,jj) 1188 1193 e_i (ji,jj,jk,jl1) = e_i(ji,jj,jk,jl1) - erdg1(ji,jj,jk) - eirft(ji,jj,jk) 1189 ! sea water heat content 1190 ztmelts = - tmut * sss_m(ji,jj) + rtt 1191 ! heat content per unit volume 1192 zdummy0 = - rcp * ( sst_m(ji,jj) + rt0 - rtt ) * vsw(ji,jj) 1193 1194 ! corrected sea water salinity 1195 zindb = MAX( 0._wp , SIGN( 1._wp , vsw(ji,jj) - epsi20 ) ) 1196 zdummy = zindb * ( srdg1(ji,jj) - srdg2(ji,jj) ) / MAX( ridge_por * vsw(ji,jj), epsi20 ) 1197 1198 ztmelts = - tmut * zdummy + rtt 1199 ersw(ji,jj,jk) = - rcp * ( ztmelts - rtt ) * vsw(ji,jj) 1200 1201 ! heat flux 1202 fheat_mec(ji,jj) = fheat_mec(ji,jj) + ( ersw(ji,jj,jk) - zdummy0 ) * r1_rdtice 1194 1195 1196 ! enthalpy of the trapped seawater (J/m2, >0) 1197 zsstK = sst_m(ji,jj) + rt0 1198 ersw(ji,jj,jk) = - rhoic * vsw(ji,jj) * rcp * ( zsstK - rt0 ) 1199 1200 ! heat flux to the ocean 1201 rdq_ice(ji,jj) = rdq_ice(ji,jj) + ersw(ji,jj,jk) !! MV HC 2014 1203 1202 1204 1203 ! Correct dimensions to avoid big values … … 1206 1205 1207 1206 ! Mutliply by ice volume, and divide by number of layers to get heat content in 10^9 J 1207 ! it is added to sea ice because the sign convention is the opposite of the sign convention for the ocean 1208 !! MV HC 2014 1208 1209 ersw (ji,jj,jk) = ersw(ji,jj,jk) * area(ji,jj) * vsw(ji,jj) / REAL( nlay_i ) 1209 1210 1210 1211 erdg2(ji,jj,jk) = erdg1(ji,jj,jk) + ersw(ji,jj,jk) 1212 1211 1213 END DO ! ij 1212 1214 END DO !jk -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90
r4345 r4506 106 106 INTEGER :: iflt, ial , iadv , ifral, ifrdv ! - - 107 107 REAL(wp) :: zinda, zemp, zemp_snow, zfmm ! local scalars 108 REAL(wp) :: zemp_snw ! - - 108 REAL(wp) :: zemp_snw, zqmass, zcd ! - - 109 REAL(wp) :: zswitch ! - - 109 110 REAL(wp) :: zfcm1 , zfcm2 ! - - 110 111 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace … … 113 114 114 115 IF( lk_cpl ) CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp ) 116 117 SELECT CASE( nn_ice_embd ) ! levitating or embedded sea-ice option 118 CASE( 0 ) ; zswitch = 1 ! (0) standard levitating sea-ice : salt exchange only 119 CASE( 1, 2 ) ; zswitch = 0 ! (1) levitating sea-ice: salt and volume exchange but no pressure effect 120 ! (2) embedded sea-ice : salt and volume fluxes and pressure 121 END SELECT ! 115 122 116 123 !------------------------------------------! … … 172 179 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice & 173 180 & + fhmec(ji,jj) & ! snow melt when ridging 174 & + fheat_mec(ji,jj) & ! ridge formation175 181 & + fheat_res(ji,jj) ! residual heat flux 176 182 ! qcmif Energy needed to bring the ocean surface layer until its freezing (ok) … … 184 190 fsbbq(ji,jj) = ( 1._wp - ( ifvt + iflt ) ) * fscmbq(ji,jj) 185 191 186 ! used to compute the oceanic heat flux at the next time step 192 ! heat flux associated with ice-ocean mass exchange 193 zqmass = ( rdq_snw(ji,jj) & 194 & + rdq_ice(ji,jj) * ( 1.- zswitch) ) * r1_rdtice ! heat flux due to snow & ice heat content 187 195 qsr(ji,jj) = zfcm1 ! solar heat flux 188 qns(ji,jj) = zfcm2 - fdtcn(ji,jj) ! non solar heat flux189 !! fdtcn : turbulent oceanic heat flux196 qns(ji,jj) = zfcm2 - fdtcn(ji,jj) + zqmass ! non solar heat flux 197 ! fdtcn : turbulent oceanic heat flux 190 198 END DO 191 199 END DO -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r4332 r4506 184 184 & + fdtcn(ji,jj) & ! turbulent ice-ocean heat 185 185 & + fsbbq(ji,jj) * ( 1.0 - zinda ) ) & ! residual heat from previous step 186 ! MV HC 2014 check that this is good. 187 ! We should remove the heat content of precip that has fallen on sea ice 186 188 & - pfrld(ji,jj)**betas * sprecip(ji,jj) * lfus ) ! latent heat of sprecip melting 189 ! MV HC 2014 partie heat content manque 187 190 ! 188 191 ! Positive heat budget is used for bottom ablation … … 282 285 CALL tab_2d_1d( nbpb, qldif_1d (1:nbpb), qldif , jpi, jpj, npb(1:nbpb) ) 283 286 CALL tab_2d_1d( nbpb, rdm_ice_1d (1:nbpb), rdm_ice , jpi, jpj, npb(1:nbpb) ) 287 CALL tab_2d_1d( nbpb, rdq_ice_1d (1:nbpb), rdq_ice , jpi, jpj, npb(1:nbpb) ) 284 288 CALL tab_2d_1d( nbpb, rdm_snw_1d (1:nbpb), rdm_snw , jpi, jpj, npb(1:nbpb) ) 289 CALL tab_2d_1d( nbpb, rdq_snw_1d (1:nbpb), rdq_snw , jpi, jpj, npb(1:nbpb) ) 285 290 CALL tab_2d_1d( nbpb, dmgwi_1d (1:nbpb), dmgwi , jpi, jpj, npb(1:nbpb) ) 286 291 CALL tab_2d_1d( nbpb, qlbbq_1d (1:nbpb), zqlbsbq , jpi, jpj, npb(1:nbpb) ) … … 349 354 CALL tab_1d_2d( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb) , jpi, jpj ) 350 355 CALL tab_1d_2d( nbpb, rdm_ice , npb, rdm_ice_1d(1:nbpb) , jpi, jpj ) 356 CALL tab_1d_2d( nbpb, rdq_ice , npb, rdq_ice_1d(1:nbpb) , jpi, jpj ) 351 357 CALL tab_1d_2d( nbpb, rdm_snw , npb, rdm_snw_1d(1:nbpb) , jpi, jpj ) 358 CALL tab_1d_2d( nbpb, rdq_snw , npb, rdq_snw_1d(1:nbpb) , jpi, jpj ) 352 359 CALL tab_1d_2d( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb) , jpi, jpj ) 353 360 CALL tab_1d_2d( nbpb, rdvosif , npb, dvsbq_1d (1:nbpb) , jpi, jpj ) -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r4332 r4506 74 74 INTEGER :: ji , jk ! dummy loop indices 75 75 INTEGER :: ii, ij ! 2D corresponding indices to ji 76 INTEGER :: i_use ! debugging flag 76 77 INTEGER :: isnow ! switch for presence (1) or absence (0) of snow 77 78 INTEGER :: isnowic ! snow ice formation not … … 85 86 REAL(wp) :: zfracs ! fractionation coefficient for bottom salt entrapment 86 87 REAL(wp) :: zcoeff ! dummy argument for snowfall partitioning over ice and leads 87 REAL(wp) :: zs m_snowice! snow-ice salinity88 REAL(wp) :: zs_snic ! snow-ice salinity 88 89 REAL(wp) :: zswi1 ! switch for computation of bottom salinity 89 90 REAL(wp) :: zswi12 ! switch for computation of bottom salinity 90 91 REAL(wp) :: zswi2 ! switch for computation of bottom salinity 91 92 REAL(wp) :: zgrr ! bottom growth rate 92 REAL(wp) :: ztform ! bottom formation temperature 93 ! 93 REAL(wp) :: zt_i_new ! bottom formation temperature 94 95 REAL(wp) :: zQm ! enthalpy exchanged with the ocean (J/m2), >0 towards the ocean 96 REAL(wp) :: zEi ! specific enthalpy of sea ice (J/kg) 97 REAL(wp) :: zEw ! specific enthalpy of exchanged water (J/kg) 98 REAL(wp) :: zdE ! specific enthalpy difference (J/kg) 99 REAL(wp) :: zfmdt ! exchange mass flux x time step (J/m2), >0 towards the ocean 100 REAL(wp) :: zsstK ! SST in Kelvin 101 94 102 REAL(wp), POINTER, DIMENSION(:) :: zh_i ! ice layer thickness 95 103 REAL(wp), POINTER, DIMENSION(:) :: zh_s ! snow layer thickness … … 119 127 120 128 ! mass and salt flux (clem) 121 REAL(wp) :: zdvres, zdvsur, zdvbot 129 REAL(wp) :: zdvres, zdvsur, zdvbot, zswitch_sal 122 130 REAL(wp), POINTER, DIMENSION(:) :: zviold, zvsold ! old ice volume... 123 131 … … 129 137 REAL(wp), POINTER, DIMENSION(:) :: zfbase, zdq_i 130 138 !!------------------------------------------------------------------ 139 140 ! Discriminate between varying salinity (num_sal=2) and prescribed cases (other values) 141 SELECT CASE( num_sal ) ! varying salinity or not 142 CASE( 1, 3, 4 ) ; zswitch_sal = 0 ! prescribed salinity profile 143 CASE( 2 ) ; zswitch_sal = 1 ! varying salinity profile 144 END SELECT 131 145 132 146 CALL wrk_alloc( jpij, zh_i, zh_s, ztfs, zhsold, zqprec, zqfont_su, zqfont_bo, z_f_surf, zhgnew, zfmass_i ) … … 222 236 DO ji = kideb, kiut 223 237 ! tatm_ice is now in K 224 zqprec (ji) = rhosn * ( cpic * ( rtt - tatm_ice_1d(ji) ) + lfus )238 zqprec (ji) = rhosn * ( cpic * ( rtt - MIN( tatm_ice_1d(ji), rt0_snow) ) + lfus ) 225 239 zqfont_su(ji) = z_f_surf(ji) * rdt_ice 226 240 zdeltah (ji,1) = MIN( 0.e0 , - zqfont_su(ji) / MAX( zqprec(ji) , epsi13 ) ) … … 272 286 DO jk = 1, nlay_i 273 287 DO ji = kideb, kiut 274 ! ! melt of layer jk 275 zdeltah (ji,jk) = - zqfont_su(ji) / q_i_b(ji,jk) 276 ! ! recompute heat available 277 zqfont_su(ji ) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk) 278 ! ! melt of layer jk cannot be higher than its thickness 279 zdeltah (ji,jk) = MAX( zdeltah(ji,jk) , - zh_i(ji) ) 280 ! ! update surface melt 281 dh_i_surf(ji ) = dh_i_surf(ji) + zdeltah(ji,jk) 282 ! ! for energy conservation 283 zdq_i (ji ) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) * r1_rdtice 284 ! 285 ! clem 286 sfx_thd_1d(ji) = sfx_thd_1d(ji) - sm_i_b(ji) * a_i_b(ji) & 288 zEi = - q_i_b(ji,jk) / rhoic ! Specific enthalpy of layer k [J/kg, <0] 289 290 ztmelts = - tmut * s_i_b(ji,jk) + rtt ! Melting point of layer k [K] 291 292 zEw = rcp * ( ztmelts - rt0 ) ! Specific enthalpy of resulting meltwater [J/kg, <0] 293 294 zdE = zEi - zEw ! Specific enthalpy difference < 0 295 296 zfmdt = - zqfont_su(ji) / zdE ! Mass flux to the ocean [kg/m2, >0] 297 298 zdeltah(ji,jk) = - zfmdt / rhoic ! Melt of layer jk [m, <0] 299 300 zqfont_su(ji) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * rhoic * ( - zdE ) 301 ! Energy remaining in case of melting of the full layer [J/m2, >0] 302 zdeltah(ji,jk) = MAX( zdeltah(ji,jk) , - zh_i(ji) ) ! Melt of layer jk cannot exceed the layer thickness [m, <0] 303 304 dh_i_surf(ji) = dh_i_surf(ji) + zdeltah(ji,jk) ! Cumulate surface melt 305 306 zfmdt = - rhoic*zdeltah(ji,jk) ! Recompute mass flux [kg/m2, >0] 307 308 zQm = MAX ( zfmdt, 0._wp ) * zEw ! Energy of the melt water sent to the ocean [J/m2, <0] 309 310 ! Contribution to salt flux 311 sfx_thd_1d(ji) = sfx_thd_1d(ji) - sm_i_b(ji) * a_i_b(ji) & 287 312 & * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic / rdt_ice 313 ! Contribution to heat flux 314 rdq_ice_1d(ji) = rdq_ice_1d(ji) + zQm * a_i_b(ji) 315 316 ! Conservation test 317 zdq_i(ji) = zdq_i(ji) + zdeltah(ji,jk) * rhoic * zdE * r1_rdtice ! ? still don't know 288 318 END DO 289 319 END DO … … 364 394 !------------------------------------------------------------------------------! 365 395 ! 366 ! Ice basal growth / melt is given by the ratio of heat budget over basal 367 ! ice heat content. Basal heat budget is given by the difference between 368 ! the inner conductive flux (fc_bo_i), from the open water heat flux 396 !------------------ 397 ! 4.1 Basal growth 398 !------------------ 399 ! Basal growth is driven by heat imbalance at the ice-ocean interface, 400 ! between the inner conductive flux (fc_bo_i), from the open water heat flux 369 401 ! (qlbbqb) and the turbulent ocean flux (fbif). 370 402 ! fc_bo_i is positive downwards. fbif and qlbbq are positive to the ice 371 403 372 !----------------------------------------------------- 373 ! 4.1 Basal growth - (a) salinity not varying in time 374 !----------------------------------------------------- 375 IF( num_sal /= 2 ) THEN ! ice salinity constant in time 404 ! If salinity varies in time, an iterative procedure is required, because 405 ! the involved quantities are inter-dependent. 406 ! Basal growth (dh_i_bott) depends upon new ice specific enthalpy (zEi), 407 ! which depends on forming ice salinity (s_i_new), which depends on dh/dt (dh_i_bott) 408 ! -> need for an iterative procedure, which converges quickly 409 410 IF ( num_sal == 2 ) THEN 411 num_iter_max = 5 412 ELSE 413 num_iter_max = 1 414 ENDIF 415 416 ! Initialize dh_i_bott 417 dh_i_bott(:) = 0.0e0 418 419 ! Iterative procedure 420 DO iter = 1, num_iter_max 376 421 DO ji = kideb, kiut 377 IF( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) < 0._wp ) THEN 378 s_i_new(ji) = sm_i_b(ji) 379 ! Melting point in K 380 ztmelts = - tmut * s_i_new(ji) + rtt 381 ! New ice heat content (Bitz and Lipscomb, 1999) 382 ztform = t_i_b(ji,nlay_i) ! t_bo_b crashes in the 383 ! Baltic 384 q_i_b(ji,nlay_i+1) = rhoic * ( cpic * ( ztmelts - ztform ) & 385 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( ztform - rtt ) ) & 386 & - rcp * ( ztmelts - rtt ) ) 387 ! Basal growth rate = - F*dt / q 388 dh_i_bott(ji) = - rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1) 389 sfx_thd_1d(ji) = sfx_thd_1d(ji) - s_i_new(ji) * a_i_b(ji) * dh_i_bott(ji) * rhoic * r1_rdtice 390 ENDIF 391 END DO 392 ENDIF 393 394 !------------------------------------------------- 395 ! 4.1 Basal growth - (b) salinity varying in time 396 !------------------------------------------------- 397 IF( num_sal == 2 ) THEN 398 ! the growth rate (dh_i_bott) is function of the new ice heat content (q_i_b(nlay_i+1)). 399 ! q_i_b depends on the new ice salinity (snewice). 400 ! snewice depends on dh_i_bott ; it converges quickly, so, no problem 401 ! See Vancoppenolle et al., OM08 for more info on this 402 403 ! Initial value (tested 1D, can be anything between 1 and 20) 404 num_iter_max = 4 405 s_i_new(:) = 4.0 406 407 ! Iterative procedure 408 DO iter = 1, num_iter_max 409 DO ji = kideb, kiut 410 IF( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) < 0.e0 ) THEN 411 ii = MOD( npb(ji) - 1, jpi ) + 1 412 ij = ( npb(ji) - 1 ) / jpi + 1 413 ! Melting point in K 414 ztmelts = - tmut * s_i_new(ji) + rtt 415 ! New ice heat content (Bitz and Lipscomb, 1999) 416 q_i_b(ji,nlay_i+1) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 417 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( t_bo_b(ji) - rtt ) ) & 418 & - rcp * ( ztmelts-rtt ) ) 419 ! Bottom growth rate = - F*dt / q 420 dh_i_bott(ji) = - rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1) 421 ! New ice salinity ( Cox and Weeks, JGR, 1988 ) 422 ! zswi2 (1) if dh_i_bott/rdt .GT. 3.6e-7 423 ! zswi12 (1) if dh_i_bott/rdt .LT. 3.6e-7 and .GT. 2.0e-8 424 ! zswi1 (1) if dh_i_bott/rdt .LT. 2.0e-8 425 zgrr = MIN( 1.0e-3, MAX ( dh_i_bott(ji) * r1_rdtice , epsi13 ) ) 426 zswi2 = MAX( zzero , SIGN( zone , zgrr - 3.6e-7 ) ) 427 zswi12 = MAX( zzero , SIGN( zone , zgrr - 2.0e-8 ) ) * ( 1.0 - zswi2 ) 428 zswi1 = 1. - zswi2 * zswi12 429 zfracs = zswi1 * 0.12 + zswi12 * ( 0.8925 + 0.0568 * LOG( 100.0 * zgrr ) ) & 430 & + zswi2 * 0.26 / ( 0.26 + 0.74 * EXP ( - 724300.0 * zgrr ) ) 431 zfracs = MIN( 0.5 , zfracs ) 432 s_i_new(ji) = zfracs * sss_m(ii,ij) 433 ENDIF ! fc_bo_i 434 END DO ! ji 435 END DO ! iter 436 437 ! Final values 438 DO ji = kideb, kiut 439 IF( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) .LT. 0.0 ) THEN 440 ! New ice salinity must not exceed 20 psu 441 s_i_new(ji) = MIN( s_i_new(ji), s_i_max ) 442 ! Metling point in K 443 ztmelts = - tmut * s_i_new(ji) + rtt 444 ! New ice heat content (Bitz and Lipscomb, 1999) 445 q_i_b(ji,nlay_i+1) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 446 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( t_bo_b(ji) - rtt ) ) & 447 & - rcp * ( ztmelts - rtt ) ) 448 ! Basal growth rate = - F*dt / q 449 dh_i_bott(ji) = - rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / q_i_b(ji,nlay_i+1) 450 ! Salinity update 451 ! entrapment during bottom growth 452 sfx_thd_1d(ji) = sfx_thd_1d(ji) - s_i_new(ji) * a_i_b(ji) * dh_i_bott(ji) * rhoic * r1_rdtice 453 ENDIF ! heat budget 454 END DO 455 ENDIF 422 423 IF( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) < 0.e0 ) THEN 424 425 ! New bottom ice salinity (Cox & Weeks, JGR88 ) 426 !--- zswi1 if dh/dt < 2.0e-8 427 !--- zswi12 if 2.0e-8 < dh/dt < 3.6e-7 428 !--- zswi2 if dh/dt > 3.6e-7 429 zgrr = MIN( 1.0e-3, MAX ( dh_i_bott(ji) * r1_rdtice , epsi13 ) ) 430 zswi2 = MAX( zzero , SIGN( zone , zgrr - 3.6e-7 ) ) 431 zswi12 = MAX( zzero , SIGN( zone , zgrr - 2.0e-8 ) ) * ( 1.0 - zswi2 ) 432 zswi1 = 1. - zswi2 * zswi12 433 zfracs = MIN ( zswi1 * 0.12 + zswi12 * ( 0.8925 + 0.0568 * LOG( 100.0 * zgrr ) ) & 434 & + zswi2 * 0.26 / ( 0.26 + 0.74 * EXP ( - 724300.0 * zgrr ) ) , 0.5 ) 435 436 ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1 437 438 s_i_new(ji) = zswitch_sal * zfracs * sss_m(ii,ij) & ! New ice salinity 439 + ( 1. - zswitch_sal ) * sm_i_b(ji) 440 ! New ice growth 441 ztmelts = - tmut * s_i_new(ji) + rtt ! New ice melting point (K) 442 443 zt_i_new = zswitch_sal * t_bo_b(ji) + ( 1. - zswitch_sal) * t_i_b(ji, nlay_i) 444 445 zEi = cpic * ( zt_i_new - ztmelts ) & ! Specific enthalpy of forming ice (J/kg, <0) 446 & - lfus * ( 1.0 - ( ztmelts - rtt ) / ( zt_i_new - rtt ) ) & 447 & + rcp * ( ztmelts-rtt ) 448 449 zEw = rcp * ( t_bo_b(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0) 450 451 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) 452 453 dh_i_bott(ji) = rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) / ( zdE * rhoic ) 454 455 !!! not sure we still need the next line... useful to keep this in memory ? check limthd_ent... 456 q_i_b(ji,nlay_i+1) = -zEi * rhoic ! New ice energy of melting (J/m3, >0) 457 458 ENDIF ! fc_bo_i 459 END DO ! ji 460 END DO ! iter 461 462 ! Contribution to Energy and Salt Fluxes 463 DO ji = kideb, kiut 464 465 zEw = rcp * ( t_bo_b(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0) 466 467 zfmdt = - rhoic * dh_i_bott(ji) ! Mass flux x time step (kg/m2, < 0) 468 469 ! Contribution to energy flux to the ocean (J/m2) 470 rdq_ice_1d(ji) = rdq_ice_1d(ji) + zEw * a_i_b(ji) * zfmdt 471 472 ! Contribution to salt flux () 473 sfx_thd_1d(ji) = sfx_thd_1d(ji) + s_i_new(ji) * a_i_b(ji) * zfmdt * r1_rdtice 474 475 END DO 456 476 457 477 !---------------- … … 465 485 ! heat convergence at the surface > 0 466 486 IF( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) >= 0._wp ) THEN 467 s_i_new(ji) = s_i_b(ji,nlay_i)487 s_i_new(ji) = s_i_b(ji,nlay_i) ! is this line still useful ? 468 488 zqfont_bo(ji) = rdt_ice * ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) 469 489 zfbase(ji) = zqfont_bo(ji) * r1_rdtice ! heat conservation test … … 476 496 DO ji = kideb, kiut 477 497 IF( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) >= 0._wp ) THEN 478 ztmelts = - tmut * s_i_b(ji,jk) + rtt 479 IF( t_i_b(ji,jk) >= ztmelts ) THEN !!gm : a comment is needed 480 zdeltah (ji,jk) = - zh_i(ji) 481 dh_i_bott (ji ) = dh_i_bott(ji) + zdeltah(ji,jk) 482 zinnermelt(ji ) = 1._wp 483 ELSE ! normal ablation 484 zdeltah (ji,jk) = - zqfont_bo(ji) / q_i_b(ji,jk) 485 zqfont_bo(ji ) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * q_i_b(ji,jk) 486 zdeltah (ji,jk) = MAX(zdeltah(ji,jk), - zh_i(ji) ) 487 dh_i_bott(ji ) = dh_i_bott(ji) + zdeltah(ji,jk) 488 zdq_i (ji ) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) * r1_rdtice 498 499 ztmelts = - tmut * s_i_b(ji,jk) + rtt ! Melting point of layer jk (K) 500 501 IF( t_i_b(ji,jk) >= ztmelts ) THEN !!! Internal melting 502 503 zdeltah (ji,jk) = - zh_i(ji) ! internal melting occurs when the internal temperature is above freezing 504 ! this should normally not happen, but sometimes, heat diffusion leads to this 505 506 dh_i_bott (ji) = dh_i_bott(ji) + zdeltah(ji,jk) 507 508 zinnermelt(ji) = 1._wp 509 510 zQm = 0. ! Not sure which specific enthalpy we should use here (MV HC 2014) 511 ! If that happens, heat is probably not well counted 512 ! Put zero by default, but more bug analysis should be done to investigate this case 513 514 ELSE !!! Basal melting 515 516 zEi = - q_i_b(ji,jk) / rhoic ! Specific enthalpy of melting ice (J/kg, <0) 517 518 zEw = rcp * ( ztmelts - rtt )! Specific enthalpy of meltwater (J/kg, <0) 519 520 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) 521 522 zfmdt = - zqfont_bo(ji) / zdE ! Mass flux x time step (kg/m2, >0) 523 524 zdeltah(ji,jk) = - zfmdt / rhoic ! Gross thickness change 525 526 zqfont_bo(ji) = MAX( 0.0 , - zh_i(ji) - zdeltah(ji,jk) ) * rhoic * ( - zdE ) ! Update heat available 527 528 zdeltah(ji,jk) = MAX( zdeltah(ji,jk), - zh_i(ji) ) ! Update thickness change 529 530 dh_i_bott(ji) = dh_i_bott(ji) + zdeltah(ji,jk) ! Update basal melt 531 532 zfmdt = - zdeltah(ji,jk) * rhoic ! Mass flux x time step 533 534 zQm = MAX ( zfmdt , 0.0 ) * zEw ! Heat exchanged with ocean 535 536 zdq_i(ji) = zdq_i(ji) + zdeltah(ji,jk) * rhoic * zdE * r1_rdtice ! for heat conservation 537 489 538 ENDIF 490 ! clem: contribution to salt flux 539 540 ! Contribution to salt flux 491 541 sfx_thd_1d(ji) = sfx_thd_1d(ji) - sm_i_b(ji) * a_i_b(ji) & 492 542 & * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic * r1_rdtice 543 544 ! Contribution to energy flux to the ocean (J/m2) 545 rdq_ice_1d(ji) = rdq_ice_1d(ji) + zQm * a_i_b(ji) 546 493 547 ENDIF 494 548 END DO ! ji … … 571 625 dvbbq_1d (ji) = a_i_b(ji) * dh_i_bott(ji) 572 626 573 ! remaining heat627 ! remaining heat 574 628 zfdt_final(ji) = ( 1.0 - zihgnew ) * ( zqfont_su(ji) + zqfont_bo(ji) ) 575 629 … … 661 715 DO ji = kideb, kiut 662 716 ht_i_b(ji) = zhgnew(ji) 663 END DO ! ji 717 END DO ! ji 718 664 719 ! 665 720 !------------------------------------------------------------------------------| … … 681 736 !clem rdm_snw_1d(ji) = rdm_snw_1d(ji) + a_i_b(ji) * ( zhnnew - ht_s_b(ji) ) * rhosn 682 737 683 ! Equivalent salt flux (1) Snow-ice formation component 684 ! ----------------------------------------------------- 685 ii = MOD( npb(ji) - 1, jpi ) + 1 686 ij = ( npb(ji) - 1 ) / jpi + 1 687 688 IF( num_sal == 2 ) THEN ; zsm_snowice = sss_m(ii,ij) * ( rhoic - rhosn ) / rhoic 689 ELSE ; zsm_snowice = sm_i_b(ji) 690 ENDIF 738 ! Salinity of snow ice 739 ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1 740 741 zs_snic = zswitch_sal * sss_m(ii,ij) * ( rhoic - rhosn ) / rhoic & 742 + ( 1. - zswitch_sal ) * sm_i_b(ji) 743 691 744 ! entrapment during snow ice formation 692 745 ! clem: new salinity difference stored (to be used in limthd_ent.F90) … … 694 747 i_ice_switch = MAX( 0._wp , SIGN( 1._wp , zhgnew(ji) - epsi10 ) ) 695 748 ! salinity dif due to snow-ice formation 696 dsm_i_si_1d(ji) = ( zs m_snowice- sm_i_b(ji) ) * dh_snowice(ji) / MAX( zhgnew(ji), epsi10 ) * i_ice_switch749 dsm_i_si_1d(ji) = ( zs_snic - sm_i_b(ji) ) * dh_snowice(ji) / MAX( zhgnew(ji), epsi10 ) * i_ice_switch 697 750 ! salinity dif due to bottom growth 698 751 IF ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) < 0._wp ) THEN … … 709 762 710 763 ! diagnostic ( snow ice growth ) 711 ii = MOD( npb(ji) - 1, jpi ) + 1 712 ij = ( npb(ji) - 1 ) / jpi + 1 764 ii = MOD( npb(ji) - 1, jpi ) + 1 !MV HC 2014 useless 765 ij = ( npb(ji) - 1 ) / jpi + 1 ! MV HC 2014 useless 713 766 diag_sni_gr(ii,ij) = diag_sni_gr(ii,ij) + dh_snowice(ji)*a_i_b(ji) * r1_rdtice 714 ! 715 ! salt flux 716 sfx_thd_1d(ji) = sfx_thd_1d(ji) - zsm_snowice * a_i_b(ji) * dh_snowice(ji) * rhoic * r1_rdtice 717 !-------------------------------- 718 ! Update mass fluxes (clem) 719 !-------------------------------- 767 768 ! Contribution to energy flux to the ocean [J/m2], >0 769 zfmdt = ( rhosn - rhoic ) * MAX( dh_snowice(ji), 0.0 ) ! <0 770 zsstK = sst_m(ii,ij) + rt0 771 zEw = rcp * ( zsstK - rt0 ) 772 zQm = zfmdt * zEw 773 rdq_ice_1d(ji) = rdq_ice_1d(ji) + zQm * a_i_b(ji) 774 775 ! Contribution to salt flux 776 sfx_thd_1d(ji) = sfx_thd_1d(ji) + sss_m(ii,ij) * a_i_b(ji) * zfmdt * r1_rdtice 777 778 ! Contribution to mass fluxes 720 779 rdm_ice_1d(ji) = rdm_ice_1d(ji) + ( a_i_b(ji) * ht_i_b(ji) - zviold(ji) ) * rhoic 721 780 rdm_snw_1d(ji) = rdm_snw_1d(ji) + ( a_i_b(ji) * ht_s_b(ji) - zvsold(ji) ) * rhosn -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limthd_ent.F90
r4332 r4506 22 22 USE domain ! 23 23 USE phycst ! physical constants 24 USE sbc_oce ! Surface boundary condition: ocean fields 24 25 USE ice ! LIM variables 25 26 USE par_ice ! LIM parameters … … 86 87 ztmelts , & ! ice melting point 87 88 zqsnic , & ! enthalpy of snow ice layer 89 zsstK , & ! SST in Kelvin 88 90 zhsnow , & ! temporary snow thickness variable 89 91 zswitch , & ! dummy switch argument … … 540 542 DO ji = kideb, kiut 541 543 ! energy of the flooding seawater 542 zqsnic = rau0 * rcp * ( rtt - t_bo_b(ji) ) * dh_snowice(ji) * & 543 (rhoic - rhosn) / rhoic * REAL(snicswi(ji)) ! generally positive 544 ii = MOD( npb(ji) - 1, jpi ) + 1 545 ij = ( npb(ji) - 1 ) / jpi + 1 546 zsstK = sst_m(ii,ij) + rt0 547 zqsnic = ( rhosn - rhoic ) * dh_snowice(ji) * rcp * ( zsstK - rt0 ) * REAL ( snicswi(ji) ) ! MV HC 2014 548 544 549 ! Heat conservation diagnostic 545 550 qt_i_in(ji,jl) = qt_i_in(ji,jl) + zqsnic 546 547 qldif_1d(ji) = qldif_1d(ji) + zqsnic * a_i_b(ji)548 551 549 552 ! enthalpy of the newly formed snow-ice layer -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90
r4332 r4506 81 81 LOGICAL :: iterate_frazil ! iterate frazil ice collection thickness 82 82 CHARACTER (len = 15) :: fieldid 83 ! 83 84 REAL(wp) :: zQm ! enthalpy exchanged with the ocean (J/m2, >0 towards ocean) 85 REAL(wp) :: zEi ! sea ice specific enthalpy (J/kg) 86 REAL(wp) :: zEw ! seawater specific enthalpy (J/kg) 87 REAL(wp) :: zfmdt ! mass flux x time step (kg/m2, >0 towards ocean) 88 84 89 INTEGER , POINTER, DIMENSION(:) :: zcatac ! indexes of categories where new ice grows 85 90 REAL(wp), POINTER, DIMENSION(:) :: zswinew ! switch for new ice or not … … 323 328 CALL tab_2d_1d( nbpac, sfx_thd_1d(1:nbpac) , sfx_thd, jpi, jpj, npac(1:nbpac) ) 324 329 CALL tab_2d_1d( nbpac, rdm_ice_1d(1:nbpac) , rdm_ice, jpi, jpj, npac(1:nbpac) ) 330 CALL tab_2d_1d( nbpac, rdq_ice_1d(1:nbpac) , rdq_ice, jpi, jpj, npac(1:nbpac) ) 325 331 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , jpi, jpj, npac(1:nbpac) ) 326 332 CALL tab_2d_1d( nbpac, zvrel_ac (1:nbpac) , zvrel , jpi, jpj, npac(1:nbpac) ) … … 365 371 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( t_bo_b(ji) - rtt ) ) & 366 372 & - rcp * ( ztmelts - rtt ) ) 373 ! MV HC 2014 comment I dont see why this line below is here... ? 374 ! This implies that ze_newice gets to rhoic*Lfus if it was negative, but this should never happen 367 375 ze_newice(ji) = MAX( ze_newice(ji) , 0._wp ) & 368 376 & + MAX( 0.0 , SIGN( 1.0 , - ze_newice(ji) ) ) * rhoic * lfus … … 386 394 !------------------- 387 395 DO ji = 1, nbpac 388 zv_newice(ji) = - zqbgow(ji) / ze_newice(ji) 396 397 zEi = - ze_newice(ji) / rhoic ! specific enthalpy of forming ice [J/kg] 398 399 zEw = rcp * ( t_bo_b(ji) - rt0 ) ! specific enthalpy of seawater at t_bo_b [J/kg] 400 401 zdE = zEi - zEw ! specific enthalpy difference [J/kg] 402 403 zfmdt = - zqbgow(ji) / zdE ! Fm.dt [kg/m2] (<0) 404 405 zv_newice(ji) = - zfmdt / rhoic 406 407 zQm = zfmdt * zEw ! heat to the ocean >0 associated with mass flux 408 409 ! Contribution to energy flux to the ocean [J/m2], >0 410 rdq_ice_1d(ji) = rdq_ice_1d(ji) + zQm 389 411 390 412 ! A fraction zfrazb of frazil ice is accreted at the ice bottom … … 539 561 END DO 540 562 !!gm ??? why the previous do loop if ocerwriten by the following one ? 563 !! MV HC 2014 just because it'snot the same index and the expression is different 541 564 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 542 565 DO ji = 1, nbpac … … 623 646 CALL tab_1d_2d( nbpac, sfx_thd, npac(1:nbpac), sfx_thd_1d(1:nbpac), jpi, jpj ) 624 647 CALL tab_1d_2d( nbpac, rdm_ice, npac(1:nbpac), rdm_ice_1d(1:nbpac), jpi, jpj ) 648 CALL tab_1d_2d( nbpac, rdq_ice, npac(1:nbpac), rdq_ice_1d(1:nbpac), jpi, jpj ) 625 649 ! 626 650 ENDIF ! nbpac > 0 -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90
r4332 r4506 448 448 CALL prt_ctl(tab2d_1=fmmec , clinfo1= ' lim_update1 : fmmec : ', tab2d_2=fhmec , clinfo2= ' fhmec : ') 449 449 CALL prt_ctl(tab2d_1=sst_m , clinfo1= ' lim_update1 : sst : ', tab2d_2=sss_m , clinfo2= ' sss : ') 450 CALL prt_ctl(tab2d_1=fhbri , clinfo1= ' lim_update1 : fhbri : ' , tab2d_2=fheat_mec , clinfo2= ' fheat_mec : ')450 CALL prt_ctl(tab2d_1=fhbri , clinfo1= ' lim_update1 : fhbri : ') 451 451 452 452 CALL prt_ctl_info(' ') -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90
r4332 r4506 598 598 CALL prt_ctl(tab2d_1=fmmec , clinfo1= ' lim_update2 : fmmec : ', tab2d_2=fhmec , clinfo2= ' fhmec : ') 599 599 CALL prt_ctl(tab2d_1=sst_m , clinfo1= ' lim_update2 : sst : ', tab2d_2=sss_m , clinfo2= ' sss : ') 600 CALL prt_ctl(tab2d_1=fhbri , clinfo1= ' lim_update2 : fhbri : ' , tab2d_2=fheat_mec , clinfo2= ' fheat_mec : ')600 CALL prt_ctl(tab2d_1=fhbri , clinfo1= ' lim_update2 : fhbri : ') 601 601 602 602 CALL prt_ctl_info(' ') -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90
r4045 r4506 66 66 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: fbif_1d !: <==> the 2D fbif 67 67 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: rdm_ice_1d !: <==> the 2D rdm_ice 68 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: rdq_ice_1d !: <==> the 2D rdq_ice 68 69 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: rdm_snw_1d !: <==> the 2D rdm_snw 70 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: rdq_snw_1d !: <==> the 2D rdq_snw 69 71 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: qlbbq_1d !: <==> the 2D qlbsbq 70 72 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: dmgwi_1d !: <==> the 2D dmgwi … … 164 166 ALLOCATE( sprecip_1d (jpij) , frld_1d (jpij) , at_i_b (jpij) , & 165 167 & fbif_1d (jpij) , rdm_ice_1d (jpij) , rdm_snw_1d (jpij) , & 168 & rdq_ice_1d (jpij) , rdq_snw_1d (jpij) , & !MV HC 2014 166 169 & qlbbq_1d (jpij) , dmgwi_1d (jpij) , dvsbq_1d (jpij) , & 167 170 & dvbbq_1d (jpij) , dvlbq_1d (jpij) , dvnbq_1d (jpij) , & -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90
r3625 r4506 54 54 REAL(wp), PUBLIC :: r1_rau0 !: = 1. / rau0 [m3/kg] 55 55 REAL(wp), PUBLIC :: rauw = 1000._wp !: volumic mass of pure water [m3/kg] 56 REAL(wp), PUBLIC :: rcp = 4.e3_wp !: ocean specific heat [J/ Kelvin]57 REAL(wp), PUBLIC :: r1_rcp !: = 1. / rcp [ Kelvin/J]56 REAL(wp), PUBLIC :: rcp = 4.e3_wp !: ocean specific heat [J/kg/K] 57 REAL(wp), PUBLIC :: r1_rcp !: = 1. / rcp [kg.K/J] 58 58 REAL(wp), PUBLIC :: r1_rau0_rcp !: = 1. / ( rau0 * rcp ) 59 59 … … 69 69 #if defined key_lim3 || defined key_cice 70 70 REAL(wp), PUBLIC :: rhoic = 917._wp !: volumic mass of sea ice [kg/m3] 71 REAL(wp), PUBLIC :: rcdic = 2.034396_wp !: thermal conductivity of fresh ice 72 REAL(wp), PUBLIC :: rcdsn = 0.31_wp !: thermal conductivity of snow 73 REAL(wp), PUBLIC :: cpic = 2067.0_wp !: specific heat for ice 71 REAL(wp), PUBLIC :: rcdic = 2.034396_wp !: thermal conductivity of fresh ice [W/m/K] 72 REAL(wp), PUBLIC :: rcdsn = 0.31_wp !: thermal conductivity of snow [W/m/K] 73 REAL(wp), PUBLIC :: cpic = 2067.0_wp !: specific heat for ice [J/kg/K] 74 74 REAL(wp), PUBLIC :: lsub = 2.834e+6_wp !: pure ice latent heat of sublimation [J/kg] 75 75 REAL(wp), PUBLIC :: lfus = 0.334e+6_wp !: latent heat of fusion of fresh ice [J/kg] 76 REAL(wp), PUBLIC :: tmut = 0.054_wp !: decrease of seawater meltpoint with salinity 76 REAL(wp), PUBLIC :: tmut = 0.054_wp !: decrease of seawater meltpoint with salinity [degC/ppt] 77 77 REAL(wp), PUBLIC :: xlsn !: = lfus*rhosn (volumetric latent heat fusion of snow) [J/m3] 78 78 #else -
branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r4358 r4506 136 136 REAL(wp) :: zcoef ! local scalar 137 137 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb_ice_os, zalb_ice_cs ! albedo of the ice under overcast/clear sky 138 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb_ice ! mean albedo of ice (for coupled)138 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb_ice ! mean albedo of ice 139 139 140 140 REAL(wp), POINTER, DIMENSION(:,:) :: zalb_ice_all ! Mean albedo over all categories … … 152 152 IF( nn_timing == 1 ) CALL timing_start('sbc_ice_lim') 153 153 154 CALL wrk_alloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs )154 CALL wrk_alloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice ) 155 155 156 156 #if defined key_coupled 157 IF ( ln_cpl .OR. ln_iceflx_ave .OR. ln_iceflx_linear ) CALL wrk_alloc( jpi,jpj,jpl, zalb_ice)158 157 IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) & 159 158 & CALL wrk_alloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) … … 168 167 ! 169 168 IF( ln_nicep ) THEN ! control print at a given point 170 jiindx = 1 77 ; jjindx = 112169 jiindx = 14 ; jjindx = 25 171 170 IF(lwp) WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx 172 171 ENDIF … … 310 309 sfx (:,:) = 0._wp ; sfx_thd (:,:) = 0._wp 311 310 sfx_bri(:,:) = 0._wp ; sfx_mec (:,:) = 0._wp ; sfx_res (:,:) = 0._wp 312 fhbri (:,:) = 0._wp ; fheat_mec(:,:) = 0._wp ; fheat_res(:,:) = 0._wp 311 fhbri (:,:) = 0._wp ; 312 fheat_res(:,:) = 0._wp 313 313 fhmec (:,:) = 0._wp ; 314 314 fmmec (:,:) = 0._wp … … 327 327 rdm_snw(:,:) = 0._wp ! variation of snow mass per unit area 328 328 rdm_ice(:,:) = 0._wp ! variation of ice mass per unit area 329 rdq_snw(:,:) = 0._wp ! heat flux associated with mass exchange (ice contribution) 330 rdq_ice(:,:) = 0._wp ! heat flux associated with mass exchange (snow contribution) 329 331 hicifp (:,:) = 0._wp ! daily thermodynamic ice production. 330 332 ! … … 421 423 422 424 !!gm remark, the ocean-ice stress is not saved in ice diag call above ..... find a solution!!! 423 ! 424 CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs ) 425 CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice ) 425 426 426 427 #if defined key_coupled 427 IF ( ln_cpl .OR. ln_iceflx_ave .OR. ln_iceflx_linear ) CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice)428 428 429 IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) & 429 430 & CALL wrk_dealloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all) … … 807 808 WRITE(numout,*) ' fhmec : ', fhmec (ji,jj) 808 809 WRITE(numout,*) ' fhbri : ', fhbri (ji,jj) 809 WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj)810 810 WRITE(numout,*) 811 811 WRITE(numout,*) ' sst : ', sst_m(ji,jj)
Note: See TracChangeset
for help on using the changeset viewer.