Changeset 14789 for NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/ICE/icethd.F90
- Timestamp:
- 2021-05-05T13:18:04+02:00 (3 years ago)
- Location:
- NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev _r12970_AGRIF_CMEMSext/AGRIF5 ^/vendors/AGRIF/dev@HEAD ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 ^/vendors/PPR@HEAD ext/PPR 8 9 9 10 # SETTE 10 ^/utils/CI/sette@1 3559sette11 ^/utils/CI/sette@14244 sette
-
- Property svn:externals
-
NEMO/branches/2021/dev_r13747_HPC-11_mcastril_HPDAonline_DiagGPU/src/ICE/icethd.F90
r13643 r14789 69 69 SUBROUTINE ice_thd( kt ) 70 70 !!------------------------------------------------------------------- 71 !! *** ROUTINE ice_thd *** 72 !! 71 !! *** ROUTINE ice_thd *** 72 !! 73 73 !! ** Purpose : This routine manages ice thermodynamics 74 !! 74 !! 75 75 !! ** Action : - computation of oceanic sensible heat flux at the ice base 76 76 !! energy budget in the leads … … 114 114 ztice_cvgerr = 0._wp ; ztice_cvgstp = 0._wp 115 115 ENDIF 116 116 117 117 !---------------------------------------------! 118 118 ! computation of friction velocity at T points … … 136 136 END_2D 137 137 ENDIF 138 CALL lbc_lnk _multi( 'icethd', zfric, 'T', 1.0_wp, zvel, 'T', 1.0_wp )138 CALL lbc_lnk( 'icethd', zfric, 'T', 1.0_wp, zvel, 'T', 1.0_wp ) 139 139 ! 140 140 !--------------------------------------------------------------------! … … 157 157 ! --- Sensible ocean-to-ice heat flux (W/m2) --- ! 158 158 ! (mostly>0 but <0 if supercooling) 159 zfric_u = MAX( SQRT( zfric(ji,jj) ), zfric_umin ) 159 zfric_u = MAX( SQRT( zfric(ji,jj) ), zfric_umin ) 160 160 qsb_ice_bot(ji,jj) = rswitch * rho0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) 161 162 ! upper bound for qsb_ice_bot: the heat retrieved from the ocean must be smaller than the heat necessary to reach 161 162 ! upper bound for qsb_ice_bot: the heat retrieved from the ocean must be smaller than the heat necessary to reach 163 163 ! the freezing point, so that we do not have SST < T_freeze 164 164 ! This implies: qsb_ice_bot(ji,jj) * at_i(ji,jj) * rtdice <= - zqfr_neg … … 166 166 qsb_ice_bot(ji,jj) = rswitch * MIN( qsb_ice_bot(ji,jj), - zqfr_neg * r1_Dt_ice / MAX( at_i(ji,jj), epsi10 ) ) 167 167 168 ! If conditions are always supercooled (such as at the mouth of ice-shelves), then ice grows continuously 169 ! ==> stop ice formation by artificially setting up the turbulent fluxes to 0 when volume > 20m (arbitrary) 170 IF( ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ) > 0._wp .AND. vt_i(ji,jj) >= 20._wp ) THEN 171 zqfr = 0._wp 172 zqfr_pos = 0._wp 173 qsb_ice_bot(ji,jj) = 0._wp 174 ENDIF 175 ! 168 176 ! --- Energy Budget of the leads (qlead, J.m-2) --- ! 169 177 ! qlead is the energy received from the atm. in the leads. … … 202 210 ! 203 211 END_2D 204 212 205 213 ! In case we bypass open-water ice formation 206 214 IF( .NOT. ln_icedO ) qlead(:,:) = 0._wp … … 219 227 npti = 0 ; nptidx(:) = 0 220 228 DO_2D( 1, 1, 1, 1 ) 221 IF ( a_i(ji,jj,jl) > epsi10 ) THEN 229 IF ( a_i(ji,jj,jl) > epsi10 ) THEN 222 230 npti = npti + 1 223 231 nptidx(npti) = (jj - 1) * jpi + ji … … 226 234 227 235 IF( npti > 0 ) THEN ! If there is no ice, do nothing. 228 ! 236 ! 229 237 CALL ice_thd_1d2d( jl, 1 ) ! --- Move to 1D arrays --- ! 230 238 ! ! --- & Change units of e_i, e_s from J/m2 to J/m3 --- ! 231 239 ! 232 s_i_new (1:npti) = 0._wp ; dh_s_tot(1:npti) = 0._wp ! --- some init --- ! (important to have them here) 233 dh_i_sum (1:npti) = 0._wp ; dh_i_bom(1:npti) = 0._wp ; dh_i_itm (1:npti) = 0._wp 240 s_i_new (1:npti) = 0._wp ; dh_s_tot(1:npti) = 0._wp ! --- some init --- ! (important to have them here) 241 dh_i_sum (1:npti) = 0._wp ; dh_i_bom(1:npti) = 0._wp ; dh_i_itm (1:npti) = 0._wp 234 242 dh_i_sub (1:npti) = 0._wp ; dh_i_bog(1:npti) = 0._wp 235 243 dh_snowice(1:npti) = 0._wp ; dh_s_mlt(1:npti) = 0._wp 236 ! 244 ! 237 245 CALL ice_thd_zdf ! --- Ice-Snow temperature --- ! 238 246 ! 239 247 IF( ln_icedH ) THEN ! --- Growing/Melting --- ! 240 CALL ice_thd_dh ! Ice-Snow thickness 241 CALL ice_thd_pnd ! Melt ponds formation 248 CALL ice_thd_dh ! Ice-Snow thickness 242 249 CALL ice_thd_ent( e_i_1d(1:npti,:) ) ! Ice enthalpy remapping 243 250 ENDIF 244 CALL ice_thd_sal( ln_icedS ) ! --- Ice salinity --- ! 251 CALL ice_thd_sal( ln_icedS ) ! --- Ice salinity --- ! 245 252 ! 246 253 CALL ice_thd_temp ! --- Temperature update --- ! … … 259 266 IF( ln_icediachk ) CALL ice_cons_hsm(1, 'icethd', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) 260 267 IF( ln_icediachk ) CALL ice_cons2D (1, 'icethd', diag_v, diag_s, diag_t, diag_fv, diag_fs, diag_ft) 261 ! 268 ! 269 IF ( ln_pnd .AND. ln_icedH ) & 270 & CALL ice_thd_pnd ! --- Melt ponds 271 ! 262 272 IF( jpl > 1 ) CALL ice_itd_rem( kt ) ! --- Transport ice between thickness categories --- ! 263 273 ! … … 266 276 CALL ice_cor( kt , 2 ) ! --- Corrections --- ! 267 277 ! 268 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * r dt_ice ! ice natural aging incrementation278 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rDt_ice ! ice natural aging incrementation 269 279 ! 270 280 ! convergence tests … … 280 290 IF( ln_timing ) CALL timing_stop('icethd') ! timing 281 291 ! 282 END SUBROUTINE ice_thd 283 284 292 END SUBROUTINE ice_thd 293 294 285 295 SUBROUTINE ice_thd_temp 286 296 !!----------------------------------------------------------------------- 287 !! *** ROUTINE ice_thd_temp *** 288 !! 297 !! *** ROUTINE ice_thd_temp *** 298 !! 289 299 !! ** Purpose : Computes sea ice temperature (Kelvin) from enthalpy 290 300 !! … … 292 302 !!------------------------------------------------------------------- 293 303 INTEGER :: ji, jk ! dummy loop indices 294 REAL(wp) :: ztmelts, zbbb, zccc ! local scalar 304 REAL(wp) :: ztmelts, zbbb, zccc ! local scalar 295 305 !!------------------------------------------------------------------- 296 306 ! Recover ice temperature … … 302 312 zccc = SQRT( MAX( zbbb * zbbb - 4._wp * rcpi * rLfus * ztmelts, 0._wp ) ) 303 313 t_i_1d(ji,jk) = rt0 - ( zbbb + zccc ) * 0.5_wp * r1_rcpi 304 314 305 315 ! mask temperature 306 rswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - h_i_1d(ji) ) ) 316 rswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - h_i_1d(ji) ) ) 307 317 t_i_1d(ji,jk) = rswitch * t_i_1d(ji,jk) + ( 1._wp - rswitch ) * rt0 308 END DO 309 END DO 318 END DO 319 END DO 310 320 ! 311 321 END SUBROUTINE ice_thd_temp … … 314 324 SUBROUTINE ice_thd_mono 315 325 !!----------------------------------------------------------------------- 316 !! *** ROUTINE ice_thd_mono *** 317 !! 326 !! *** ROUTINE ice_thd_mono *** 327 !! 318 328 !! ** Purpose : Lateral melting in case virtual_itd 319 329 !! ( dA = A/2h dh ) … … 322 332 REAL(wp) :: zhi_bef ! ice thickness before thermo 323 333 REAL(wp) :: zdh_mel, zda_mel ! net melting 324 REAL(wp) :: zvi, zvs ! ice/snow volumes 334 REAL(wp) :: zvi, zvs ! ice/snow volumes 325 335 !!----------------------------------------------------------------------- 326 336 ! … … 334 344 rswitch = MAX( 0._wp , SIGN( 1._wp , zhi_bef - epsi20 ) ) 335 345 zda_mel = rswitch * a_i_1d(ji) * zdh_mel / ( 2._wp * MAX( zhi_bef, epsi20 ) ) 336 a_i_1d(ji) = MAX( epsi20, a_i_1d(ji) + zda_mel ) 346 a_i_1d(ji) = MAX( epsi20, a_i_1d(ji) + zda_mel ) 337 347 ! adjust thickness 338 h_i_1d(ji) = zvi / a_i_1d(ji) 339 h_s_1d(ji) = zvs / a_i_1d(ji) 348 h_i_1d(ji) = zvi / a_i_1d(ji) 349 h_s_1d(ji) = zvs / a_i_1d(ji) 340 350 ! retrieve total concentration 341 351 at_i_1d(ji) = a_i_1d(ji) … … 348 358 SUBROUTINE ice_thd_1d2d( kl, kn ) 349 359 !!----------------------------------------------------------------------- 350 !! *** ROUTINE ice_thd_1d2d *** 351 !! 360 !! *** ROUTINE ice_thd_1d2d *** 361 !! 352 362 !! ** Purpose : move arrays from 1d to 2d and the reverse 353 363 !!----------------------------------------------------------------------- 354 INTEGER, INTENT(in) :: kl ! index of the ice category 364 INTEGER, INTENT(in) :: kl ! index of the ice category 355 365 INTEGER, INTENT(in) :: kn ! 1= from 2D to 1D ; 2= from 1D to 2D 356 366 ! … … 377 387 CALL tab_2d_1d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl) ) 378 388 END DO 379 CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip (:,:,kl) )380 CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,kl) )381 CALL tab_2d_1d( npti, nptidx(1:npti), h_il_1d (1:npti), h_il (:,:,kl) )382 389 ! 383 390 CALL tab_2d_1d( npti, nptidx(1:npti), qprec_ice_1d (1:npti), qprec_ice ) … … 387 394 CALL tab_2d_1d( npti, nptidx(1:npti), dqns_ice_1d (1:npti), dqns_ice(:,:,kl) ) 388 395 CALL tab_2d_1d( npti, nptidx(1:npti), t_bo_1d (1:npti), t_bo ) 389 CALL tab_2d_1d( npti, nptidx(1:npti), sprecip_1d (1:npti), sprecip ) 396 CALL tab_2d_1d( npti, nptidx(1:npti), sprecip_1d (1:npti), sprecip ) 390 397 CALL tab_2d_1d( npti, nptidx(1:npti), qsb_ice_bot_1d(1:npti), qsb_ice_bot ) 391 398 CALL tab_2d_1d( npti, nptidx(1:npti), fhld_1d (1:npti), fhld ) 392 399 393 400 CALL tab_2d_1d( npti, nptidx(1:npti), qml_ice_1d (1:npti), qml_ice (:,:,kl) ) 394 401 CALL tab_2d_1d( npti, nptidx(1:npti), qcn_ice_1d (1:npti), qcn_ice (:,:,kl) ) … … 409 416 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr ) 410 417 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam ) 411 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd )412 418 ! 413 419 CALL tab_2d_1d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog ) … … 464 470 v_s_1d (1:npti) = h_s_1d (1:npti) * a_i_1d (1:npti) 465 471 sv_i_1d(1:npti) = s_i_1d (1:npti) * v_i_1d (1:npti) 466 v_ip_1d(1:npti) = h_ip_1d(1:npti) * a_ip_1d(1:npti)467 v_il_1d(1:npti) = h_il_1d(1:npti) * a_ip_1d(1:npti)468 472 oa_i_1d(1:npti) = o_i_1d (1:npti) * a_i_1d (1:npti) 469 473 470 474 CALL tab_1d_2d( npti, nptidx(1:npti), at_i_1d(1:npti), at_i ) 471 475 CALL tab_1d_2d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,kl) ) … … 483 487 CALL tab_1d_2d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,kl) ) 484 488 END DO 485 CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip (:,:,kl) )486 CALL tab_1d_2d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,kl) )487 CALL tab_1d_2d( npti, nptidx(1:npti), h_il_1d (1:npti), h_il (:,:,kl) )488 489 ! 489 490 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sni_1d(1:npti), wfx_snw_sni ) … … 501 502 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_spr_1d (1:npti), wfx_spr ) 502 503 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_lam_1d (1:npti), wfx_lam ) 503 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd )504 504 ! 505 505 CALL tab_1d_2d( npti, nptidx(1:npti), sfx_bog_1d (1:npti), sfx_bog ) … … 529 529 CALL tab_1d_2d( npti, nptidx(1:npti), cnd_ice_1d(1:npti), cnd_ice(:,:,kl) ) 530 530 CALL tab_1d_2d( npti, nptidx(1:npti), t1_ice_1d (1:npti), t1_ice (:,:,kl) ) 531 ! SIMIP diagnostics 531 ! Melt ponds 532 CALL tab_1d_2d( npti, nptidx(1:npti), dh_i_sum (1:npti) , dh_i_sum_2d(:,:,kl) ) 533 CALL tab_1d_2d( npti, nptidx(1:npti), dh_s_mlt (1:npti) , dh_s_mlt_2d(:,:,kl) ) 534 ! SIMIP diagnostics 532 535 CALL tab_1d_2d( npti, nptidx(1:npti), t_si_1d (1:npti), t_si (:,:,kl) ) 533 536 CALL tab_1d_2d( npti, nptidx(1:npti), qcn_ice_bot_1d(1:npti), qcn_ice_bot(:,:,kl) ) … … 537 540 CALL tab_1d_2d( npti, nptidx(1:npti), v_s_1d (1:npti), v_s (:,:,kl) ) 538 541 CALL tab_1d_2d( npti, nptidx(1:npti), sv_i_1d(1:npti), sv_i(:,:,kl) ) 539 CALL tab_1d_2d( npti, nptidx(1:npti), v_ip_1d(1:npti), v_ip(:,:,kl) )540 CALL tab_1d_2d( npti, nptidx(1:npti), v_il_1d(1:npti), v_il(:,:,kl) )541 542 CALL tab_1d_2d( npti, nptidx(1:npti), oa_i_1d(1:npti), oa_i(:,:,kl) ) 542 543 ! check convergence of heat diffusion scheme … … 553 554 SUBROUTINE ice_thd_init 554 555 !!------------------------------------------------------------------- 555 !! *** ROUTINE ice_thd_init *** 556 !! 556 !! *** ROUTINE ice_thd_init *** 557 !! 557 558 !! ** Purpose : Physical constants and parameters associated with 558 559 !! ice thermodynamics
Note: See TracChangeset
for help on using the changeset viewer.