- Timestamp:
- 2020-06-24T14:38:26+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r12527_Gurvan_ShallowWater/src/ICE/iceistate.F90
r12489 r13151 18 18 USE oce ! dynamics and tracers variables 19 19 USE dom_oce ! ocean domain 20 USE sbc_oce , ONLY : sst_m, sss_m, ln_ice_embd 20 USE sbc_oce , ONLY : sst_m, sss_m, ln_ice_embd 21 21 USE sbc_ice , ONLY : tn_ice, snwice_mass, snwice_mass_b 22 22 USE eosbn2 ! equation of state … … 60 60 INTEGER , PARAMETER :: jp_hpd = 9 ! index of pnd depth (m) 61 61 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: si ! structure of input fields (file informations, fields read) 62 ! 62 63 63 !! * Substitutions 64 64 # include "do_loop_substitute.h90" … … 77 77 !! 78 78 !! ** Method : This routine will put some ice where ocean 79 !! is at the freezing point, then fill in ice 80 !! state variables using prescribed initial 81 !! values in the namelist 79 !! is at the freezing point, then fill in ice 80 !! state variables using prescribed initial 81 !! values in the namelist 82 82 !! 83 83 !! ** Steps : 1) Set initial surface and basal temperatures … … 91 91 !! where there is no ice 92 92 !!-------------------------------------------------------------------- 93 INTEGER, INTENT(in) :: kt ! time step 93 INTEGER, INTENT(in) :: kt ! time step 94 94 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 95 95 ! … … 102 102 REAL(wp), DIMENSION(jpi,jpj) :: zt_su_ini, zht_s_ini, zsm_i_ini, ztm_i_ini !data from namelist or nc file 103 103 REAL(wp), DIMENSION(jpi,jpj) :: zapnd_ini, zhpnd_ini !data from namelist or nc file 104 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zti_3d , zts_3d ! temporaryarrays104 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zti_3d , zts_3d !locak arrays 105 105 !! 106 106 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zhi_2d, zhs_2d, zai_2d, zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d … … 117 117 ! basal temperature (considered at freezing point) [Kelvin] 118 118 CALL eos_fzp( sss_m(:,:), t_bo(:,:) ) 119 t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) 119 t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) 120 120 ! 121 121 ! surface temperature and conductivity … … 142 142 e_i (:,:,:,:) = 0._wp 143 143 e_s (:,:,:,:) = 0._wp 144 144 145 145 ! general fields 146 146 a_i (:,:,:) = 0._wp … … 213 213 IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) & 214 214 & si(jp_apd)%fnow(:,:,1) = ( rn_apd_ini_n * zswitch + rn_apd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) & ! rn_apd = pond fraction => rn_apnd * a_i = pond conc. 215 & * si(jp_ati)%fnow(:,:,1) 215 & * si(jp_ati)%fnow(:,:,1) 216 216 ! 217 217 ! pond depth … … 227 227 ! 228 228 ! change the switch for the following 229 WHERE( zat_i_ini(:,:) > 0._wp ) ; zswitch(:,:) = tmask(:,:,1) 229 WHERE( zat_i_ini(:,:) > 0._wp ) ; zswitch(:,:) = tmask(:,:,1) 230 230 ELSEWHERE ; zswitch(:,:) = 0._wp 231 231 END WHERE … … 234 234 ! !---------------! 235 235 ! no ice if (sst - Tfreez) >= thresold 236 WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst ) ; zswitch(:,:) = 0._wp 236 WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst ) ; zswitch(:,:) = 0._wp 237 237 ELSEWHERE ; zswitch(:,:) = tmask(:,:,1) 238 238 END WHERE … … 247 247 zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:) 248 248 ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:) 249 zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 249 zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 250 250 zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:) 251 251 ELSEWHERE … … 268 268 zhpnd_ini(:,:) = 0._wp 269 269 ENDIF 270 270 271 271 !-------------! 272 272 ! fill fields ! … … 295 295 ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), & 296 296 & zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) ) 297 297 298 298 ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl) 299 299 CALL ice_var_itd( h_i_1d(1:npti) , h_s_1d(1:npti) , at_i_1d(1:npti), & … … 341 341 DO jl = 1, jpl 342 342 DO_3D_11_11( 1, nlay_i ) 343 t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl) 343 t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl) 344 344 ztmelts = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K 345 345 e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * & … … 357 357 END WHERE 358 358 v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 359 359 360 360 ! specific temperatures for coupled runs 361 361 tn_ice(:,:,:) = t_su(:,:,:) … … 377 377 ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rho0 378 378 ! 379 IF( .NOT.ln_linssh ) THEN 380 ! 381 WHERE( ht_0(:,:) > 0 ) ; z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 382 ELSEWHERE ; z2d(:,:) = 1._wp ; END WHERE 383 ! 384 DO jk = 1,jpkm1 ! adjust initial vertical scale factors 385 e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 386 e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 387 e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 388 END DO 389 ! 390 ! Reconstruction of all vertical scale factors at now and before time-steps 391 ! ========================================================================= 392 ! Horizontal scale factor interpolations 393 ! -------------------------------------- 394 CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 395 CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 396 CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 397 CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 398 CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 399 ! Vertical scale factor interpolations 400 ! ------------------------------------ 401 CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W' ) 402 CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 403 CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 404 CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 405 CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 406 ! t- and w- points depth 407 ! ---------------------- 408 !!gm not sure of that.... 409 gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 410 gdepw(:,:,1,Kmm) = 0.0_wp 411 gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 412 DO jk = 2, jpk 413 gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk ,Kmm) 414 gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 415 gde3w(:,:,jk) = gdept(:,:,jk ,Kmm) - ssh (:,:,Kmm) 416 END DO 417 ENDIF 379 IF( .NOT.ln_linssh ) CALL dom_vvl_zgr( Kbb, Kmm, Kaa ) ! interpolation scale factor, depth and water column 380 ! !!st 381 ! IF( .NOT.ln_linssh ) THEN 382 ! ! 383 ! WHERE( ht_0(:,:) > 0 ) ; z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 384 ! ELSEWHERE ; z2d(:,:) = 1._wp ; END WHERE 385 ! ! 386 ! DO jk = 1,jpkm1 ! adjust initial vertical scale factors 387 ! e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 388 ! e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 389 ! e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 390 ! END DO 391 ! ! 392 ! ! Reconstruction of all vertical scale factors at now and before time-steps 393 ! ! ========================================================================= 394 ! ! Horizontal scale factor interpolations 395 ! ! -------------------------------------- 396 ! CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 397 ! CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 398 ! CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 399 ! CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 400 ! CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 401 ! ! Vertical scale factor interpolations 402 ! ! ------------------------------------ 403 ! CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W' ) 404 ! CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 405 ! CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 406 ! CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 407 ! CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 408 ! ! t- and w- points depth 409 ! ! ---------------------- 410 ! !!gm not sure of that.... 411 ! gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 412 ! gdepw(:,:,1,Kmm) = 0.0_wp 413 ! gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 414 ! DO jk = 2, jpk 415 ! gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk ,Kmm) 416 ! gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 417 ! gde3w(:,:,jk) = gdept(:,:,jk ,Kmm) - ssh (:,:,Kmm) 418 ! END DO 419 ! ENDIF 418 420 ENDIF 419 421 420 422 !------------------------------------ 421 423 ! 4) store fields at before time-step … … 432 434 v_ice_b(:,:) = v_ice(:,:) 433 435 ! total concentration is needed for Lupkes parameterizations 434 at_i_b (:,:) = at_i (:,:) 436 at_i_b (:,:) = at_i (:,:) 435 437 436 438 !!clem: output of initial state should be written here but it is impossible because 437 439 !! the ocean and ice are in the same file 438 !! CALL dia_wri_state( 'output.init' )440 !! CALL dia_wri_state( Kmm, 'output.init' ) 439 441 ! 440 442 END SUBROUTINE ice_istate … … 444 446 !!------------------------------------------------------------------- 445 447 !! *** ROUTINE ice_istate_init *** 446 !! 447 !! ** Purpose : Definition of initial state of the ice 448 !! 449 !! ** Method : Read the namini namelist and check the parameter 448 !! 449 !! ** Purpose : Definition of initial state of the ice 450 !! 451 !! ** Method : Read the namini namelist and check the parameter 450 452 !! values called at the first timestep (nit000) 451 453 !! … … 453 455 !! 454 456 !!----------------------------------------------------------------------------- 455 INTEGER :: ios ! Local integer output status for namelist read456 INTEGER :: ifpr, ierror 457 INTEGER :: ios, ifpr, ierror ! Local integers 458 457 459 ! 458 460 CHARACTER(len=256) :: cn_dir ! Root directory for location of ice files … … 488 490 WRITE(numout,*) ' max ocean temp. above Tfreeze with initial ice rn_thres_sst = ', rn_thres_sst 489 491 IF( ln_iceini .AND. .NOT.ln_iceini_file ) THEN 490 WRITE(numout,*) ' initial snw thickness in the north-south rn_hts_ini = ', rn_hts_ini_n,rn_hts_ini_s 492 WRITE(numout,*) ' initial snw thickness in the north-south rn_hts_ini = ', rn_hts_ini_n,rn_hts_ini_s 491 493 WRITE(numout,*) ' initial ice thickness in the north-south rn_hti_ini = ', rn_hti_ini_n,rn_hti_ini_s 492 494 WRITE(numout,*) ' initial ice concentr in the north-south rn_ati_ini = ', rn_ati_ini_n,rn_ati_ini_s
Note: See TracChangeset
for help on using the changeset viewer.
