- Timestamp:
- 2018-11-07T18:25:49+01:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2018/dev_r9866_HPC_03_globcom/src/OCE/ICB/icbthm.F90
r9598 r10288 33 33 !! NEMO/OCE 4.0 , NEMO Consortium (2018) 34 34 !! $Id$ 35 !! Software governed by the CeCILL licen ce (./LICENSE)35 !! Software governed by the CeCILL license (see ./LICENSE) 36 36 !!---------------------------------------------------------------------- 37 37 CONTAINS … … 50 50 REAL(wp) :: zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn 51 51 REAL(wp) :: zMv, zMe, zMb, zmelt, zdvo, zdva, zdM, zSs, zdMe, zdMb, zdMv 52 REAL(wp) :: zMnew, zMnew1, zMnew2, zheat 52 REAL(wp) :: zMnew, zMnew1, zMnew2, zheat_hcflux, zheat_latent, z1_12 53 53 REAL(wp) :: zMbits, znMbits, zdMbitsE, zdMbitsM, zLbits, zAbits, zMbb 54 54 REAL(wp) :: zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e2 … … 58 58 ! 59 59 z1_rday = 1._wp / rday 60 z1_12 = 1._wp / 12._wp 61 zdt = berg_dt 62 z1_dt = 1._wp / zdt 60 63 ! 61 64 ! we're either going to ignore berg fresh water melt flux and associated heat … … 65 68 ! 66 69 berg_grid%floating_melt(:,:) = 0._wp 70 ! calving_hflx re-used here as temporary workspace for the heat flux associated with melting 67 71 berg_grid%calving_hflx(:,:) = 0._wp 68 72 ! … … 91 95 ij = mj1( ij ) 92 96 zVol = zT * zW * zL 93 zdt = berg_dt ; z1_dt = 1._wp / zdt94 97 95 98 ! Environment 96 99 zdvo = SQRT( (pt%uvel-pt%uo)**2 + (pt%vvel-pt%vo)**2 ) 97 100 zdva = SQRT( (pt%ua -pt%uo)**2 + (pt%va -pt%vo)**2 ) 98 zSs = 1.5 * SQRT( zdva ) + 0.1* zdva ! Sea state (eqn M.A9)101 zSs = 1.5_wp * SQRT( zdva ) + 0.1_wp * zdva ! Sea state (eqn M.A9) 99 102 100 103 ! Melt rates in m/s (i.e. division by rday) 101 zMv = MAX( 7.62 e-3*zSST+1.29e-3*(zSST**2), 0._wp ) * z1_rday ! Buoyant convection at sides (eqn M.A10)102 zMb = MAX( 0.58 *(zdvo**0.8)*(zSST+4.0)/(zL**0.2), 0._wp ) * z1_rday ! Basal turbulent melting (eqn M.A7 )103 zMe = MAX( 1./12.*(zSST+2.)*zSs*(1+COS(rpi*(zIC**3))), 0._wp ) * z1_rday ! Wave erosion (eqn M.A8 )104 zMv = MAX( 7.62d-3*zSST+1.29d-3*(zSST**2) , 0._wp ) * z1_rday ! Buoyant convection at sides (eqn M.A10) 105 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 ) 106 zMe = MAX( z1_12*(zSST+2.)*zSs*(1._wp+COS(rpi*(zIC**3))) , 0._wp ) * z1_rday ! Wave erosion (eqn M.A8 ) 104 107 105 108 IF( ln_operator_splitting ) THEN ! Operator split update of volume/mass 106 109 zTn = MAX( zT - zMb*zdt , 0._wp ) ! new total thickness (m) 107 110 znVol = zTn * zW * zL ! new volume (m^3) 108 zMnew1 = ( znVol/zVol) * zM! new mass (kg)111 zMnew1 = ( znVol / zVol ) * zM ! new mass (kg) 109 112 zdMb = zM - zMnew1 ! mass lost to basal melting (>0) (kg) 110 113 ! … … 112 115 zWn = MAX( zW - zMv*zdt , 0._wp ) ! new width (m) 113 116 znVol = zTn * zWn * zLn ! new volume (m^3) 114 zMnew2 = ( znVol/zVol) * zM! new mass (kg)117 zMnew2 = ( znVol / zVol ) * zM ! new mass (kg) 115 118 zdMv = zMnew1 - zMnew2 ! mass lost to buoyant convection (>0) (kg) 116 119 ! … … 142 145 zLbits = MIN( zL, zW, zT, 40._wp ) ! assume bergy bits are smallest dimension or 40 meters 143 146 zAbits = ( zMbits / rn_rho_bergs ) / zLbits ! Effective bottom area (assuming T=Lbits) 144 zMbb = MAX( 0.58*(zdvo**0.8)*(zSST+2.0)/(zLbits**0.2), 0.) * z1_rday ! Basal turbulent melting (for bits) 147 zMbb = MAX( 0.58_wp*(zdvo**0.8_wp)*(zSST+2._wp) / & 148 & ( zLbits**0.2_wp ) , 0._wp ) * z1_rday ! Basal turbulent melting (for bits) 145 149 zMbb = rn_rho_bergs * zAbits * zMbb ! in kg/s 146 150 zdMbitsM = MIN( zMbb*zdt , znMbits ) ! bergy bits mass lost to melting (kg) … … 163 167 zmelt = ( zdM - ( zdMbitsE - zdMbitsM ) ) * z1_dt ! kg/s 164 168 berg_grid%floating_melt(ii,ij) = berg_grid%floating_melt(ii,ij) + zmelt * z1_e1e2 ! kg/m2/s 165 zheat = zmelt * pt%heat_density ! kg/s x J/kg = J/s 166 berg_grid%calving_hflx (ii,ij) = berg_grid%calving_hflx (ii,ij) + zheat * z1_e1e2 ! W/m2 167 CALL icb_dia_melt( ii, ij, zMnew, zheat, this%mass_scaling, & 169 !! NB. The src_calving_hflx field is currently hardwired to zero in icb_stp, which means that the 170 !! heat density of the icebergs is zero and the heat content flux to the ocean from iceberg 171 !! melting is always zero. Leaving the term in the code until such a time as this is fixed. DS. 172 zheat_hcflux = zmelt * pt%heat_density ! heat content flux : kg/s x J/kg = J/s 173 zheat_latent = - zmelt * rLfus ! latent heat flux: kg/s x J/kg = J/s 174 berg_grid%calving_hflx (ii,ij) = berg_grid%calving_hflx (ii,ij) + ( zheat_hcflux + zheat_latent ) * z1_e1e2 ! W/m2 175 CALL icb_dia_melt( ii, ij, zMnew, zheat_hcflux, zheat_latent, this%mass_scaling, & 168 176 & zdM, zdMbitsE, zdMbitsM, zdMb, zdMe, & 169 177 & zdMv, z1_dt_e1e2 ) … … 211 219 IF(.NOT. ln_passive_mode ) THEN 212 220 emp (:,:) = emp (:,:) - berg_grid%floating_melt(:,:) 213 !! qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:) !!gm heat flux not yet properly coded ==>> need it, SOLVE that! 221 qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:) 214 222 ENDIF 215 223 !
Note: See TracChangeset
for help on using the changeset viewer.