- Timestamp:
- 2017-07-15T17:27:14+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r8341 r8342 42 42 CONTAINS 43 43 44 SUBROUTINE lim_thd_dh ( kideb, kiut )44 SUBROUTINE lim_thd_dh 45 45 !!------------------------------------------------------------------ 46 46 !! *** ROUTINE lim_thd_dh *** … … 66 66 !! Vancoppenolle et al.,2009, Ocean Modelling 67 67 !!------------------------------------------------------------------ 68 INTEGER , INTENT(in) :: kideb, kiut ! Start/End point on which the the computation is applied69 !!70 68 INTEGER :: ji , jk ! dummy loop indices 71 69 INTEGER :: iter … … 130 128 131 129 ! Initialize enthalpy at nlay_i+1 132 DO ji = kideb, kiut130 DO ji = 1, nidx 133 131 e_i_1d(ji,nlay_i+1) = 0._wp 134 132 END DO … … 138 136 eh_i_old(:,0:nlay_i+1) = 0._wp 139 137 DO jk = 1, nlay_i 140 DO ji = kideb, kiut138 DO ji = 1, nidx 141 139 h_i_old (ji,jk) = ht_i_1d(ji) * r1_nlay_i 142 140 eh_i_old(ji,jk) = e_i_1d(ji,jk) * h_i_old(ji,jk) … … 148 146 !------------------------------------------------------------------------------! 149 147 ! 150 DO ji = kideb, kiut148 DO ji = 1, nidx 151 149 zdum = qns_ice_1d(ji) + ( 1._wp - i0(ji) ) * qsr_ice_1d(ji) - fc_su(ji) 152 150 zf_tt(ji) = fc_bo_i(ji) + fhtur_1d(ji) + fhld_1d(ji) … … 161 159 ! (should not happen but sometimes it does) 162 160 !------------------------------------------------------------------------------! 163 DO ji = kideb, kiut161 DO ji = 1, nidx 164 162 IF( t_s_1d(ji,1) > rt0 ) THEN !!! Internal melting 165 163 ! Contribution to heat flux to the ocean [W.m-2], < 0 … … 179 177 ! 180 178 DO jk = 1, nlay_i 181 DO ji = kideb, kiut179 DO ji = 1, nidx 182 180 zh_i(ji,jk) = ht_i_1d(ji) * r1_nlay_i 183 181 zeh_i(ji) = zeh_i(ji) + e_i_1d(ji,jk) * zh_i(ji,jk) … … 203 201 ! Martin Vancoppenolle, December 2006 204 202 205 CALL lim_thd_snwblow( 1. - at_i_1d( kideb:kiut), zsnw(kideb:kiut) ) ! snow distribution over ice after wind blowing203 CALL lim_thd_snwblow( 1. - at_i_1d(1:nidx), zsnw(1:nidx) ) ! snow distribution over ice after wind blowing 206 204 207 205 zdeltah(:,:) = 0._wp 208 DO ji = kideb, kiut206 DO ji = 1, nidx 209 207 !----------- 210 208 ! Snow fall … … 242 240 zdeltah(:,:) = 0._wp 243 241 DO jk = 1, nlay_s 244 DO ji = kideb, kiut242 DO ji = 1, nidx 245 243 ! thickness change 246 244 rswitch = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) ) … … 265 263 ! clem comment: not counted in mass/heat exchange in limsbc since this is an exchange with atm. (not ocean) 266 264 zdeltah(:,:) = 0._wp 267 DO ji = kideb, kiut265 DO ji = 1, nidx 268 266 zdh_s_sub(ji) = MAX( - ht_s_1d(ji) , - evap_ice_1d(ji) * r1_rhosn * rdt_ice ) 269 267 ! remaining evap in kg.m-2 (used for ice melting later on) … … 284 282 285 283 ! --- Update snow diags --- ! 286 DO ji = kideb, kiut284 DO ji = 1, nidx 287 285 dh_s_tot(ji) = zdh_s_mel(ji) + zdh_s_pre(ji) + zdh_s_sub(ji) 288 286 END DO … … 293 291 ! new temp and enthalpy of the snow (remaining snow precip + remaining pre-existing snow) 294 292 DO jk = 1, nlay_s 295 DO ji = kideb,kiut293 DO ji = 1,nidx 296 294 rswitch = MAX( 0._wp , SIGN( 1._wp, ht_s_1d(ji) - epsi20 ) ) 297 295 e_s_1d(ji,jk) = rswitch / MAX( ht_s_1d(ji), epsi20 ) * & … … 306 304 zdeltah(:,:) = 0._wp ! important 307 305 DO jk = 1, nlay_i 308 DO ji = kideb, kiut306 DO ji = 1, nidx 309 307 ztmelts = - tmut * s_i_1d(ji,jk) + rt0 ! Melting point of layer k [K] 310 308 … … 394 392 END DO 395 393 ! update ice thickness 396 DO ji = kideb, kiut394 DO ji = 1, nidx 397 395 ht_i_1d(ji) = MAX( 0._wp , ht_i_1d(ji) + dh_i_surf(ji) + dh_i_sub(ji) ) 398 396 END DO 399 397 400 398 ! remaining "potential" evap is sent to ocean 401 DO ji = kideb, kiut399 DO ji = 1, nidx 402 400 wfx_err_sub_1d(ji) = wfx_err_sub_1d(ji) - zevap_rema(ji) * a_i_1d(ji) * r1_rdtice ! <=0 (net evap for the ocean in kg.m-2.s-1) 403 401 END DO … … 426 424 427 425 ! Iterative procedure 428 DO ji = kideb, kiut426 DO ji = 1, nidx 429 427 IF( zf_tt(ji) < 0._wp ) THEN 430 428 DO iter = 1, num_iter_max … … 501 499 zdeltah(:,:) = 0._wp ! important 502 500 DO jk = nlay_i, 1, -1 503 DO ji = kideb, kiut501 DO ji = 1, nidx 504 502 IF( zf_tt(ji) > 0._wp .AND. jk > icount(ji,jk) ) THEN ! do not calculate where layer has already disappeared by surface melting 505 503 … … 575 573 ! Update temperature, energy 576 574 !------------------------------------------- 577 DO ji = kideb, kiut575 DO ji = 1, nidx 578 576 ht_i_1d(ji) = MAX( 0._wp , ht_i_1d(ji) + dh_i_bott(ji) ) 579 577 END DO … … 585 583 !------------------------------------------- 586 584 zdeltah(:,:) = 0._wp ! important 587 DO ji = kideb, kiut585 DO ji = 1, nidx 588 586 zq_rema(ji) = zq_su(ji) + zq_bo(ji) 589 587 rswitch = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) ) ! =1 if snow … … 612 610 ! When snow load excesses Archimede's limit, snow-ice interface goes down under sea-level, 613 611 ! flooding of seawater transforms snow into ice dh_snowice is positive for the ice 614 DO ji = kideb, kiut612 DO ji = 1, nidx 615 613 ! 616 614 dh_snowice(ji) = MAX( 0._wp , ( rhosn * ht_s_1d(ji) + (rhoic-rau0) * ht_i_1d(ji) ) / ( rhosn+rau0-rhoic ) ) … … 651 649 ! Update temperature, energy 652 650 !------------------------------------------- 653 DO ji = kideb, kiut651 DO ji = 1, nidx 654 652 rswitch = 1.0 - MAX( 0._wp , SIGN( 1._wp , - ht_i_1d(ji) ) ) 655 653 t_su_1d(ji) = rswitch * t_su_1d(ji) + ( 1.0 - rswitch ) * rt0 … … 657 655 658 656 DO jk = 1, nlay_s 659 DO ji = kideb,kiut657 DO ji = 1,nidx 660 658 ! mask enthalpy 661 659 rswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp, - ht_s_1d(ji) ) )
Note: See TracChangeset
for help on using the changeset viewer.
