- Timestamp:
- 2016-05-09T16:42:28+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r6399 r6515 140 140 !----------------------------------------------------------------- 141 141 SELECT CASE( kblk ) 142 CASE( jp_clio ) ;CALL blk_ice_clio_tau ! CLIO bulk formulation143 CASE( jp_core ) ;CALL blk_ice_core_tau ! CORE bulk formulation144 CASE( jp_purecpl ) ;CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) ! Coupled formulation142 CASE( jp_clio ) ; CALL blk_ice_clio_tau ! CLIO bulk formulation 143 CASE( jp_core ) ; CALL blk_ice_core_tau ! CORE bulk formulation 144 CASE( jp_purecpl ) ; CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) ! Coupled formulation 145 145 END SELECT 146 146 147 IF( ln_mixcpl) THEN ! Case of a mixed Bulk/Coupled formulation147 IF( ln_mixcpl) THEN ! Case of a mixed Bulk/Coupled formulation 148 148 CALL wrk_alloc( jpi,jpj , zutau_ice, zvtau_ice) 149 CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice )149 CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice ) 150 150 utau_ice(:,:) = utau_ice(:,:) * xcplmask(:,:,0) + zutau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) 151 151 vtau_ice(:,:) = vtau_ice(:,:) * xcplmask(:,:,0) + zvtau_ice(:,:) * ( 1. - xcplmask(:,:,0) ) … … 158 158 numit = numit + nn_fsbc ! Ice model time step 159 159 ! 160 CALL sbc_lim_bef ! Store previous ice values 161 CALL sbc_lim_diag0 ! set diag of mass, heat and salt fluxes to 0 162 CALL lim_rst_opn( kt ) ! Open Ice restart file 163 ! 164 IF( .NOT. lk_c1d ) THEN 160 CALL sbc_lim_bef ! Store previous ice values 161 CALL sbc_lim_diag0 ! set diag of mass, heat and salt fluxes to 0 162 CALL lim_rst_opn( kt ) ! Open Ice restart file 163 ! 164 ! --- zap this if no ice dynamics --- ! 165 IF( .NOT. lk_c1d .AND. ln_limdyn ) THEN 165 166 ! 166 CALL lim_dyn( kt ) ! Ice dynamics ( rheology/dynamics ) 167 ! 168 CALL lim_trp( kt ) ! Ice transport ( Advection/diffusion ) 169 ! 170 IF( nn_monocat /= 2 ) CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) 171 ! 172 #if defined key_bdy 173 CALL bdy_ice_lim( kt ) ! bdy ice thermo 174 IF( ln_icectl ) CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 175 #endif 176 ! 177 CALL lim_update1( kt ) ! Corrections 167 IF( nn_limdyn /= 0 ) THEN ! -- Ice dynamics 168 CALL lim_dyn( kt ) ! rheology 169 ELSE 170 u_ice(:,:) = rn_uice * umask(:,:,1) ! or prescribed velocity 171 v_ice(:,:) = rn_vice * vmask(:,:,1) 172 ENDIF 173 CALL lim_trp( kt ) ! -- Ice transport (Advection/diffusion) 174 IF( nn_limdyn == 2 .AND. nn_monocat /= 2 ) & ! -- Mechanical redistribution (ridging/rafting) 175 & CALL lim_itd_me 176 IF( nn_limdyn == 2 ) CALL lim_update1( kt ) ! -- Corrections 178 177 ! 179 178 ENDIF 180 179 ! --- 180 #if defined key_bdy 181 IF( ln_limthd ) CALL bdy_ice_lim( kt ) ! bdy ice thermo 182 IF( ln_icectl ) CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 183 #endif 184 181 185 ! previous lead fraction and ice volume for flux calculations 182 CALL sbc_lim_bef 183 CALL lim_var_glo2eqv ! ht_i and ht_s for ice albedo calculation 184 CALL lim_var_agg(1) ! at_i for coupling (via pfrld) 186 CALL sbc_lim_bef 187 CALL lim_var_glo2eqv ! ht_i and ht_s for ice albedo calculation 188 CALL lim_var_agg(1) ! at_i for coupling (via pfrld) 189 ! 185 190 pfrld(:,:) = 1._wp - at_i(:,:) 186 191 phicif(:,:) = vt_i(:,:) … … 197 202 !---------------------------------------------------------------------------------------- 198 203 CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs ) 199 CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos200 204 205 CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 201 206 SELECT CASE( kblk ) 202 CASE( jp_clio ) ! CLIO bulk formulation 203 ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 204 ! (alb_ice) is computed within the bulk routine 205 CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, alb_ice ) 206 IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 207 IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 208 CASE( jp_core ) ! CORE bulk formulation 209 ! albedo depends on cloud fraction because of non-linear spectral effects 210 alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 211 CALL blk_ice_core_flx( t_su, alb_ice ) 212 IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 213 IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 214 CASE ( jp_purecpl ) 215 ! albedo depends on cloud fraction because of non-linear spectral effects 216 alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 217 CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 218 IF( nn_limflx == 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 207 208 CASE( jp_clio ) ! CLIO bulk formulation 209 ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 210 ! (alb_ice) is computed within the bulk routine 211 CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, alb_ice ) 212 IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 213 IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 214 215 CASE( jp_core ) ! CORE bulk formulation 216 ! albedo depends on cloud fraction because of non-linear spectral effects 217 alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 218 CALL blk_ice_core_flx( t_su, alb_ice ) 219 IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 220 IF( nn_limflx /= 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 221 222 CASE ( jp_purecpl ) ! Coupled formulation 223 ! albedo depends on cloud fraction because of non-linear spectral effects 224 alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 225 CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=alb_ice, psst=sst_m, pist=t_su ) 226 IF( nn_limflx == 2 ) CALL ice_lim_flx( t_su, alb_ice, qns_ice, qsr_ice, dqns_ice, evap_ice, devap_ice, nn_limflx ) 227 219 228 END SELECT 229 220 230 CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs ) 221 231 … … 223 233 ! --- ice thermodynamics --- ! 224 234 !----------------------------! 225 CALL lim_thd( kt ) ! Ice thermodynamics 226 ! 227 CALL lim_update2( kt ) ! Corrections 228 ! 229 CALL lim_sbc_flx( kt ) ! Update surface ocean mass, heat and salt fluxes 230 ! 231 IF(ln_limdiaout) CALL lim_diahsb ! Diagnostics and outputs 232 ! 233 CALL lim_wri( 1 ) ! Ice outputs 235 ! --- zap this if no ice thermo --- ! 236 IF( ln_limthd ) CALL lim_thd( kt ) ! -- Ice thermodynamics 237 IF( ln_limthd ) CALL lim_update2( kt ) ! -- Corrections 238 ! --- 239 CALL lim_var_glo2eqv ! necessary calls (at least for coupling) 240 CALL lim_var_agg( 2 ) ! necessary calls (at least for coupling) 241 ! 242 CALL lim_sbc_flx( kt ) ! -- Update surface ocean mass, heat and salt fluxes 243 ! 244 IF(ln_limdiaout) CALL lim_diahsb ! -- Diagnostics and outputs 245 ! 246 CALL lim_wri( 1 ) ! -- Ice outputs 234 247 ! 235 248 IF( kt == nit000 .AND. ln_rstart ) & 236 & CALL iom_close( numrir ) ! close input ice restart file237 ! 238 IF( lrst_ice ) CALL lim_rst_write( kt ) !Ice restart file239 ! 240 IF( ln_icectl ) CALL lim_ctl( kt ) ! alerts in case of model crash249 & CALL iom_close( numrir ) ! close input ice restart file 250 ! 251 IF( lrst_ice ) CALL lim_rst_write( kt ) ! -- Ice restart file 252 ! 253 IF( ln_icectl ) CALL lim_ctl( kt ) ! alerts in case of model crash 241 254 ! 242 255 ENDIF ! End sea-ice time step only … … 246 259 !-------------------------! 247 260 ! Update surface ocean stresses (only in ice-dynamic case) otherwise the atm.-ocean stresses are used everywhere 248 IF( ln_limdyn ) CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) ) ! using before instantaneous surf. currents 261 ! using before instantaneous surf. currents 262 IF( ln_limdyn ) CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) ) 249 263 !!gm remark, the ocean-ice stress is not saved in ice diag call above ..... find a solution!!! 250 264 ! … … 264 278 !!---------------------------------------------------------------------- 265 279 IF(lwp) WRITE(numout,*) 266 IF(lwp) WRITE(numout,*) 'sbc_ ice_lim: update ocean surface boudary condition'280 IF(lwp) WRITE(numout,*) 'sbc_lim_init : update ocean surface boudary condition' 267 281 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ via Louvain la Neuve Ice Model (LIM-3) time stepping' 268 282 ! … … 279 293 ierr = ierr + sbc_ice_alloc () ! surface forcing 280 294 ierr = ierr + thd_ice_alloc () ! thermodynamics 281 ierr = ierr + lim_itd_me_alloc () ! ice thickness distribution - mechanics295 IF( ln_limdyn ) ierr = ierr + lim_itd_me_alloc () ! ice thickness distribution - mechanics 282 296 ! 283 297 IF( lk_mpp ) CALL mpp_sum( ierr ) … … 300 314 CALL lim_msh ! ice mesh initialization 301 315 ! 302 CALL lim_itd_me_init ! ice thickness distribution initialization for mecanical deformation316 IF( ln_limdyn ) CALL lim_itd_me_init ! ice thickness distribution initialization for mecanical deformation 303 317 ! ! Initial sea-ice state 304 318 IF( .NOT. ln_rstart ) THEN ! start from rest: sea-ice deduced from sst … … 347 361 !!------------------------------------------------------------------- 348 362 INTEGER :: ios ! Local integer output status for namelist read 349 NAMELIST/namicerun/ jpl, nlay_i, nlay_s, cn_icerst_in, cn_icerst_indir, cn_icerst_out, cn_icerst_outdir, & 350 & ln_limdyn, rn_amax_n, rn_amax_s, ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt 363 NAMELIST/namicerun/ jpl, nlay_i, nlay_s, rn_amax_n, rn_amax_s, cn_icerst_in, cn_icerst_indir, & 364 & cn_icerst_out, cn_icerst_outdir, ln_limthd, ln_limdyn, nn_limdyn, rn_uice, rn_vice 365 NAMELIST/namicediag/ ln_limdiahsb, ln_limdiaout, ln_icectl, iiceprt, jiceprt 351 366 !!------------------------------------------------------------------- 352 367 ! … … 360 375 IF(lwm) WRITE ( numoni, namicerun ) 361 376 ! 377 REWIND( numnam_ice_ref ) ! Namelist namicediag in reference namelist : Parameters for ice 378 READ ( numnam_ice_ref, namicediag, IOSTAT = ios, ERR = 903) 379 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicediag in reference namelist', lwp ) 380 381 REWIND( numnam_ice_cfg ) ! Namelist namicediag in configuration namelist : Parameters for ice 382 READ ( numnam_ice_cfg, namicediag, IOSTAT = ios, ERR = 904 ) 383 904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicediag in configuration namelist', lwp ) 384 IF(lwm) WRITE ( numoni, namicediag ) 362 385 ! 363 386 IF(lwp) THEN ! control print 364 387 WRITE(numout,*) 365 WRITE(numout,*) 'ice_run : ice share parameters for dynamics/advection/thermo of sea-ice'388 WRITE(numout,*) 'ice_run : ice share~d parameters for dynamics/advection/thermo of sea-ice' 366 389 WRITE(numout,*) ' ~~~~~~' 367 390 WRITE(numout,*) ' number of ice categories = ', jpl 368 391 WRITE(numout,*) ' number of ice layers = ', nlay_i 369 392 WRITE(numout,*) ' number of snow layers = ', nlay_s 370 WRITE(numout,*) ' switch for ice dynamics (1) or not (0) ln_limdyn = ', ln_limdyn371 393 WRITE(numout,*) ' maximum ice concentration for NH = ', rn_amax_n 372 394 WRITE(numout,*) ' maximum ice concentration for SH = ', rn_amax_s 373 WRITE(numout,*) ' Diagnose heat/salt budget or not ln_limdiahsb = ', ln_limdiahsb 374 WRITE(numout,*) ' Output heat/salt budget or not ln_limdiaout = ', ln_limdiaout 395 WRITE(numout,*) ' Ice thermodynamics (T) or not (F) ln_limthd = ', ln_limthd 396 WRITE(numout,*) ' Ice dynamics (T) or not (F) ln_limdyn = ', ln_limdyn 397 WRITE(numout,*) ' (ln_limdyn=T) Ice dynamics switch nn_limdyn = ', nn_limdyn 398 WRITE(numout,*) ' 2: total' 399 WRITE(numout,*) ' 1: advection only (no diffusion, no ridging/rafting)' 400 WRITE(numout,*) ' 0: advection only (as 1 + prescribed velocity, bypass rheology)' 401 WRITE(numout,*) ' (ln_limdyn=T) prescribed u-vel (case nn_limdyn=0) = ', rn_uice 402 WRITE(numout,*) ' (ln_limdyn=T) prescribed v-vel (case nn_limdyn=0) = ', rn_vice 403 WRITE(numout,*) 404 WRITE(numout,*) '...and ice diagnostics' 405 WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~' 406 WRITE(numout,*) ' Diagnose heat/mass/salt budget or not ln_limdiahsb = ', ln_limdiahsb 407 WRITE(numout,*) ' Output heat/mass/salt budget or not ln_limdiaout = ', ln_limdiaout 375 408 WRITE(numout,*) ' control prints in ocean.out for (i,j)=(iiceprt,jiceprt) = ', ln_icectl 376 409 WRITE(numout,*) ' i-index for control prints (ln_icectl=true) = ', iiceprt … … 410 443 ! 411 444 REWIND( numnam_ice_ref ) ! Namelist namiceitd in reference namelist : Parameters for ice 412 READ ( numnam_ice_ref, namiceitd, IOSTAT = ios, ERR = 90 3)413 90 3IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in reference namelist', lwp )445 READ ( numnam_ice_ref, namiceitd, IOSTAT = ios, ERR = 905) 446 905 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in reference namelist', lwp ) 414 447 415 448 REWIND( numnam_ice_cfg ) ! Namelist namiceitd in configuration namelist : Parameters for ice 416 READ ( numnam_ice_cfg, namiceitd, IOSTAT = ios, ERR = 90 4)417 90 4IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in configuration namelist', lwp )449 READ ( numnam_ice_cfg, namiceitd, IOSTAT = ios, ERR = 906 ) 450 906 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namiceitd in configuration namelist', lwp ) 418 451 IF(lwm) WRITE ( numoni, namiceitd ) 419 !420 452 ! 421 453 IF(lwp) THEN ! control print 422 454 WRITE(numout,*) 423 WRITE(numout,*) ' ice_itd : ice cat distribution'424 WRITE(numout,*) ' 455 WRITE(numout,*) 'lim_itd_init : Initialization of ice cat distribution ' 456 WRITE(numout,*) '~~~~~~~~~~~~' 425 457 WRITE(numout,*) ' shape of ice categories distribution nn_catbnd = ', nn_catbnd 426 458 WRITE(numout,*) ' mean ice thickness in the domain (only active if nn_catbnd=2) rn_himean = ', rn_himean … … 430 462 !- Thickness categories boundaries 431 463 !---------------------------------- 432 IF(lwp) WRITE(numout,*)433 IF(lwp) WRITE(numout,*) 'lim_itd_init : Initialization of ice cat distribution '434 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'435 436 464 hi_max(:) = 0._wp 437 465 … … 581 609 !!---------------------------------------------------------------------- 582 610 sfx (:,:) = 0._wp ; 583 sfx_bri(:,:) = 0._wp ; 611 sfx_bri(:,:) = 0._wp ; sfx_lam(:,:) = 0._wp 584 612 sfx_sni(:,:) = 0._wp ; sfx_opw(:,:) = 0._wp 585 613 sfx_bog(:,:) = 0._wp ; sfx_dyn(:,:) = 0._wp … … 592 620 wfx_bom(:,:) = 0._wp ; wfx_sum(:,:) = 0._wp 593 621 wfx_res(:,:) = 0._wp ; wfx_sub(:,:) = 0._wp 594 wfx_spr(:,:) = 0._wp ; 622 wfx_spr(:,:) = 0._wp ; wfx_lam(:,:) = 0._wp 595 623 596 624 hfx_thd(:,:) = 0._wp ;
Note: See TracChangeset
for help on using the changeset viewer.