- Timestamp:
- 2020-07-30T13:20:57+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/tickets_icb_1900/src/OCE/ICB/icbthm.F90
r13265 r13357 48 48 INTEGER, INTENT(in) :: kt ! timestep number, just passed to icb_utl_print_berg 49 49 ! 50 INTEGER :: ii, ij 51 REAL(wp) :: zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn 52 REAL(wp) :: zMv, zMe, zMb, zmelt, zdvo, zdv a, zdM, zSs, zdMe, zdMb, zdMv50 INTEGER :: ii, ij, jk, ikb 51 REAL(wp) :: zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn, zD, zvb, zub, ztb 52 REAL(wp) :: zMv, zMe, zMb, zmelt, zdvo, zdvob, zdva, zdM, zSs, zdMe, zdMb, zdMv 53 53 REAL(wp) :: zMnew, zMnew1, zMnew2, zheat_hcflux, zheat_latent, z1_12 54 54 REAL(wp) :: zMbits, znMbits, zdMbitsE, zdMbitsM, zLbits, zAbits, zMbb 55 REAL(wp) :: zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e2 55 REAL(wp) :: zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e2, zdepw 56 REAL(wp), DIMENSION(jpk) :: ztoce, zuoce, zvoce, ze3t, ztmp 56 57 TYPE(iceberg), POINTER :: this, next 57 58 TYPE(point) , POINTER :: pt … … 83 84 pt => this%current_point 84 85 nknberg = this%number(1) 86 87 CALL icb_utl_interp( pt%xi, pt%yj, & ! position 88 & pssu=pt%ssu, pua=pt%ua, & ! oce/atm velocities 89 & pssv=pt%ssv, pva=pt%va, & ! oce/atm velocities 90 & psst=pt%sst, pcn=pt%cn ) 91 85 92 IF ( nn_sample_rate > 0 .AND. MOD(kt-1,nn_sample_rate) == 0 ) THEN 86 93 CALL icb_utl_interp( pt%xi, pt%yj, pe1=pt%e1, pe2=pt%e2, & 87 & puo=pt%uo, pui=pt%ui, pua=pt%ua, pssh_i=pt%ssh_x, & 88 & pvo=pt%vo, pvi=pt%vi, pva=pt%va, pssh_j=pt%ssh_y, & 89 & psst=pt%sst, pcn=pt%cn, phi=pt%hi, & 90 & plat=pt%lat, plon=pt%lon ) 91 ELSE 92 CALL icb_utl_interp( pt%xi, pt%yj, & ! position 93 & puo=pt%uo, pua=pt%ua, & ! oce/atm velocities 94 & pvo=pt%vo, pva=pt%va, & ! oce/atm velocities 95 & psst=pt%sst, pcn=pt%cn ) ! sst and ice concentration 96 ! preparation Nacho add t3d and uo, vo, to basal 94 & pui=pt%ui, pssh_i=pt%ssh_x, & 95 & pvi=pt%vi, pssh_j=pt%ssh_y, & 96 & phi=pt%hi, & 97 & plat=pt%lat, plon=pt%lon ) 97 98 END IF 98 99 ! … … 101 102 zM = pt%mass 102 103 zT = pt%thickness ! total thickness 103 ! D = (rn_rho_bergs/pp_rho_seawater)*zT! draught (keel depth)104 zD = rho_berg_1_oce * zT ! draught (keel depth) 104 105 ! F = zT - D ! freeboard 105 106 zW = pt%width … … 114 115 115 116 ! Environment 116 zdvo = SQRT( (pt%uvel-pt%uo)**2 + (pt%vvel-pt%vo)**2 ) 117 zdva = SQRT( (pt%ua -pt%uo)**2 + (pt%va -pt%vo)**2 ) 118 zSs = 1.5_wp * SQRT( zdva ) + 0.1_wp * zdva ! Sea state (eqn M.A9) 119 117 ! default sst, ssu and ssv 118 ! ln_M2016: use temp, u and v profile 119 IF ( ln_M2016 ) THEN 120 121 ! load t, u, v and e3 profile at icb position 122 CALL icb_utl_interp( pt%xi, pt%yj, ptoce=ztoce, puoce=zuoce, pvoce=zvoce, pe3t=ze3t ) 123 124 !compute bottom level 125 CALL icb_utl_getkb( pt%kb, ze3t, zD ) 126 127 ikb = MIN(pt%kb,mbkt(ii,ij)) 128 ztb = ztoce(ikb) ! basal temperature 129 zub = zuoce(ikb) 130 zvb = zvoce(ikb) 131 ELSE 132 ztb = pt%sst 133 zub = pt%ssu 134 zvb = pt%ssv 135 END IF 136 137 zdvob = SQRT( (pt%uvel-zub)**2 + (pt%vvel-zvb)**2 ) ! relative basal velocity 138 zdva = SQRT( (pt%ua -pt%ssu)**2 + (pt%va -pt%ssv)**2 ) ! relative wind 139 zSs = 1.5_wp * SQRT( zdva ) + 0.1_wp * zdva ! Sea state (eqn M.A9) 140 ! 120 141 ! Melt rates in m/s (i.e. division by rday) 121 zMv = MAX( 7.62d-3*zSST+1.29d-3*(zSST**2) , 0._wp ) * z1_rday ! Buoyant convection at sides (eqn M.A10) 122 zMb = MAX( 0.58_wp*(zdvo**0.8_wp)*(zSST+4.0_wp)/(zL**0.2_wp) , 0._wp ) * z1_rday ! Basal turbulent melting (eqn M.A7 ) 123 zMe = MAX( z1_12*(zSST+2.)*zSs*(1._wp+COS(rpi*(zIC**3))) , 0._wp ) * z1_rday ! Wave erosion (eqn M.A8 ) 142 IF ( ln_M2016 ) THEN 143 ! Buoyant convection at sides (eqn M.A10) but averaging along all the iceberg draft 144 ztmp(:) = MAX( 7.62d-3*ztoce(:)+1.29d-3*(ztoce(:)**2), 0._wp ) * z1_rday 145 zMv = 0.0; zdepw = 0.0 146 DO jk = 1,ikb-1 147 zMv = zMv + ze3t(jk)*ztmp(jk) 148 zdepw = zdepw + ze3t(jk) 149 END DO 150 zMv = (zMv + (zD - zdepw)*ztmp(ikb)) / zD 151 ELSE 152 zMv = MAX( 7.62d-3*zSST+1.29d-3*(zSST**2), 0._wp ) * z1_rday 153 END IF 154 155 ! Basal turbulent melting (eqn M.A7 ) but using the basal temperature 156 zMb = MAX( 0.58_wp*(zdvob**0.8_wp)*(ztb+4.0_wp)/(zL**0.2_wp) , 0._wp ) * z1_rday 157 158 ! Wave erosion (eqn M.A8 ) 159 zMe = MAX( z1_12*(zSST+2.)*zSs*(1._wp+COS(rpi*(zIC**3))) , 0._wp ) * z1_rday 124 160 125 161 IF( ln_operator_splitting ) THEN ! Operator split update of volume/mass … … 209 245 210 246 ! Rolling 211 zDn = ( rn_rho_bergs / pp_rho_seawater )* zTn ! draught (keel depth)247 zDn = rho_berg_1_oce * zTn ! draught (keel depth) 212 248 IF( zDn > 0._wp .AND. MAX(zWn,zLn) < SQRT( 0.92*(zDn**2) + 58.32*zDn ) ) THEN 213 249 zT = zTn
Note: See TracChangeset
for help on using the changeset viewer.