Changeset 9888
- Timestamp:
- 2018-07-06T12:33:01+02:00 (7 years ago)
- Location:
- NEMO/trunk/src
- Files:
-
- 3 edited
-
ICE/icedyn_adv.F90 (modified) (1 diff)
-
ICE/icevar.F90 (modified) (1 diff)
-
OCE/BDY/bdyice.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icedyn_adv.F90
r9885 r9888 105 105 ! ==> 1) remove negative ice areas and volumes (conservation is ensure) 106 106 CALL ice_var_zapsmall 107 ! ==> 2) remove remaining negative advected fields (conservation is not preserved) 107 ! ==> 2) remove remaining negative advected fields (conservation is not preserved) => conservation issue 108 108 WHERE( v_s (:,:,:) < 0._wp ) v_s (:,:,:) = 0._wp 109 109 WHERE( sv_i(:,:,:) < 0._wp ) sv_i(:,:,:) = 0._wp -
NEMO/trunk/src/ICE/icevar.F90
r9725 r9888 853 853 DO ji = 1, npti 854 854 ztmelts = - tmut * sz_i_1d(ji,jk) 855 t_i_1d(ji,jk) = MIN( t_i_1d(ji,jk), ztmelts + rt0 ) ! Force t_i_1d to be lower than melting point 855 t_i_1d(ji,jk) = MIN( t_i_1d(ji,jk), ztmelts + rt0 ) ! Force t_i_1d to be lower than melting point => likely conservation issue 856 856 ! (sometimes zdf scheme produces abnormally high temperatures) 857 857 e_i_1d(ji,jk) = rhoic * ( cpic * ( ztmelts - ( t_i_1d(ji,jk) - rt0 ) ) & -
NEMO/trunk/src/OCE/BDY/bdyice.F90
r9885 r9888 139 139 CALL lbc_bdy_lnk( h_s(:,:,jl), 'T', 1., ib_bdy ) 140 140 ENDDO 141 ! retrieve at_i142 at_i(:,:) = 0._wp143 DO jl = 1, jpl144 at_i(:,:) = a_i(:,:,jl) + at_i(:,:)145 END DO146 141 147 142 DO jl = 1, jpl … … 162 157 ! ! do not make state variables dependent on velocity 163 158 ! 164 rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ii,ij) - 0.01 ) ) ! 0 if no ice 165 ! 166 ! concentration and thickness 167 a_i(ji,jj,jl) = a_i(ii,ij,jl) * rswitch 168 h_i(ji,jj,jl) = h_i(ii,ij,jl) * rswitch 169 h_s(ji,jj,jl) = h_s(ii,ij,jl) * rswitch 170 ! 171 ! Ice and snow volumes 172 v_i(ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl) 173 v_s(ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl) 174 ! 175 SELECT CASE( jpbound ) 176 ! 177 CASE( 0 ) ! velocity is inward 178 ! 179 ! Ice salinity, age, temperature 180 s_i (ji,jj,jl) = rswitch * rn_ice_sal(ib_bdy) + ( 1.0 - rswitch ) * rn_simin 181 oa_i(ji,jj,jl) = rswitch * rn_ice_age(ib_bdy) * a_i(ji,jj,jl) 182 t_su(ji,jj,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rn_ice_tem(ib_bdy) 159 IF( a_i(ii,ij,jl) > 0._wp ) THEN ! there is ice at the boundary 160 ! 161 a_i(ji,jj,jl) = a_i(ii,ij,jl) ! concentration 162 h_i(ji,jj,jl) = h_i(ii,ij,jl) ! thickness ice 163 h_s(ji,jj,jl) = h_s(ii,ij,jl) ! thickness snw 164 ! 165 SELECT CASE( jpbound ) 166 ! 167 CASE( 0 ) ! velocity is inward 168 ! 169 oa_i(ji,jj, jl) = rn_ice_age(ib_bdy) * a_i(ji,jj,jl) ! age 170 a_ip(ji,jj, jl) = 0._wp ! pond concentration 171 v_ip(ji,jj, jl) = 0._wp ! pond volume 172 t_su(ji,jj, jl) = rn_ice_tem(ib_bdy) ! temperature surface 173 t_s (ji,jj,:,jl) = rn_ice_tem(ib_bdy) ! temperature snw 174 t_i (ji,jj,:,jl) = rn_ice_tem(ib_bdy) ! temperature ice 175 s_i (ji,jj, jl) = rn_ice_sal(ib_bdy) ! salinity 176 sz_i(ji,jj,:,jl) = rn_ice_sal(ib_bdy) ! salinity profile 177 ! 178 CASE( 1 ) ! velocity is outward 179 ! 180 oa_i(ji,jj, jl) = oa_i(ii,ij, jl) ! age 181 a_ip(ji,jj, jl) = a_ip(ii,ij, jl) ! pond concentration 182 v_ip(ji,jj, jl) = v_ip(ii,ij, jl) ! pond volume 183 t_su(ji,jj, jl) = t_su(ii,ij, jl) ! temperature surface 184 t_s (ji,jj,:,jl) = t_s (ii,ij,:,jl) ! temperature snw 185 t_i (ji,jj,:,jl) = t_i (ii,ij,:,jl) ! temperature ice 186 s_i (ji,jj, jl) = s_i (ii,ij, jl) ! salinity 187 sz_i(ji,jj,:,jl) = sz_i(ii,ij,:,jl) ! salinity profile 188 ! 189 END SELECT 190 ! 191 IF( nn_icesal == 1 ) THEN ! if constant salinity 192 s_i (ji,jj ,jl) = rn_icesal 193 sz_i(ji,jj,:,jl) = rn_icesal 194 ENDIF 195 ! 196 ! global fields 197 v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl) ! volume ice 198 v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl) ! volume snw 199 sv_i(ji,jj,jl) = MIN( s_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content 183 200 DO jk = 1, nlay_s 184 t_s(ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0 201 e_s(ji,jj,jk,jl) = rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) ! enthalpy in J/m3 202 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s ! enthalpy in J/m2 203 END DO 204 DO jk = 1, nlay_i 205 ztmelts = - tmut * sz_i(ji,jj,jk,jl) ! Melting temperature in C 206 t_i(ji,jj,jk,jl) = MIN( t_i(ji,jj,jk,jl), ztmelts + rt0 ) ! Force t_i to be lower than melting point => likely conservation issue 207 ! 208 e_i(ji,jj,jk,jl) = rhoic * ( cpic * ( ztmelts - ( t_i(ji,jj,jk,jl) - rt0 ) ) & ! enthalpy in J/m3 209 & + lfus * ( 1._wp - ztmelts / ( t_i(ji,jj,jk,jl) - rt0 ) ) & 210 & - rcp * ztmelts ) 211 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * v_i(ji,jj,jl) * r1_nlay_i ! enthalpy in J/m2 185 212 END DO 186 DO jk = 1, nlay_i 187 t_i (ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0 188 sz_i(ji,jj,jk,jl) = rswitch * rn_ice_sal(ib_bdy) + ( 1.0 - rswitch ) * rn_simin 189 END DO 190 ! 191 ! Ice ponds 192 a_ip(ji,jj,jl) = 0._wp 193 v_ip(ji,jj,jl) = 0._wp 194 ! 195 CASE( 1 ) ! velocity is outward 196 ! 197 ! Ice salinity, age, temperature 198 s_i (ji,jj,jl) = rswitch * s_i (ii,ij,jl) + ( 1.0 - rswitch ) * rn_simin 199 oa_i(ji,jj,jl) = rswitch * oa_i(ii,ij,jl) 200 t_su(ji,jj,jl) = rswitch * t_su(ii,ij,jl) + ( 1.0 - rswitch ) * rt0 201 DO jk = 1, nlay_s 202 t_s(ji,jj,jk,jl) = rswitch * t_s(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0 203 END DO 204 DO jk = 1, nlay_i 205 t_i (ji,jj,jk,jl) = rswitch * t_i (ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0 206 sz_i(ji,jj,jk,jl) = rswitch * sz_i(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rn_simin 207 END DO 208 ! 209 ! Ice ponds 210 a_ip(ji,jj,jl) = rswitch * a_ip(ii,ij,jl) 211 v_ip(ji,jj,jl) = rswitch * v_ip(ii,ij,jl) 212 ! 213 END SELECT 214 ! 215 IF( nn_icesal == 1 ) THEN ! constant salinity : overwrite rn_icesal 216 s_i (ji,jj ,jl) = rn_icesal 217 sz_i(ji,jj,:,jl) = rn_icesal 213 ! 214 ELSE ! no ice at the boundary 215 ! 216 a_i (ji,jj, jl) = 0._wp 217 h_i (ji,jj, jl) = 0._wp 218 h_s (ji,jj, jl) = 0._wp 219 oa_i(ji,jj, jl) = 0._wp 220 a_ip(ji,jj, jl) = 0._wp 221 v_ip(ji,jj, jl) = 0._wp 222 t_su(ji,jj, jl) = rt0 223 t_s (ji,jj,:,jl) = rt0 224 t_i (ji,jj,:,jl) = rt0 225 226 IF( nn_icesal == 1 ) THEN ! if constant salinity 227 s_i (ji,jj ,jl) = rn_icesal 228 sz_i(ji,jj,:,jl) = rn_icesal 229 ELSE ! if variable salinity 230 s_i (ji,jj,jl) = rn_simin 231 sz_i(ji,jj,:,jl) = rn_simin 232 ENDIF 233 ! 234 ! global fields 235 v_i (ji,jj, jl) = 0._wp 236 v_s (ji,jj, jl) = 0._wp 237 sv_i(ji,jj, jl) = 0._wp 238 e_s (ji,jj,:,jl) = 0._wp 239 e_i (ji,jj,:,jl) = 0._wp 240 218 241 ENDIF 219 ! 220 ! contents 221 sv_i(ji,jj,jl) = MIN( s_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) 222 DO jk = 1, nlay_s 223 ! Snow energy of melting 224 e_s(ji,jj,jk,jl) = rswitch * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) 225 ! Multiply by volume, so that heat content in J/m2 226 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s 227 END DO 228 DO jk = 1, nlay_i 229 ztmelts = - tmut * sz_i(ji,jj,jk,jl) + rt0 !Melting temperature in K 230 ! heat content per unit volume 231 e_i(ji,jj,jk,jl) = rswitch * rhoic * & 232 ( cpic * ( ztmelts - t_i(ji,jj,jk,jl) ) & 233 + lfus * ( 1.0 - (ztmelts-rt0) / MIN((t_i(ji,jj,jk,jl)-rt0),-epsi20) ) & 234 - rcp * ( ztmelts - rt0 ) ) 235 ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2 236 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * h_i(ji,jj,jl) * r1_nlay_i 237 END DO 238 ! 242 239 243 END DO 240 244 ! 241 CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy ) 242 CALL lbc_bdy_lnk( h_i(:,:,jl), 'T', 1., ib_bdy ) 243 CALL lbc_bdy_lnk( h_s(:,:,jl), 'T', 1., ib_bdy ) 244 CALL lbc_bdy_lnk( v_i(:,:,jl), 'T', 1., ib_bdy ) 245 CALL lbc_bdy_lnk( v_s(:,:,jl), 'T', 1., ib_bdy ) 246 ! 245 CALL lbc_bdy_lnk( a_i (:,:,jl), 'T', 1., ib_bdy ) 246 CALL lbc_bdy_lnk( h_i (:,:,jl), 'T', 1., ib_bdy ) 247 CALL lbc_bdy_lnk( h_s (:,:,jl), 'T', 1., ib_bdy ) 248 CALL lbc_bdy_lnk( oa_i(:,:,jl), 'T', 1., ib_bdy ) 249 CALL lbc_bdy_lnk( a_ip(:,:,jl), 'T', 1., ib_bdy ) 250 CALL lbc_bdy_lnk( v_ip(:,:,jl), 'T', 1., ib_bdy ) 251 CALL lbc_bdy_lnk( s_i (:,:,jl), 'T', 1., ib_bdy ) 252 CALL lbc_bdy_lnk( t_su(:,:,jl), 'T', 1., ib_bdy ) 253 CALL lbc_bdy_lnk( v_i (:,:,jl), 'T', 1., ib_bdy ) 254 CALL lbc_bdy_lnk( v_s (:,:,jl), 'T', 1., ib_bdy ) 247 255 CALL lbc_bdy_lnk( sv_i(:,:,jl), 'T', 1., ib_bdy ) 248 CALL lbc_bdy_lnk( s_i(:,:,jl), 'T', 1., ib_bdy ) 249 CALL lbc_bdy_lnk( oa_i(:,:,jl), 'T', 1., ib_bdy ) 250 CALL lbc_bdy_lnk( t_su(:,:,jl), 'T', 1., ib_bdy ) 251 DO jk = 1, nlay_s 256 DO jk = 1, nlay_s 252 257 CALL lbc_bdy_lnk(t_s(:,:,jk,jl), 'T', 1., ib_bdy ) 253 258 CALL lbc_bdy_lnk(e_s(:,:,jk,jl), 'T', 1., ib_bdy ) … … 258 263 END DO 259 264 ! 260 CALL lbc_bdy_lnk( a_ip(:,:,jl), 'T', 1., ib_bdy ) 261 CALL lbc_bdy_lnk( v_ip(:,:,jl), 'T', 1., ib_bdy ) 262 ! 263 END DO !jl 265 END DO ! jl 264 266 ! 265 267 END SUBROUTINE bdy_ice_frs … … 309 311 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice 310 312 ! 311 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce)313 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then do not change u_ice) 312 314 u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + & 313 315 & u_ice(ji-1,jj) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + & 314 & u_ oce(ji ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) )316 & u_ice(ji ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 315 317 ELSE ! everywhere else 316 !u_ice(ji,jj) = u_oce(ji,jj)317 318 u_ice(ji,jj) = 0._wp 318 319 ENDIF 319 ! mask ice velocities320 rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01_wp ) ) ! 0 if no ice321 u_ice(ji,jj) = rswitch * u_ice(ji,jj)322 320 ! 323 321 END DO … … 336 334 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj-1) ) ) ! 0 if no ice 337 335 ! 338 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce)336 ! v_ice = v_ice of the adjacent grid point except if this grid point is ice-free (then do not change v_ice) 339 337 v_ice (ji,jj) = v_ice(ji,jj+1) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + & 340 338 & v_ice(ji,jj-1) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + & 341 & v_ oce(ji,jj ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) )339 & v_ice(ji,jj ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 342 340 ELSE ! everywhere else 343 !v_ice(ji,jj) = v_oce(ji,jj)344 341 v_ice(ji,jj) = 0._wp 345 342 ENDIF 346 ! mask ice velocities347 rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01 ) ) ! 0 if no ice348 v_ice(ji,jj) = rswitch * v_ice(ji,jj)349 343 ! 350 344 END DO
Note: See TracChangeset
for help on using the changeset viewer.
