Changeset 2715 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
- Timestamp:
- 2011-03-30T17:58:35+02:00 (13 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90 (modified) (14 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
r2528 r2715 1 1 MODULE limvar 2 !!----------------------------------------------------------------------3 !! 'key_lim3' LIM3 sea-ice model4 !!----------------------------------------------------------------------5 2 !!====================================================================== 6 3 !! *** MODULE limvar *** … … 32 29 !! - ot_i(jpi,jpj) !average ice age 33 30 !!====================================================================== 31 !! History : - ! 2006-01 (M. Vancoppenolle) Original code 32 !! 4.0 ! 2011-02 (G. Madec) dynamical allocation 33 !!---------------------------------------------------------------------- 34 34 #if defined key_lim3 35 35 !!---------------------------------------------------------------------- 36 !! * Modules used 37 USE dom_ice 36 !! 'key_lim3' LIM3 sea-ice model 37 !!---------------------------------------------------------------------- 38 !! lim_var_agg : 39 !! lim_var_glo2eqv : 40 !! lim_var_eqv2glo : 41 !! lim_var_salprof : 42 !! lim_var_salprof1d : 43 !! lim_var_bv : 44 !!---------------------------------------------------------------------- 38 45 USE par_oce ! ocean parameters 39 46 USE phycst ! physical constants (ocean directory) 40 47 USE sbc_oce ! Surface boundary condition: ocean fields 41 USE thd_ice 42 USE in_out_manager 43 USE ice 44 USE par_ice 48 USE ice ! LIM variables 49 USE par_ice ! LIM parameters 50 USE dom_ice ! LIM domain 51 USE thd_ice ! LIM thermodynamics 52 USE wrk_nemo ! workspace manager 53 USE in_out_manager ! I/O manager 54 USE lib_mpp ! MPP library 45 55 46 56 IMPLICIT NONE 47 57 PRIVATE 48 58 49 !! * Routine accessibility50 PUBLIC lim_var_agg51 PUBLIC lim_var_glo2eqv52 PUBLIC lim_var_eqv2glo53 PUBLIC lim_var_salprof54 PUBLIC lim_var_bv55 PUBLIC lim_var_salprof1d 56 57 !! * Module variables58 REAL(wp) :: & ! constant values59 epsi20 = 1e-20 , &60 epsi13 = 1e-13 , &61 zzero = 0.e0 , &62 zone = 1.e063 64 !!---------------------------------------------------------------------- 65 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010)59 PUBLIC lim_var_agg ! 60 PUBLIC lim_var_glo2eqv ! 61 PUBLIC lim_var_eqv2glo ! 62 PUBLIC lim_var_salprof ! 63 PUBLIC lim_var_bv ! 64 PUBLIC lim_var_salprof1d ! 65 66 REAL(wp) :: eps20 = 1.e-20_wp ! module constants 67 REAL(wp) :: eps16 = 1.e-16_wp ! - - 68 REAL(wp) :: eps13 = 1.e-13_wp ! - - 69 REAL(wp) :: eps10 = 1.e-10_wp ! - - 70 REAL(wp) :: eps06 = 1.e-06_wp ! - - 71 REAL(wp) :: zzero = 0.e0 ! - - 72 REAL(wp) :: zone = 1.e0 ! - - 73 74 !!---------------------------------------------------------------------- 75 !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011) 66 76 !! $Id$ 67 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 68 !!---------------------------------------------------------------------- 69 70 77 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 78 !!---------------------------------------------------------------------- 71 79 CONTAINS 72 80 73 SUBROUTINE lim_var_agg( n)81 SUBROUTINE lim_var_agg( kn ) 74 82 !!------------------------------------------------------------------ 75 83 !! *** ROUTINE lim_var_agg *** 76 !! ** Purpose : 77 !! This routine aggregates ice-thickness-category variables to 78 !! all-ice variables 79 !! i.e. it turns VGLO into VAGG 84 !! 85 !! ** Purpose : aggregates ice-thickness-category variables to all-ice variables 86 !! i.e. it turns VGLO into VAGG 80 87 !! ** Method : 81 88 !! 82 !! ** Arguments :83 !! kideb , kiut : Starting and ending points on which the84 !! the computation is applied85 !!86 !! ** Inputs / Ouputs : (global commons)87 89 !! ** Arguments : n = 1, at_i vt_i only 88 90 !! n = 2 everything 89 91 !! 90 !! ** External :91 !!92 !! ** References :93 !!94 !! ** History :95 !! (01-2006) Martin Vancoppenolle, UCL-ASTR96 !!97 92 !! note : you could add an argument when you need only at_i, vt_i 98 93 !! and when you need everything 99 94 !!------------------------------------------------------------------ 100 !! * Arguments 101 102 !! * Local variables 103 INTEGER :: ji, & ! spatial dummy loop index 104 jj, & ! spatial dummy loop index 105 jk, & ! vertical layering dummy loop index 106 jl ! ice category dummy loop index 107 108 REAL :: zeps, epsi16, zinda, epsi06 109 110 INTEGER, INTENT( in ) :: n ! describes what is needed 111 112 !!-- End of declarations 113 !!---------------------------------------------------------------------------------------------- 114 zeps = 1.0e-13 115 epsi16 = 1.0e-16 116 epsi06 = 1.0e-6 117 118 !------------------ 119 ! Zero everything 120 !------------------ 121 122 vt_i(:,:) = 0.0 123 vt_s(:,:) = 0.0 124 at_i(:,:) = 0.0 125 ato_i(:,:) = 1.0 126 127 IF ( n .GT. 1 ) THEN 128 et_s(:,:) = 0.0 129 ot_i(:,:) = 0.0 130 smt_i(:,:) = 0.0 131 et_i(:,:) = 0.0 132 ENDIF 95 INTEGER, INTENT( in ) :: kn ! =1 at_i & vt only ; = what is needed 96 ! 97 INTEGER :: ji, jj, jk, jl ! dummy loop indices 98 REAL(wp) :: zinda 99 !!------------------------------------------------------------------ 133 100 134 101 !-------------------- 135 102 ! Compute variables 136 103 !-------------------- 137 104 vt_i (:,:) = 0._wp 105 vt_s (:,:) = 0._wp 106 at_i (:,:) = 0._wp 107 ato_i(:,:) = 1._wp 108 ! 138 109 DO jl = 1, jpl 139 110 DO jj = 1, jpj 140 111 DO ji = 1, jpi 141 112 ! 142 113 vt_i(ji,jj) = vt_i(ji,jj) + v_i(ji,jj,jl) ! ice volume 143 114 vt_s(ji,jj) = vt_s(ji,jj) + v_s(ji,jj,jl) ! snow volume 144 115 at_i(ji,jj) = at_i(ji,jj) + a_i(ji,jj,jl) ! ice concentration 145 116 ! 146 117 zinda = MAX( zzero , SIGN( zone , at_i(ji,jj) - 0.10 ) ) 147 icethi(ji,jj) = vt_i(ji,jj) / MAX(at_i(ji,jj),epsi16)*zinda 148 ! ice thickness 118 icethi(ji,jj) = vt_i(ji,jj) / MAX( at_i(ji,jj) , eps16 ) * zinda ! ice thickness 149 119 END DO 150 120 END DO … … 153 123 DO jj = 1, jpj 154 124 DO ji = 1, jpi 155 ato_i(ji,jj) = MAX(1.0 - at_i(ji,jj), 0.0) ! open water fraction 156 END DO 157 END DO 158 159 IF ( n .GT. 1 ) THEN 160 125 ato_i(ji,jj) = MAX( 1._wp - at_i(ji,jj), 0._wp ) ! open water fraction 126 END DO 127 END DO 128 129 IF( kn > 1 ) THEN 130 et_s (:,:) = 0._wp 131 ot_i (:,:) = 0._wp 132 smt_i(:,:) = 0._wp 133 et_i (:,:) = 0._wp 134 ! 161 135 DO jl = 1, jpl 162 136 DO jj = 1, jpj 163 137 DO ji = 1, jpi 164 et_s(ji,jj) = et_s(ji,jj) + & ! snow heat content 165 e_s(ji,jj,1,jl) 138 et_s(ji,jj) = et_s(ji,jj) + e_s(ji,jj,1,jl) ! snow heat content 166 139 zinda = MAX( zzero , SIGN( zone , vt_i(ji,jj) - 0.10 ) ) 167 smt_i(ji,jj) = smt_i(ji,jj) + & ! ice salinity 168 smv_i(ji,jj,jl) / MAX( vt_i(ji,jj) , zeps ) * & 169 zinda 140 smt_i(ji,jj) = smt_i(ji,jj) + smv_i(ji,jj,jl) / MAX( vt_i(ji,jj) , eps13 ) * zinda ! ice salinity 170 141 zinda = MAX( zzero , SIGN( zone , at_i(ji,jj) - 0.10 ) ) 171 ot_i(ji,jj) = ot_i(ji,jj) + & ! ice age 172 oa_i(ji,jj,jl) / MAX( at_i(ji,jj) , zeps ) * & 173 zinda 174 END DO 175 END DO 176 END DO 177 142 ot_i(ji,jj) = ot_i(ji,jj) + oa_i(ji,jj,jl) / MAX( at_i(ji,jj) , eps13 ) * zinda ! ice age 143 END DO 144 END DO 145 END DO 146 ! 178 147 DO jl = 1, jpl 179 148 DO jk = 1, nlay_i 180 DO jj = 1, jpj 181 DO ji = 1, jpi 182 et_i(ji,jj) = et_i(ji,jj) + e_i(ji,jj,jk,jl) ! ice heat 183 ! content 184 END DO 185 END DO 186 END DO 187 END DO 188 189 ENDIF ! n .GT. 1 190 149 et_i(:,:) = et_i(:,:) + e_i(:,:,jk,jl) ! ice heat content 150 END DO 151 END DO 152 ! 153 ENDIF 154 ! 191 155 END SUBROUTINE lim_var_agg 192 156 193 !==============================================================================194 157 195 158 SUBROUTINE lim_var_glo2eqv 196 159 !!------------------------------------------------------------------ 197 !! *** ROUTINE lim_var_glo2eqv ***' 198 !! ** Purpose : 199 !! This routine computes equivalent variables as function of 200 !! global variables 201 !! i.e. it turns VGLO into VEQV 202 !! ** Method : 203 !! 204 !! ** Arguments : 205 !! kideb , kiut : Starting and ending points on which the 206 !! the computation is applied 207 !! 208 !! ** Inputs / Ouputs : 209 !! 210 !! ** External : 211 !! 212 !! ** References : 213 !! 214 !! ** History : 215 !! (01-2006) Martin Vancoppenolle, UCL-ASTR 216 !! 217 !!------------------------------------------------------------------ 218 219 !! * Local variables 220 INTEGER :: ji, & ! spatial dummy loop index 221 jj, & ! spatial dummy loop index 222 jk, & ! vertical layering dummy loop index 223 jl ! ice category dummy loop index 224 225 REAL :: zq_i, zaaa, zbbb, zccc, zdiscrim, & 226 ztmelts, zindb, zq_s, zfac1, zfac2 227 228 REAL :: zeps, epsi06 229 230 zeps = 1.0e-10 231 epsi06 = 1.0e-06 232 233 !!-- End of declarations 234 !!------------------------------------------------------------------------------ 160 !! *** ROUTINE lim_var_glo2eqv *** 161 !! 162 !! ** Purpose : computes equivalent variables as function of global variables 163 !! i.e. it turns VGLO into VEQV 164 !!------------------------------------------------------------------ 165 INTEGER :: ji, jj, jk, jl ! dummy loop indices 166 REAL(wp) :: zq_i, zaaa, zbbb, zccc, zdiscrim ! local scalars 167 REAL(wp) :: ztmelts, zindb, zq_s, zfac1, zfac2 ! - - 168 !!------------------------------------------------------------------ 235 169 236 170 !------------------------------------------------------- 237 171 ! Ice thickness, snow thickness, ice salinity, ice age 238 172 !------------------------------------------------------- 239 !CDIR NOVERRCHK240 173 DO jl = 1, jpl 241 !CDIR NOVERRCHK242 174 DO jj = 1, jpj 243 !CDIR NOVERRCHK244 175 DO ji = 1, jpi 245 zindb = 1.0-MAX(0.0,SIGN(1.0,- a_i(ji,jj,jl))) !0 if no ice and 1 if yes 246 ht_i(ji,jj,jl) = v_i(ji,jj,jl) / MAX( a_i(ji,jj,jl) , zeps ) * zindb 247 ht_s(ji,jj,jl) = v_s(ji,jj,jl) / MAX( a_i(ji,jj,jl) , zeps ) * zindb 248 o_i(ji,jj,jl) = oa_i(ji,jj,jl) / MAX( a_i(ji,jj,jl) , zeps ) * zindb 249 END DO 250 END DO 251 END DO 252 253 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )THEN 254 255 !CDIR NOVERRCHK 176 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp,- a_i(ji,jj,jl) ) ) !0 if no ice and 1 if yes 177 ht_i(ji,jj,jl) = v_i (ji,jj,jl) / MAX( a_i(ji,jj,jl) , eps10 ) * zindb 178 ht_s(ji,jj,jl) = v_s (ji,jj,jl) / MAX( a_i(ji,jj,jl) , eps10 ) * zindb 179 o_i(ji,jj,jl) = oa_i(ji,jj,jl) / MAX( a_i(ji,jj,jl) , eps10 ) * zindb 180 END DO 181 END DO 182 END DO 183 184 IF( num_sal == 2 .OR. num_sal == 4 )THEN 256 185 DO jl = 1, jpl 257 !CDIR NOVERRCHK 258 DO jj = 1, jpj 259 !CDIR NOVERRCHK 260 DO ji = 1, jpi 261 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes 262 sm_i(ji,jj,jl) = smv_i(ji,jj,jl) / MAX(v_i(ji,jj,jl),zeps) * zindb 263 END DO 264 END DO 265 END DO 266 186 DO jj = 1, jpj 187 DO ji = 1, jpi 188 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp,- a_i(ji,jj,jl) ) ) !0 if no ice and 1 if yes 189 sm_i(ji,jj,jl) = smv_i(ji,jj,jl) / MAX( v_i(ji,jj,jl) , eps10 ) * zindb 190 END DO 191 END DO 192 END DO 267 193 ENDIF 268 194 269 ! salinity profile 270 CALL lim_var_salprof 195 CALL lim_var_salprof ! salinity profile 271 196 272 197 !------------------- … … 281 206 !CDIR NOVERRCHK 282 207 DO ji = 1, jpi 283 !Energy of melting q(S,T) [J.m-3] 284 zq_i = e_i(ji,jj,jk,jl) / area(ji,jj) / & 285 MAX( v_i(ji,jj,jl) , epsi06 ) * nlay_i 286 ! zindb = 0 if no ice and 1 if yes 287 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_i(ji,jj,jl) ) ) 288 !convert units ! very important that this line is here 289 zq_i = zq_i * unit_fac * zindb 290 !Ice layer melt temperature 291 ztmelts = -tmut*s_i(ji,jj,jk,jl) + rtt 292 !Conversion q(S,T) -> T (second order equation) 293 zaaa = cpic 294 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + & 295 zq_i / rhoic - lfus 208 ! ! Energy of melting q(S,T) [J.m-3] 209 zq_i = e_i(ji,jj,jk,jl) / area(ji,jj) / MAX( v_i(ji,jj,jl) , eps06 ) * REAL(nlay_i,wp) 210 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_i(ji,jj,jl) ) ) ! zindb = 0 if no ice and 1 if yes 211 zq_i = zq_i * unit_fac * zindb !convert units 212 ztmelts = -tmut * s_i(ji,jj,jk,jl) + rtt ! Ice layer melt temperature 213 ! 214 zaaa = cpic ! Conversion q(S,T) -> T (second order equation) 215 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + zq_i / rhoic - lfus 296 216 zccc = lfus * (ztmelts-rtt) 297 zdiscrim = SQRT( MAX(zbbb*zbbb - 4.0*zaaa*zccc,0.0) ) 298 t_i(ji,jj,jk,jl) = rtt + zindb *( - zbbb - zdiscrim ) / & 299 ( 2.0 *zaaa ) 300 t_i(ji,jj,jk,jl) = MIN( rtt, MAX(173.15, t_i(ji,jj,jk,jl) ) ) 217 zdiscrim = SQRT( MAX(zbbb*zbbb - 4._wp*zaaa*zccc , 0._wp) ) 218 t_i(ji,jj,jk,jl) = rtt + zindb *( - zbbb - zdiscrim ) / ( 2.0 *zaaa ) 219 t_i(ji,jj,jk,jl) = MIN( rtt, MAX( 173.15_wp, t_i(ji,jj,jk,jl) ) ) ! 100-rtt < t_i < rtt 301 220 END DO 302 221 END DO … … 307 226 ! Snow temperatures 308 227 !-------------------- 309 zfac1 = 1. / ( rhosn * cpic )228 zfac1 = 1._wp / ( rhosn * cpic ) 310 229 zfac2 = lfus / cpic 311 !CDIR NOVERRCHK312 230 DO jl = 1, jpl 313 !CDIR NOVERRCHK314 231 DO jk = 1, nlay_s 315 !CDIR NOVERRCHK 316 DO jj = 1, jpj 317 !CDIR NOVERRCHK 232 DO jj = 1, jpj 318 233 DO ji = 1, jpi 319 234 !Energy of melting q(S,T) [J.m-3] 320 zq_s = e_s(ji,jj,jk,jl) / area(ji,jj) / & 321 MAX( v_s(ji,jj,jl) , epsi06 ) * nlay_s 322 ! zindb = 0 if no ice and 1 if yes 323 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_s(ji,jj,jl) ) ) 324 !convert units ! very important that this line is here 325 zq_s = zq_s * unit_fac * zindb 235 zq_s = e_s(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_s(ji,jj,jl) , eps06 ) ) * REAL(nlay_s,wp) 236 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - v_s(ji,jj,jl) ) ) ! zindb = 0 if no ice and 1 if yes 237 zq_s = zq_s * unit_fac * zindb ! convert units 238 ! 326 239 t_s(ji,jj,jk,jl) = rtt + zindb * ( - zfac1 * zq_s + zfac2 ) 327 t_s(ji,jj,jk,jl) = MIN( rtt, MAX(173.15, t_s(ji,jj,jk,jl) ) ) 328 240 t_s(ji,jj,jk,jl) = MIN( rtt, MAX( 173.15, t_s(ji,jj,jk,jl) ) ) ! 100-rtt < t_i < rtt 329 241 END DO 330 242 END DO … … 335 247 ! Mean temperature 336 248 !------------------- 337 tm_i(:,:) = 0.0 338 !CDIR NOVERRCHK 249 tm_i(:,:) = 0._wp 339 250 DO jl = 1, jpl 340 !CDIR NOVERRCHK341 251 DO jk = 1, nlay_i 342 !CDIR NOVERRCHK 343 DO jj = 1, jpj 344 !CDIR NOVERRCHK 345 DO ji = 1, jpi 346 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) 347 zindb = zindb*1.0-MAX(0.0,SIGN(1.0,-v_i(ji,jj,jl))) 348 tm_i(ji,jj) = tm_i(ji,jj) + t_i(ji,jj,jk,jl)*v_i(ji,jj,jl) / & 349 REAL(nlay_i) / MAX( vt_i(ji,jj) , zeps ) 350 END DO 351 END DO 352 END DO 353 END DO 354 252 DO jj = 1, jpj 253 DO ji = 1, jpi 254 zindb = ( 1._wp - MAX( 0._wp , SIGN( 1._wp , -a_i(ji,jj,jl) ) ) ) & 255 & * ( 1._wp - MAX( 0._wp , SIGN( 1._wp , -v_i(ji,jj,jl) ) ) ) 256 tm_i(ji,jj) = tm_i(ji,jj) + t_i(ji,jj,jk,jl) * v_i(ji,jj,jl) & 257 & / ( REAL(nlay_i,wp) * MAX( vt_i(ji,jj) , eps10 ) ) 258 END DO 259 END DO 260 END DO 261 END DO 262 ! 355 263 END SUBROUTINE lim_var_glo2eqv 356 264 357 !===============================================================================358 265 359 266 SUBROUTINE lim_var_eqv2glo 360 267 !!------------------------------------------------------------------ 361 !! *** ROUTINE lim_var_eqv2glo ***' 362 !! ** Purpose : 363 !! This routine computes global variables as function of 364 !! equivalent variables 365 !! i.e. it turns VEQV into VGLO 268 !! *** ROUTINE lim_var_eqv2glo *** 269 !! 270 !! ** Purpose : computes global variables as function of equivalent variables 271 !! i.e. it turns VEQV into VGLO 366 272 !! ** Method : 367 273 !! 368 !! ** Arguments : 369 !! 370 !! ** Inputs / Ouputs : (global commons) 371 !! 372 !! ** External : 373 !! 374 !! ** References : 375 !! 376 !! ** History : 377 !! (01-2006) Martin Vancoppenolle, UCL-ASTR 378 !! Take it easy man 379 !! Life is just a simple game, between 380 !! ups / and downs \ :@) 381 !! 382 !!------------------------------------------------------------------ 383 274 !! ** History : (01-2006) Martin Vancoppenolle, UCL-ASTR 275 !!------------------------------------------------------------------ 276 ! 384 277 v_i(:,:,:) = ht_i(:,:,:) * a_i(:,:,:) 385 278 v_s(:,:,:) = ht_s(:,:,:) * a_i(:,:,:) 386 279 smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:) 387 280 oa_i (:,:,:) = o_i (:,:,:) * a_i(:,:,:) 388 281 ! 389 282 END SUBROUTINE lim_var_eqv2glo 390 283 391 !===============================================================================392 284 393 285 SUBROUTINE lim_var_salprof 394 286 !!------------------------------------------------------------------ 395 !! *** ROUTINE lim_var_salprof ***' 396 !! ** Purpose : 397 !! This routine computes salinity profile in function of 398 !! bulk salinity 287 !! *** ROUTINE lim_var_salprof *** 288 !! 289 !! ** Purpose : computes salinity profile in function of bulk salinity 399 290 !! 400 291 !! ** Method : If bulk salinity greater than s_i_1, … … 406 297 !! 407 298 !! ** References : Vancoppenolle et al., 2007 (in preparation) 408 !! 409 !! ** History : 410 !! (08-2006) Martin Vancoppenolle, UCL-ASTR 411 !! Take it easy man 412 !! Life is just a simple game, between ups 413 !! / and downs \ :@) 414 !! 415 !!------------------------------------------------------------------ 416 !! * Arguments 417 418 !! * Local variables 419 INTEGER :: & 420 ji , & !: spatial dummy loop index 421 jj , & !: spatial dummy loop index 422 jk , & !: vertical layering dummy loop index 423 jl !: ice category dummy loop index 424 425 REAL(wp) :: & 426 dummy_fac0 , & !: dummy factor used in computations 427 dummy_fac1 , & !: dummy factor used in computations 428 dummy_fac , & !: dummy factor used in computations 429 zind0 , & !: switch, = 1 if sm_i lt s_i_0 430 zind01 , & !: switch, = 1 if sm_i between s_i_0 and s_i_1 431 zindbal , & !: switch, = 1, if 2*sm_i gt sss_m 432 zargtemp !: dummy factor 433 434 REAL(wp), DIMENSION(nlay_i) :: & 435 zs_zero !: linear salinity profile for salinities under 436 !: s_i_0 437 438 REAL(wp), DIMENSION(jpi,jpj,jpl) :: & 439 z_slope_s , & !: slope of the salinity profile 440 zalpha !: weight factor for s between s_i_0 and s_i_1 441 442 !!-- End of declarations 443 !!------------------------------------------------------------------------------ 299 !!------------------------------------------------------------------ 300 INTEGER :: ji, jj, jk, jl ! dummy loop index 301 REAL(wp) :: dummy_fac0, dummy_fac1, dummy_fac, zsal ! local scalar 302 REAL(wp) :: zind0, zind01, zindbal, zargtemp , zs_zero ! - - 303 ! 304 REAL(wp), POINTER, DIMENSION(:,:,:) :: z_slope_s, zalpha ! 3D pointer 305 !!------------------------------------------------------------------ 306 307 IF( wrk_in_use( 2, 1,2 ) ) THEN 308 CALL ctl_stop( 'lim_var_salprof: requested workspace arrays unavailable' ) ; RETURN 309 END IF 310 311 z_slope_s => wrk_3d_1(:,:,1:jpl) ! slope of the salinity profile 312 zalpha => wrk_3d_2(:,:,1:jpl) ! weight factor for s between s_i_0 and s_i_1 444 313 445 314 !--------------------------------------- 446 315 ! Vertically constant, constant in time 447 316 !--------------------------------------- 448 449 IF ( num_sal .EQ. 1 ) THEN 450 451 s_i(:,:,:,:) = bulk_sal 452 453 ENDIF 317 IF( num_sal == 1 ) s_i(:,:,:,:) = bulk_sal 454 318 455 319 !----------------------------------- … … 457 321 !----------------------------------- 458 322 459 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )THEN460 323 IF( num_sal == 2 .OR. num_sal == 4 ) THEN 324 ! 461 325 DO jk = 1, nlay_i 462 326 s_i(:,:,jk,:) = sm_i(:,:,:) 463 END DO ! jk 464 465 ! Slope of the linear profile zs_zero 466 !------------------------------------- 327 END DO 328 ! 329 DO jl = 1, jpl ! Slope of the linear profile 330 DO jj = 1, jpj 331 DO ji = 1, jpi 332 z_slope_s(ji,jj,jl) = 2._wp * sm_i(ji,jj,jl) / MAX( 0.01 , ht_i(ji,jj,jl) ) 333 END DO 334 END DO 335 END DO 336 ! 337 dummy_fac0 = 1._wp / ( s_i_0 - s_i_1 ) ! Weighting factor between zs_zero and zs_inf 338 dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 ) 339 340 zalpha(:,:,:) = 0._wp 467 341 DO jl = 1, jpl 468 342 DO jj = 1, jpj 469 DO ji = 1, jpi470 z_slope_s(ji,jj,jl) = 2.0 * sm_i(ji,jj,jl) / MAX( 0.01 &471 , ht_i(ji,jj,jl) )472 END DO ! ji473 END DO ! jj474 END DO ! jl475 476 ! Weighting factor between zs_zero and zs_inf477 !---------------------------------------------478 dummy_fac0 = 1. / ( ( s_i_0 - s_i_1 ) )479 dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 )480 481 zalpha(:,:,:) = 0.0482 483 !CDIR NOVERRCHK484 DO jl = 1, jpl485 !CDIR NOVERRCHK486 DO jj = 1, jpj487 !CDIR NOVERRCHK488 343 DO ji = 1, jpi 489 344 ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise 490 345 zind0 = MAX( 0.0 , SIGN( 1.0 , s_i_0 - sm_i(ji,jj,jl) ) ) 491 346 ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws 492 zind01 = ( 1.0 - zind0 ) * & 493 MAX( 0.0 , SIGN( 1.0 , s_i_1 - sm_i(ji,jj,jl) ) ) 347 zind01 = ( 1.0 - zind0 ) * MAX( 0.0 , SIGN( 1.0 , s_i_1 - sm_i(ji,jj,jl) ) ) 494 348 ! If 2.sm_i GE sss_m then zindbal = 1 495 zindbal = MAX( 0.0 , SIGN( 1.0 , 2. * sm_i(ji,jj,jl) - & 496 sss_m(ji,jj) ) ) 497 zalpha(ji,jj,jl) = zind0 * 1.0 & 498 + zind01 * ( sm_i(ji,jj,jl) * dummy_fac0 + & 499 dummy_fac1 ) 349 zindbal = MAX( 0.0 , SIGN( 1.0 , 2. * sm_i(ji,jj,jl) - sss_m(ji,jj) ) ) 350 zalpha(ji,jj,jl) = zind0 * 1.0 + zind01 * ( sm_i(ji,jj,jl) * dummy_fac0 + dummy_fac1 ) 500 351 zalpha(ji,jj,jl) = zalpha(ji,jj,jl) * ( 1.0 - zindbal ) 501 352 END DO … … 503 354 END DO 504 355 505 ! Computation of the profile 506 !---------------------------- 507 dummy_fac = 1. / nlay_i 508 356 dummy_fac = 1._wp / nlay_i ! Computation of the profile 509 357 DO jl = 1, jpl 510 358 DO jk = 1, nlay_i 511 359 DO jj = 1, jpj 512 360 DO ji = 1, jpi 513 ! linear profile with 0 at the surface 514 zs_zero(jk) = z_slope_s(ji,jj,jl) * ( jk - 1./2. ) * & 515 ht_i(ji,jj,jl) * dummy_fac 516 ! weighting the profile 517 s_i(ji,jj,jk,jl) = zalpha(ji,jj,jl) * zs_zero(jk) + & 518 ( 1.0 - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl) 361 ! ! linear profile with 0 at the surface 362 zs_zero = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * ht_i(ji,jj,jl) * dummy_fac 363 ! ! weighting the profile 364 s_i(ji,jj,jk,jl) = zalpha(ji,jj,jl) * zs_zero + ( 1._wp - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl) 519 365 END DO ! ji 520 366 END DO ! jj … … 527 373 ! Vertically varying salinity profile, constant in time 528 374 !------------------------------------------------------- 529 ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30) 530 531 IF ( num_sal .EQ. 3 ) THEN 532 533 sm_i(:,:,:) = 2.30 534 535 !CDIR NOVERRCHK 375 376 IF( num_sal == 3 ) THEN ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30) 377 ! 378 sm_i(:,:,:) = 2.30_wp 379 ! 536 380 DO jl = 1, jpl 537 381 !CDIR NOVERRCHK 538 382 DO jk = 1, nlay_i 539 !CDIR NOVERRCHK 540 DO jj = 1, jpj 541 !CDIR NOVERRCHK 542 DO ji = 1, jpi 543 zargtemp = ( jk - 0.5 ) / nlay_i 544 s_i(ji,jj,jk,jl) = 1.6 - 1.6 * COS( 3.14169265 * & 545 ( zargtemp**(0.407/ & 546 ( 0.573 + zargtemp ) ) ) ) 547 END DO ! ji 548 END DO ! jj 549 END DO ! jk 550 END DO ! jl 383 zargtemp = ( REAL(jk,wp) - 0.5_wp ) / REAL(nlay_i,wp) 384 zsal = 1.6_wp * ( 1._wp - COS( rpi * zargtemp**(0.407_wp/(0.573_wp+zargtemp)) ) ) 385 s_i(:,:,jk,jl) = zsal 386 END DO 387 END DO 551 388 552 389 ENDIF ! num_sal 553 390 ! 391 IF( wrk_not_released(2, 1,2) ) CALL ctl_stop('lim_var_salprof: failed to release workspace arrays.') 392 ! 554 393 END SUBROUTINE lim_var_salprof 555 394 556 !===============================================================================557 395 558 396 SUBROUTINE lim_var_bv 559 397 !!------------------------------------------------------------------ 560 !! *** ROUTINE lim_var_bv *** '561 !! ** Purpose :562 !! This routinecomputes mean brine volume (%) in sea ice398 !! *** ROUTINE lim_var_bv *** 399 !! 400 !! ** Purpose : computes mean brine volume (%) in sea ice 563 401 !! 564 402 !! ** Method : e = - 0.054 * S (ppt) / T (C) 565 403 !! 566 !! ** Arguments : 567 !! 568 !! ** Inputs / Ouputs : (global commons) 569 !! 570 !! ** External : 571 !! 572 !! ** References : Vancoppenolle et al., JGR, 2007 573 !! 574 !! ** History : 575 !! (08-2006) Martin Vancoppenolle, UCL-ASTR 576 !! 577 !!------------------------------------------------------------------ 578 !! * Arguments 579 580 !! * Local variables 581 INTEGER :: ji, & ! spatial dummy loop index 582 jj, & ! spatial dummy loop index 583 jk, & ! vertical layering dummy loop index 584 jl ! ice category dummy loop index 585 586 REAL :: zbvi, & ! brine volume for a single ice category 587 zeps, & ! very small value 588 zindb ! is there ice or not 589 590 !!-- End of declarations 591 !!------------------------------------------------------------------------------ 592 593 zeps = 1.0e-13 594 bv_i(:,:) = 0.0 595 !CDIR NOVERRCHK 404 !! References : Vancoppenolle et al., JGR, 2007 405 !!------------------------------------------------------------------ 406 INTEGER :: ji, jj, jk, jl ! dummy loop indices 407 REAL(wp) :: zbvi, zindb ! local scalars 408 !!------------------------------------------------------------------ 409 ! 410 bv_i(:,:) = 0._wp 596 411 DO jl = 1, jpl 597 !CDIR NOVERRCHK598 412 DO jk = 1, nlay_i 599 !CDIR NOVERRCHK 600 DO jj = 1, jpj 601 !CDIR NOVERRCHK 602 DO ji = 1, jpi 603 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes 604 zbvi = - zindb * tmut *s_i(ji,jj,jk,jl) / & 605 MIN( t_i(ji,jj,jk,jl) - 273.15 , zeps ) & 606 * v_i(ji,jj,jl) / REAL(nlay_i) 607 bv_i(ji,jj) = bv_i(ji,jj) + zbvi & 608 / MAX( vt_i(ji,jj) , zeps ) 609 END DO 610 END DO 611 END DO 612 END DO 613 413 DO jj = 1, jpj 414 DO ji = 1, jpi 415 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes 416 zbvi = - zindb * tmut * s_i(ji,jj,jk,jl) / MIN( t_i(ji,jj,jk,jl) - 273.15 , eps13 ) & 417 & * v_i(ji,jj,jl) / REAL(nlay_i,wp) 418 bv_i(ji,jj) = bv_i(ji,jj) + zbvi / MAX( vt_i(ji,jj) , eps13 ) 419 END DO 420 END DO 421 END DO 422 END DO 423 ! 614 424 END SUBROUTINE lim_var_bv 615 425 616 !=============================================================================== 617 618 SUBROUTINE lim_var_salprof1d(kideb,kiut) 426 427 SUBROUTINE lim_var_salprof1d( kideb, kiut ) 619 428 !!------------------------------------------------------------------- 620 429 !! *** ROUTINE lim_thd_salprof1d *** 621 430 !! 622 431 !! ** Purpose : 1d computation of the sea ice salinity profile 623 !! Works with 1d vectors and is used by thermodynamic 624 !! modules 625 !! 626 !! history : 627 !! 3.0 ! May 2007 M. Vancoppenolle Original code 432 !! Works with 1d vectors and is used by thermodynamic modules 628 433 !!------------------------------------------------------------------- 629 INTEGER, INTENT(in) :: & 630 kideb, kiut ! thickness category index 631 632 INTEGER :: & 633 ji, jk, & ! geographic and layer index 634 zji, zjj 635 636 REAL(wp) :: & 637 dummy_fac0, & ! dummy factors 638 dummy_fac1, & 639 dummy_fac2, & 640 zalpha , & ! weighting factor 641 zind0 , & ! switches as in limvar 642 zind01 , & ! switch 643 zindbal , & ! switch if in freshwater area 644 zargtemp 645 646 REAL(wp), DIMENSION(jpij) :: & 647 z_slope_s 648 649 REAL(wp), DIMENSION(jpij,jkmax) :: & 650 zs_zero 651 !!------------------------------------------------------------------- 434 INTEGER, INTENT(in) :: kideb, kiut ! thickness category index 435 ! 436 INTEGER :: ji, jk ! dummy loop indices 437 INTEGER :: zji, zjj ! local integers 438 REAL(wp) :: dummy_fac0, dummy_fac1, dummy_fac2, zargtemp, zsal ! local scalars 439 REAL(wp) :: zalpha, zind0, zind01, zindbal, zs_zero ! - - 440 ! 441 REAL(wp), POINTER, DIMENSION(:) :: z_slope_s 442 !!--------------------------------------------------------------------- 443 444 IF( wrk_in_use(1, 1) ) THEN 445 CALL ctl_stop('lim_var_salprof1d : requestead workspace arrays unavailable.') ; RETURN 446 END IF 447 ! Set-up pointers to sub-arrays of workspace arrays 448 z_slope_s => wrk_1d_1 (1:jpij) 652 449 653 450 !--------------------------------------- 654 451 ! Vertically constant, constant in time 655 452 !--------------------------------------- 656 657 IF ( num_sal .EQ. 1 ) THEN 658 659 s_i_b(:,:) = bulk_sal 660 661 ENDIF 453 IF( num_sal == 1 ) s_i_b(:,:) = bulk_sal 662 454 663 455 !------------------------------------------------------ … … 665 457 !------------------------------------------------------ 666 458 667 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN 668 669 ! Slope of the linear profile zs_zero 670 !------------------------------------- 671 !CDIR NOVERRCHK 672 DO ji = kideb, kiut 673 z_slope_s(ji) = 2.0 * sm_i_b(ji) / MAX( 0.01 & 674 , ht_i_b(ji) ) 675 END DO ! ji 459 IF( num_sal == 2 .OR. num_sal == 4 ) THEN 460 ! 461 DO ji = kideb, kiut ! Slope of the linear profile zs_zero 462 z_slope_s(ji) = 2._wp * sm_i_b(ji) / MAX( 0.01 , ht_i_b(ji) ) 463 END DO 676 464 677 465 ! Weighting factor between zs_zero and zs_inf 678 466 !--------------------------------------------- 679 dummy_fac0 = 1. / ( ( s_i_0 - s_i_1 ))467 dummy_fac0 = 1._wp / ( s_i_0 - s_i_1 ) 680 468 dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 ) 681 dummy_fac2 = 1. / nlay_i469 dummy_fac2 = 1._wp / REAL(nlay_i,wp) 682 470 683 471 !CDIR NOVERRCHK … … 685 473 !CDIR NOVERRCHK 686 474 DO ji = kideb, kiut 687 zji = MOD( npb(ji) - 1, jpi ) + 1 688 zjj = ( npb(ji) - 1 ) / jpi + 1 689 zalpha = 0.0 475 zji = MOD( npb(ji) - 1 , jpi ) + 1 476 zjj = ( npb(ji) - 1 ) / jpi + 1 690 477 ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise 691 zind0 = MAX( 0. 0 , SIGN( 1.0, s_i_0 - sm_i_b(ji) ) )478 zind0 = MAX( 0._wp , SIGN( 1._wp , s_i_0 - sm_i_b(ji) ) ) 692 479 ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws 693 zind01 = ( 1.0 - zind0 ) * & 694 MAX( 0.0 , SIGN( 1.0 , s_i_1 - sm_i_b(ji) ) ) 480 zind01 = ( 1._wp - zind0 ) * MAX( 0._wp , SIGN( 1._wp , s_i_1 - sm_i_b(ji) ) ) 695 481 ! if 2.sm_i GE sss_m then zindbal = 1 696 zindbal = MAX( 0.0 , SIGN( 1.0 , 2. * sm_i_b(ji) - & 697 sss_m(zji,zjj) ) ) 698 699 zalpha = zind0 * 1.0 & 700 + zind01 * ( sm_i_b(ji) * dummy_fac0 + & 701 dummy_fac1 ) 702 zalpha = zalpha * ( 1.0 - zindbal ) 703 704 zs_zero(ji,jk) = z_slope_s(ji) * ( jk - 1./2. ) * & 705 ht_i_b(ji) * dummy_fac2 482 zindbal = MAX( 0._wp , SIGN( 1._wp , 2._wp * sm_i_b(ji) - sss_m(zji,zjj) ) ) 483 ! 484 zalpha = ( zind0 + zind01 * ( sm_i_b(ji) * dummy_fac0 + dummy_fac1 ) ) * ( 1.0 - zindbal ) 485 ! 486 zs_zero = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * ht_i_b(ji) * dummy_fac2 706 487 ! weighting the profile 707 s_i_b(ji,jk) = zalpha * zs_zero(ji,jk) + & 708 ( 1.0 - zalpha ) * sm_i_b(ji) 488 s_i_b(ji,jk) = zalpha * zs_zero + ( 1._wp - zalpha ) * sm_i_b(ji) 709 489 END DO ! ji 710 490 END DO ! jk … … 715 495 ! Vertically varying salinity profile, constant in time 716 496 !------------------------------------------------------- 717 ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30) 718 719 IF ( num_sal .EQ. 3 ) THEN720 721 sm_i_b(:) = 2.30722 723 !CDIR NOVERRCHK 724 DO ji = kideb, kiut725 !CDIR NOVERRCHK 726 DO j k = 1, nlay_i727 zargtemp = ( jk - 0.5 ) / nlay_i728 s_i_b(ji,jk) = 1.6 - 1.6*cos(3.14169265*(zargtemp**(0.407/ &729 (0.573+zargtemp))))730 END DO ! jk731 END DO ! ji732 733 ENDIF ! num_sal734 497 498 IF( num_sal == 3 ) THEN ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30) 499 ! 500 sm_i_b(:) = 2.30_wp 501 ! 502 !CDIR NOVERRCHK 503 DO jk = 1, nlay_i 504 zargtemp = ( REAL(jk,wp) - 0.5_wp ) / REAL(nlay_i,wp) 505 zsal = 1.6_wp * ( 1._wp - COS( rpi * zargtemp**(0.407_wp/(0.573_wp+zargtemp)) ) ) 506 DO ji = kideb, kiut 507 s_i_b(ji,jk) = zsal 508 END DO 509 END DO 510 ! 511 ENDIF 512 ! 513 IF( wrk_not_released(1, 1) ) CALL ctl_stop( 'lim_var_salprof1d : failed to release workspace arrays' ) 514 ! 735 515 END SUBROUTINE lim_var_salprof1d 736 516 737 !===============================================================================738 739 517 #else 740 !!====================================================================== 741 !! *** MODULE limvar *** 742 !! no sea ice model 743 !!====================================================================== 518 !!---------------------------------------------------------------------- 519 !! Default option Dummy module NO LIM3 sea-ice model 520 !!---------------------------------------------------------------------- 744 521 CONTAINS 745 522 SUBROUTINE lim_var_agg ! Empty routines … … 755 532 SUBROUTINE lim_var_salprof1d ! Emtpy routines 756 533 END SUBROUTINE lim_var_salprof1d 757 758 534 #endif 535 536 !!====================================================================== 759 537 END MODULE limvar
Note: See TracChangeset
for help on using the changeset viewer.
