Changeset 2528 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90
- Timestamp:
- 2010-12-27T18:33:53+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90
r2477 r2528 1 1 MODULE limthd_sal 2 !!----------------------------------------------------------------------3 !! 'key_lim3' LIM3 sea-ice model4 !!----------------------------------------------------------------------5 2 !!====================================================================== 6 3 !! *** MODULE limthd_sal *** 7 !! computation salinity variations in 8 !! the ice 4 !! LIM-3 sea-ice : computation of salinity variations in the ice 9 5 !!====================================================================== 6 !! History : - ! 2003-05 (M. Vancoppenolle) UCL-ASTR first coding for LIM3-1D 7 !! 3.0 ! 2005-12 (M. Vancoppenolle) adapted to the 3-D version 8 !!--------------------------------------------------------------------- 10 9 #if defined key_lim3 11 10 !!---------------------------------------------------------------------- 11 !! 'key_lim3' LIM-3 sea-ice model 12 !!---------------------------------------------------------------------- 12 13 !! lim_thd_sal : salinity variations in the ice 13 !! * Modules used14 !!---------------------------------------------------------------------- 14 15 USE par_oce ! ocean parameters 15 16 USE phycst ! physical constants (ocean directory) 16 17 USE sbc_oce ! Surface boundary condition: ocean fields 17 USE thd_ice18 USE ice19 USE in_out_manager20 USE limvar 21 USE par_ice18 USE ice ! LIM: sea-ice variables 19 USE par_ice ! LIM: sea-ice parameters 20 USE thd_ice ! LIM: sea-ice thermodynamics 21 USE limvar ! LIM: sea-ice variables 22 USE in_out_manager ! I/O manager 22 23 23 24 IMPLICIT NONE 24 25 PRIVATE 25 26 26 !! * Routine accessibility 27 PUBLIC lim_thd_sal ! called by lim_thd 28 PUBLIC lim_thd_sal_init ! called by lim_thd 29 30 !! * Module variables 31 REAL(wp) :: & ! constant values 32 epsi20 = 1e-20 , & 33 epsi13 = 1e-13 , & 34 zzero = 0.e0 , & 35 zone = 1.e0 36 37 !!---------------------------------------------------------------------- 38 !! LIM 3.0, UCL-ASTR-LOCEAN-IPSL (2008) 27 PUBLIC lim_thd_sal ! called by limthd module 28 PUBLIC lim_thd_sal_init ! called by iceini module 29 30 !!---------------------------------------------------------------------- 31 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 39 32 !! $Id$ 40 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 41 !!---------------------------------------------------------------------- 42 33 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 34 !!---------------------------------------------------------------------- 43 35 CONTAINS 44 36 45 SUBROUTINE lim_thd_sal( kideb,kiut)37 SUBROUTINE lim_thd_sal( kideb, kiut ) 46 38 !!------------------------------------------------------------------- 47 !! *** ROUTINE lim_thd_sal *** 48 !! ** Purpose :49 !! This routinecomputes new salinities in the ice39 !! *** ROUTINE lim_thd_sal *** 40 !! 41 !! ** Purpose : computes new salinities in the ice 50 42 !! 51 43 !! ** Method : 4 possibilities … … 54 46 !! -> num_sal = 3 -> S = S(z) [multiyear ice] 55 47 !! -> num_sal = 4 -> S = S(h) [Cox and Weeks 74] 56 !!57 !! ** Steps58 !!59 !! ** Arguments60 !!61 !! ** Inputs / Outputs62 !!63 !! ** External64 !!65 !! ** References66 !!67 !! ** History :68 !!69 !! "Je ne suis pour l'instant qu'a 80% de ma condition, mais c'est70 !! les 30% qui restent qui seront les plus difficiles"71 !! E. Mpenza72 !!73 !!-------------------------------------------------------------------74 !! History :75 !! ori ! 03-05 M. Vancoppenolle UCL-ASTR first coding for LIM-1D76 !! 3.0 ! 05-12 Routine rewritten for the 3-D version77 48 !!--------------------------------------------------------------------- 78 !! 79 !! * Local variables 80 INTEGER, INTENT(in) :: & 81 kideb, kiut !: thickness category index 82 83 INTEGER :: & 84 ji, jk , & !: geographic and layer index 85 zji, zjj 86 87 REAL(wp) :: & 88 zsold, & !: old salinity 89 zeps=1.0e-06 , & !: very small 90 iflush , & !: flushing (1) or not (0) 91 iaccrbo , & !: bottom accretion (1) or not (0) 92 igravdr , & !: gravity drainage or not 93 isnowic , & !: gravity drainage or not 94 i_ice_switch , & !: ice thickness above a certain treshold or not 95 ztmelts , & !: freezing point of sea ice 96 zaaa , & !: dummy factor 97 zbbb , & !: dummy factor 98 zccc , & !: dummy factor 99 zdiscrim !: dummy factor 100 101 REAL(wp), DIMENSION(jpij) :: & 102 ze_init , & !initial total enthalpy 103 zhiold , & 104 zsiold 105 49 INTEGER, INTENT(in) :: kideb, kiut ! thickness category index 50 ! 51 INTEGER :: ji, jk ! dummy loop indices 52 INTEGER :: zji, zjj ! local integers 53 REAL(wp) :: zsold, zeps, iflush, iaccrbo, igravdr, isnowic, i_ice_switch, ztmelts ! local scalars 54 REAL(wp) :: zaaa, zbbb, zccc, zdiscrim ! local scalars 55 REAL(wp), DIMENSION(jpij) :: ze_init, zhiold, zsiold ! 1D workspace 106 56 !!--------------------------------------------------------------------- 107 57 58 zeps=1.0e-06_wp 59 108 60 !------------------------------------------------------------------------------| 109 61 ! 1) Constant salinity, constant in time | 110 62 !------------------------------------------------------------------------------| 111 63 112 IF (num_sal.eq.1) THEN 113 114 ! WRITE(numout,*) 115 ! WRITE(numout,*) 'lim_thd_sal : Ice salinity computation module ', & 116 ! num_sal 117 ! WRITE(numout,*) '~~~~~~~~~~~~' 118 64 IF( num_sal == 1 ) THEN 119 65 DO jk = 1, nlay_i 120 66 DO ji = kideb, kiut … … 122 68 END DO ! ji 123 69 END DO ! jk 124 70 ! 125 71 DO ji = kideb, kiut 126 72 sm_i_b(ji) = bulk_sal 127 73 END DO ! ji 128 129 ENDIF ! num_sal .EQ. 174 ! 75 ENDIF 130 76 131 77 !------------------------------------------------------------------------------| … … 174 120 ! isnowic : 1 if snow ice formation 175 121 i_ice_switch = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_i_b(ji) + 1.0e-2 ) ) 176 isnowic = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - dh_snowice(ji) ) ) * & 177 i_ice_switch 122 isnowic = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - dh_snowice(ji) ) ) * i_ice_switch 178 123 179 124 !--------------------- … … 182 127 183 128 ! drainage by gravity drainage 184 dsm_i_gd_1d(ji) = - igravdr * & 185 MAX( sm_i_b(ji) - sal_G , 0.0 ) / & 186 time_G * rdt_ice 129 dsm_i_gd_1d(ji) = - igravdr * MAX( sm_i_b(ji) - sal_G , 0._wp ) / time_G * rdt_ice 187 130 188 131 ! drainage by flushing 189 dsm_i_fl_1d(ji) = - iflush * & 190 MAX( sm_i_b(ji) - sal_F , 0.0 ) / & 191 time_F * rdt_ice 132 dsm_i_fl_1d(ji) = - iflush * MAX( sm_i_b(ji) - sal_F , 0._wp ) / time_F * rdt_ice 192 133 193 134 !----------------- … … 199 140 ! to conserve energy 200 141 zsiold(ji) = sm_i_b(ji) 201 sm_i_b(ji) = sm_i_b(ji) & 202 + dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji) ! & 142 sm_i_b(ji) = sm_i_b(ji) + dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji) 203 143 204 144 ! if no ice, salinity eq 0.1 205 i_ice_switch = 1.0 - MAX ( 0.0, SIGN (1.0 , - ht_i_b(ji) ) ) 206 sm_i_b(ji) = i_ice_switch*sm_i_b(ji) + s_i_min * ( 1.0 - & 207 i_ice_switch ) 145 i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 146 sm_i_b(ji) = i_ice_switch*sm_i_b(ji) + s_i_min * ( 1._wp - i_ice_switch ) 208 147 END DO ! ji 209 148 210 149 ! Salinity profile 211 CALL lim_var_salprof1d( kideb,kiut)150 CALL lim_var_salprof1d( kideb, kiut ) 212 151 213 152 !---------------------------- … … 217 156 DO ji = kideb, kiut 218 157 ! iflush : 1 if summer 219 iflush = MAX( 0.0 , SIGN ( 1.0, t_su_b(ji) - rtt ) )158 iflush = MAX( 0._wp , SIGN ( 1._wp , t_su_b(ji) - rtt ) ) 220 159 ! igravdr : 1 if t_su lt t_bo 221 igravdr = MAX( 0.0 , SIGN ( 1.0, t_bo_b(ji) - t_su_b(ji) ) )160 igravdr = MAX( 0._wp , SIGN ( 1._wp , t_bo_b(ji) - t_su_b(ji) ) ) 222 161 ! iaccrbo : 1 if bottom accretion 223 iaccrbo = MAX( 0.0 , SIGN ( 1.0, dh_i_bott(ji) ) )224 225 fhbri_1d(ji) = 0. 0162 iaccrbo = MAX( 0._wp , SIGN ( 1._wp , dh_i_bott(ji) ) ) 163 ! 164 fhbri_1d(ji) = 0._wp 226 165 END DO ! ji 227 166 … … 230 169 !---------------------------- 231 170 DO ji = kideb, kiut 232 i_ice_switch = 1.0 - MAX ( 0.0, SIGN (1.0 , - ht_i_b(ji) ) ) 233 fsbri_1d(ji) = fsbri_1d(ji) - & 234 i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) * & 235 ( MAX(dsm_i_gd_1d(ji) + dsm_i_fl_1d(ji), & 236 sm_i_b(ji) - zsiold(ji) ) ) / rdt_ice 237 IF ( num_sal .EQ. 4 ) fsbri_1d(ji) = 0.0 238 171 i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 172 fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) & 173 & * ( MAX(dsm_i_gd_1d(ji) + dsm_i_fl_1d(ji), sm_i_b(ji) - zsiold(ji) ) ) / rdt_ice 174 IF( num_sal == 4 ) fsbri_1d(ji) = 0._wp 239 175 END DO ! ji 240 176 … … 244 180 !-------------------- 245 181 DO jk = 1, nlay_i 246 247 182 DO ji = kideb, kiut 248 249 183 ztmelts = -tmut*s_i_b(ji,jk) + rtt 250 184 !Conversion q(S,T) -> T (second order equation) 251 185 zaaa = cpic 252 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + & 253 q_i_b(ji,jk) / rhoic - lfus 186 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk) / rhoic - lfus 254 187 zccc = lfus * ( ztmelts - rtt ) 255 188 zdiscrim = SQRT( MAX(zbbb*zbbb - 4.0*zaaa*zccc,0.0) ) 256 t_i_b(ji,jk) = rtt - ( zbbb + zdiscrim ) / & 257 ( 2.0 *zaaa ) 189 t_i_b(ji,jk) = rtt - ( zbbb + zdiscrim ) / ( 2.0 *zaaa ) 258 190 END DO !ji 259 260 191 END DO !jk 261 192 ! 262 193 ENDIF ! num_sal .EQ. 2 263 194 … … 266 197 !------------------------------------------------------------------------------| 267 198 268 IF 199 IF( num_sal .EQ. 3 ) THEN 269 200 270 201 WRITE(numout,*) … … 320 251 zji = MOD( npb(ji) - 1, jpi ) + 1 321 252 zjj = ( npb(ji) - 1 ) / jpi + 1 322 fseqv_1d(ji) = fseqv_1d(ji) + & 323 ( sss_m(zji,zjj) - bulk_sal ) * & 324 rhoic * a_i_b(ji) * & 325 MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 253 fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - bulk_sal ) & 254 & * rhoic * a_i_b(ji) * MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 326 255 END DO 327 256 ELSE … … 329 258 zji = MOD( npb(ji) - 1, jpi ) + 1 330 259 zjj = ( npb(ji) - 1 ) / jpi + 1 331 fseqv_1d(ji) = fseqv_1d(ji) + & 332 ( sss_m(zji,zjj) - s_i_new(ji) ) * & 333 rhoic * a_i_b(ji) * & 334 MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 260 fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - s_i_new(ji) ) & 261 & * rhoic * a_i_b(ji) * MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 335 262 END DO ! ji 336 263 ENDIF 337 338 !-- End of salinity computations 264 ! 339 265 END SUBROUTINE lim_thd_sal 340 !============================================================================== 266 341 267 342 268 SUBROUTINE lim_thd_sal_init … … 346 272 !! ** Purpose : initialization of ice salinity parameters 347 273 !! 348 !! ** Method : Read the namicesal namelist and check the parameter349 !! values called at the first timestep (nit000)274 !! ** Method : Read the namicesal namelist and check the parameter 275 !! values called at the first timestep (nit000) 350 276 !! 351 277 !! ** input : Namelist namicesal 352 !!353 !! history :354 !! 3.0 ! July 2005 M. Vancoppenolle Original code355 278 !!------------------------------------------------------------------- 356 NAMELIST/namicesal/ num_sal, bulk_sal, sal_G, time_G, sal_F, time_F, &357 s_i_max, s_i_min, s_i_0, s_i_1279 NAMELIST/namicesal/ num_sal, bulk_sal, sal_G, time_G, sal_F, time_F, & 280 & s_i_max, s_i_min, s_i_0, s_i_1 358 281 !!------------------------------------------------------------------- 359 360 ! Read Namelist namicesal361 RE WIND ( numnam_ice)362 READ ( numnam_ice , namicesal )363 IF(lwp) THEN 282 ! 283 REWIND( numnam_ice ) ! Read Namelist namicesal 284 READ ( numnam_ice , namicesal ) 285 ! 286 IF(lwp) THEN ! control print 364 287 WRITE(numout,*) 365 288 WRITE(numout,*) 'lim_thd_sal_init : Ice parameters for salinity ' … … 376 299 WRITE(numout,*) ' 2nd salinity for salinity profile : ', s_i_1 377 300 ENDIF 378 301 ! 379 302 END SUBROUTINE lim_thd_sal_init 380 303 381 304 #else 382 305 !!---------------------------------------------------------------------- 383 !! Default option Empty Module No sea-ice model 384 !!---------------------------------------------------------------------- 385 CONTAINS 386 SUBROUTINE lim_thd_sal ! Empty routine 387 END SUBROUTINE lim_thd_sal 306 !! Default option Dummy Module No LIM-3 sea-ice model 307 !!---------------------------------------------------------------------- 388 308 #endif 389 309 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.