Changeset 11334
- Timestamp:
- 2019-07-24T11:20:42+02:00 (6 years ago)
- Location:
- NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D
- Files:
-
- 1 added
- 8 edited
-
cfgs/SHARED/namelist_ref (modified) (2 diffs)
-
src/ABL/ablmod.F90 (modified) (11 diffs)
-
src/ABL/sbcabl.F90 (modified) (6 diffs)
-
src/ICE/icesbc.F90 (modified) (2 diffs)
-
src/ICE/icevar.F90 (modified) (4 diffs)
-
src/OCE/DIA/diawri.F90 (modified) (13 diffs)
-
src/OCE/SBC/abl.F90 (added)
-
src/OCE/SBC/sbcabl.F90 (modified) (2 diffs)
-
src/OCE/SBC/sbcblk.F90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/cfgs/SHARED/namelist_ref
r11322 r11334 301 301 rn_Rod = 0.15 ! c0 in RMCA17 mixing length formulation (not yet implemented) 302 302 rn_Ric = 0.139 ! Critical Richardson number (to compute PBL height and diffusivities) 303 ln_smth_pblh = . false. ! Smoothing of PBL height with a 2d Hanning filter303 ln_smth_pblh = .true. ! Smoothing of PBL height with a 2d Hanning filter 304 304 / 305 305 !----------------------------------------------------------------------- … … 589 589 !----------------------------------------------------------------------- 590 590 ln_tide = .false. ! Activate tides 591 ln_tide_pot = .true. ! use tidal potential forcing 592 ln_scal_load = .false. ! Use scalar approximation for 593 rn_scal_load = 0.094 ! load potential 594 ln_read_load = .false. ! Or read load potential from file 595 cn_tide_load = 'tide_LOAD_grid_T.nc' ! filename for load potential 596 ! 597 ln_tide_ramp = .false. ! Use linear ramp for tides at startup 598 rdttideramp = 0. ! ramp duration in days 599 clname(1) = 'DUMMY' ! name of constituent - all tidal components must be set in namelist_cfg 591 ln_tide_pot = .true. ! use tidal potential forcing 592 ln_scal_load = .false. ! Use scalar approximation for 593 ln_read_load = .false. ! Or read load potential from file 594 cn_tide_load = 'tide_LOAD_grid_T.nc' ! filename for load potential 595 ln_tide_ramp = .false. ! Use linear ramp for tides at startup 596 rn_scal_load = 0.094 ! load potential 597 rdttideramp = 0. ! ramp duration in days 598 clname = 'FOUR' ! name of constituent - all tidal components must be set in namelist_cfg 600 599 / 601 600 !----------------------------------------------------------------------- -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/ABL/ablmod.F90
r11322 r11334 38 38 !=================================================================================================== 39 39 SUBROUTINE abl_stp( kt, psst, pssu, pssv, pssq, & ! in 40 & pu_dta, pv_dta, pt_dta, pq_dta, & 41 & pslp_dta, pgu_dta, pgv_dta, & 42 & pcd_du, psen, pevp, & ! in/out 43 & pwndm, ptaui, ptauj, ptaum ) 40 & pu_dta, pv_dta, pt_dta, pq_dta, & 41 & pslp_dta, pgu_dta, pgv_dta, & 42 & pcd_du, psen, pevp, & ! in/out 43 & pwndm, ptaui, ptauj, ptaum & 44 #if defined key_si3 45 & , ptm_su,pssu_ice,pssv_ice,pssq_ice,pcd_du_ice & 46 & , psen_ice, pevp_ice, pwndm_ice, pfrac_oce & 47 & , ptaui_ice, ptauj_ice & 48 #endif 49 & ) 44 50 !--------------------------------------------------------------------------------------------------- 45 51 … … 80 86 REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaum ! ||tau|| 81 87 ! 88 #if defined key_si3 89 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: ptm_su ! ice-surface temperature [K] 90 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssu_ice ! ice-surface u (U-point) 91 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssv_ice ! ice-surface v (V-point) 92 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssq_ice ! ice-surface humidity 93 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pcd_du_ice ! Cd x Du over ice (T-point) 94 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: psen_ice ! Ch x Du over ice (T-point) 95 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pevp_ice ! Ce x Du over ice (T-point) 96 REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pwndm_ice ! ||uwnd - uice|| 97 REAL(wp) , INTENT(inout), DIMENSION(:,: ) :: pfrac_oce ! 98 REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaui_ice ! ice-surface taux stress (U-point) 99 REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptauj_ice ! ice-surface tauy stress (V-point) 100 #endif 101 ! 82 102 REAL(wp), DIMENSION(1:jpi,1:jpj ) :: zrhoa, zwnd_i, zwnd_j 83 ! REAL(wp), DIMENSION(1:jpi,1:jpka ) :: zFC84 103 REAL(wp), DIMENSION(1:jpi,2:jpka ) :: zCF 85 ! REAL(wp), DIMENSION(1:jpi, jptq ) :: zBC 86 REAL(wp), DIMENSION(1:jpi,1:jpj,1:jpka) :: z_cft !--FL--to be removed after the test phase 104 REAL(wp), DIMENSION(1:jpi,1:jpj,1:jpka) :: z_cft !--FL--to be removed after the test phase 87 105 ! 88 106 REAL(wp), DIMENSION(1:jpi,1:jpka ) :: z_elem_a … … 91 109 ! 92 110 INTEGER :: ji, jj, jk, jtra, jbak ! dummy loop indices 93 REAL(wp) :: zztmp, zcff, ztemp, zhumi, zcff1 94 REAL(wp) :: zcff2, zfcor, zmsk, zsig, zcffu, zcffv 111 REAL(wp) :: zztmp, zcff, ztemp, zhumi, zcff1, zztmp1, zztmp2 112 REAL(wp) :: zcff2, zfcor, zmsk, zsig, zcffu, zcffv, zzice,zzoce 95 113 ! 96 114 !!--------------------------------------------------------------------- … … 108 126 DO jj = 1,jpj 109 127 DO ji = 1,jpi 110 ustar2(ji,jj) = pCd_du(ji,jj) * pwndm(ji,jj) 128 zzoce = pCd_du (ji,jj) * pwndm (ji,jj) 129 #if defined key_si3 130 zzice = pCd_du_ice(ji,jj) * pwndm_ice(ji,jj) 131 ustar2(ji,jj) = zzoce * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * zzice 132 #else 133 ustar2(ji,jj) = zzoce 134 #endif 111 135 END DO 112 136 END DO … … 154 178 IF(jtra == jp_ta) THEN 155 179 DO ji = 1,jpi ! surface boundary condition for temperature 156 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * psen(ji, jj) 157 tq_abl ( ji, jj, 2, nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) & 158 & + rdt * psen(ji, jj) * ( psst(ji, jj) + rt0 ) 159 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) 180 zztmp1 = psen(ji, jj) 181 zztmp2 = psen(ji, jj) * ( psst(ji, jj) + rt0 ) 182 #if defined key_si3 183 zztmp1 = zztmp1 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * psen_ice(ji,jj) 184 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * psen_ice(ji,jj) * ptm_su(ji,jj) 185 #endif 186 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * zztmp1 187 tq_abl ( ji, jj, 2, nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) + rdt * zztmp2 188 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) 160 189 END DO 161 190 ELSE 162 191 DO ji = 1,jpi ! surface boundary condition for humidity 163 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * pevp(ji, jj) 164 tq_abl ( ji, jj, 2, nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) & 165 & + rdt * pevp(ji, jj) * pssq(ji, jj) 192 zztmp1 = pevp(ji, jj) 193 zztmp2 = pevp(ji, jj) * pssq(ji, jj) 194 #if defined key_si3 195 zztmp1 = zztmp1 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pevp_ice(ji,jj) 196 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pevp_ice(ji, jj) * pssq_ice(ji, jj) 197 #endif 198 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * zztmp1 199 tq_abl ( ji, jj, 2, nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) + rdt * zztmp2 166 200 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) 167 201 END DO … … 251 285 ENDDO 252 286 END IF 253 287 254 288 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 255 289 ! ! 4 *** Advance u,v to time n+1 … … 273 307 z_elem_a( ji, 2 ) = 0._wp 274 308 z_elem_c( ji, 2 ) = - rdt * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) 275 zcff = rdt * pcd_du(ji, jj) 276 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + zcff 277 u_abl( ji, jj, 2, nt_a ) = u_abl( ji, jj, 2, nt_a ) + 0.5_wp * zcff * ( pssu(ji-1, jj) + pssu(ji,jj) ) 309 ! 310 zztmp1 = pcd_du(ji, jj) 311 zztmp2 = 0.5_wp * pcd_du(ji, jj) * ( pssu(ji-1, jj) + pssu(ji,jj) ) 312 #if defined key_si3 313 zztmp1 = zztmp1 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) 314 zzice = 0.5_wp * ( pssu_ice(ji-1, jj) + pssu_ice(ji,jj) ) 315 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) * zzice 316 #endif 317 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * zztmp1 318 u_abl( ji, jj, 2, nt_a ) = u_abl( ji, jj, 2, nt_a ) + rdt * zztmp2 319 278 320 !++ Top Neumann B.C. 279 321 !z_elem_a( ji, jpka ) = - 0.5_wp * rdt * ( Avm_abl( ji, jj, jpka )+ Avm_abl( ji+1, jj, jpka ) ) / e3w_abl( jpka ) … … 331 373 !++ Surface boundary condition 332 374 z_elem_a( ji, 2 ) = 0._wp 333 z_elem_c( ji, 2 ) = - rdt * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) 334 zcff = rdt * pcd_du(ji, jj) 335 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + zcff 336 v_abl( ji, jj, 2, nt_a ) = v_abl( ji, jj, 2, nt_a ) + 0.5_wp * zcff * ( pssv(ji, jj) + pssv(ji, jj-1) ) 375 z_elem_c( ji, 2 ) = - rdt * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) 376 ! 377 zztmp1 = pcd_du(ji, jj) 378 zztmp2 = 0.5_wp * pcd_du(ji, jj) * ( pssv(ji, jj) + pssv(ji, jj-1) ) 379 #if defined key_si3 380 zztmp1 = zztmp1 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) 381 zzice = 0.5_wp * ( pssv_ice(ji, jj) + pssv_ice(ji,jj-1) ) 382 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) * zzice 383 #endif 384 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rdt * zztmp1 385 v_abl( ji, jj, 2, nt_a ) = v_abl( ji, jj, 2, nt_a ) + rdt * zztmp2 337 386 !++ Top Neumann B.C. 338 387 !z_elem_a( ji, jpka ) = -rdt * Avm_abl( ji, jj, jpka ) / e3w_abl( jpka ) … … 372 421 END DO ! end outer loop 373 422 !------------- 374 423 375 424 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 376 425 ! ! 5 *** Apply nudging on the dynamics and the tracers … … 428 477 !------------- 429 478 END DO ! end outer loop 430 !------------- 479 !------------- 431 480 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 432 481 ! ! 6 *** MPI exchanges 433 482 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 434 435 CALL lbc_lnk_multi( 'ablmod', u_abl(:,:,:,nt_a), 'T', -1., v_abl(:,:,:,nt_a), kfillmode = jpfillnothing )483 ! 484 CALL lbc_lnk_multi( 'ablmod', u_abl(:,:,:,nt_a), 'T', -1., v_abl(:,:,:,nt_a), 'T', -1., kfillmode = jpfillnothing ) 436 485 CALL lbc_lnk_multi( 'ablmod', tq_abl(:,:,:,nt_a,jp_ta), 'T', 1., tq_abl(:,:,:,nt_a,jp_qa), 'T', 1., kfillmode = jpfillnothing ) ! ++++ this should not be needed... 437 486 ! 438 487 CALL iom_put ( 'u_abl', u_abl(:,:,2:jpka,nt_a )) 439 488 CALL iom_put ( 'v_abl', v_abl(:,:,2:jpka,nt_a )) … … 504 553 CALL prt_ctl( tab2d_2=ptauj , clinfo2= 'vtau : ' ) 505 554 ENDIF 506 555 556 #if defined key_si3 557 ! ------------------------------------------------------------ ! 558 ! Wind stress relative to the moving ice ( U10m - U_ice ) ! 559 ! ------------------------------------------------------------ ! 560 DO jj = 2, jpjm1 561 DO ji = 2, jpim1 562 563 zztmp1 = 0.5_wp * ( u_abl(ji+1,jj,2,nt_a) + u_abl(ji,jj,2,nt_a) ) 564 zztmp2 = 0.5_wp * ( v_abl(ji,jj+1,2,nt_a) + v_abl(ji,jj,2,nt_a) ) 565 566 ptaui_ice(ji,jj) = 0.5_wp * ( zrhoa(ji+1,jj) * pCd_du_ice(ji+1,jj) & 567 & + zrhoa(ji ,jj) * pCd_du_ice(ji ,jj) ) & 568 & * ( zztmp1 - rn_vfac * pssu_ice(ji,jj) ) 569 ptauj_ice(ji,jj) = 0.5_wp * ( zrhoa(ji,jj+1) * pCd_du_ice(ji,jj+1) & 570 & + zrhoa(ji,jj ) * pCd_du_ice(ji,jj ) ) & 571 & * ( zztmp2 - rn_vfac * pssv_ice(ji,jj) ) 572 END DO 573 END DO 574 CALL lbc_lnk_multi( 'ablmod', ptaui_ice, 'U', -1., ptauj_ice, 'V', -1. ) 575 ! 576 IF(ln_ctl) CALL prt_ctl( tab2d_1=ptaui_ice , clinfo1=' abl_stp: putaui : ' & 577 & , tab2d_2=ptauj_ice , clinfo2=' pvtaui : ' ) 578 #endif 507 579 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 508 580 ! ! 8 *** Swap time indices for the next timestep -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/ABL/sbcabl.F90
r11322 r11334 30 30 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 31 31 USE prtctl ! Print control 32 32 #if defined key_si3 33 USE ice , ONLY : u_ice, v_ice, tm_su, ato_i ! ato_i = total open water fractional area 34 USE sbc_ice, ONLY : wndm_ice, utau_ice, vtau_ice 35 #endif 36 #if ! defined key_iomput 37 USE diawri , ONLY : dia_wri_alloc_abl 38 #endif 33 39 IMPLICIT NONE 34 40 PRIVATE … … 141 147 CALL iom_get( inum, jpdom_unknown, 'ghw_abl', ghw_abl(:) ) 142 148 CALL iom_close( inum ) 149 150 #if ! defined key_iomput 151 IF( dia_wri_alloc_abl() /= 0 ) CALL ctl_stop( 'STOP', 'abl_init : unable to allocate arrays' ) 152 #endif 143 153 144 154 IF(lwp) THEN … … 303 313 !! 304 314 REAL(wp), DIMENSION(jpi,jpj) :: zssq, zcd_du, zsen, zevp 315 #if defined key_si3 316 REAL(wp), DIMENSION(jpi,jpj) :: zssqi, zcd_dui, zseni, zevpi 317 #endif 305 318 INTEGER :: jbak, jbak_dta, ji, jj 306 319 !!--------------------------------------------------------------------- … … 319 332 & sf(jp_slp )%fnow(:,:,1), sst_m, ssu_m, ssv_m, & 320 333 & zssq, zcd_du, zsen, zevp ) 321 334 335 #if defined key_si3 336 CALL blk_ice_1( u_abl(:,:,2,nt_n), v_abl(:,:,2,nt_n), tq_abl(:,:,2,nt_n,jp_ta), tq_abl(:,:,2,nt_n,jp_qa), & 337 & sf(jp_slp)%fnow(:,:,1), u_ice, v_ice, & 338 & pseni=zseni, pevpi=zevpi, ptsui=tm_su, pssqi=zssqi, pcd_dui=zcd_dui ) ! outputs 339 #endif 340 322 341 !!------------------------------------------------------------------------------------------- 323 342 !! 3 - Advance ABL variables from now (n) to after (n+1) … … 329 348 & sf(jp_hpgi)%fnow(:,:,:), sf(jp_hpgj)%fnow(:,:,:), & 330 349 & zcd_du, zsen, zevp, & ! in/out 331 & wndm, utau, vtau, taum ) ! out 332 350 & wndm, utau, vtau, taum & ! out 351 #if defined key_si3 352 & , tm_su, u_ice, v_ice, zssqi, zcd_dui & 353 & , zseni, zevpi, wndm_ice, ato_i, & 354 & utau_ice, vtau_ice & 355 #endif 356 & ) 333 357 !!------------------------------------------------------------------------------------------- 334 358 !! 4 - Finalize flux computation using ABL variables at (n+1), nt_n corresponds to (n+1) since … … 340 364 & sst_m , zsen, zevp ) 341 365 366 #if defined key_si3 367 ! Avoid a USE abl in icesbc module 368 sf(jp_tair)%fnow(:,:,1) = tq_abl(:,:,2,nt_n,jp_ta); sf(jp_humi)%fnow(:,:,1) = tq_abl(:,:,2,nt_n,jp_qa) 369 #endif 370 342 371 END SUBROUTINE sbc_abl 343 372 -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/ICE/icesbc.F90
r11322 r11334 75 75 & u_ice, v_ice, & ! inputs 76 76 & putaui = utau_ice, pvtaui = vtau_ice ) ! outputs 77 CASE( jp_abl ) ; CALL blk_ice_1( sf(jp_wndi)%fnow(:,:,1), sf(jp_wndj)%fnow(:,:,1), & 78 & sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1), & 79 & u_ice, v_ice, & ! inputs 80 & putaui = utau_ice, pvtaui = vtau_ice ) ! outputs 77 ! CASE( jp_abl ) utau_ice & vtau_ice are computed in ablmod 81 78 CASE( jp_purecpl ) ; CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) ! Coupled formulation 82 79 END SELECT … … 146 143 CASE( jp_usr ) !--- user defined formulation 147 144 CALL usrdef_sbc_ice_flx( kt, h_s, h_i ) 148 CASE( jp_blk, jp_abl ) !--- bulk formulation (temporary version for ABL option)145 CASE( jp_blk, jp_abl ) !--- bulk formulation & ABL formulation 149 146 CALL blk_ice_2 ( t_su, h_s, h_i, alb_ice, sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), & 150 147 & sf(jp_slp)%fnow(:,:,1), sf(jp_qlw)%fnow(:,:,1), sf(jp_prec)%fnow(:,:,1), sf(jp_snow)%fnow(:,:,1) ) ! -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/ICE/icevar.F90
r11223 r11334 114 114 ! 115 115 ato_i(:,:) = 1._wp - at_i(:,:) ! open water fraction 116 116 ! 117 ALLOCATE( z1_at_i(jpi,jpj) ) 118 WHERE( at_i(:,:) > epsi20 ) ; z1_at_i(:,:) = 1._wp / at_i(:,:) 119 ELSEWHERE ; z1_at_i(:,:) = 0._wp 120 END WHERE 121 ! 122 tm_su(:,:) = SUM( t_su(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:) 123 WHERE( at_i(:,:)<=epsi20 ); tm_su(:,:) = rt0; END WHERE 124 ! 117 125 ! The following fields are calculated for diagnostics and outputs only 118 126 ! ==> Do not use them for other purposes 119 127 IF( kn > 1 ) THEN 120 128 ! 121 ALLOCATE( z1_at_i(jpi,jpj) , z1_vt_i(jpi,jpj) , z1_vt_s(jpi,jpj) ) 122 WHERE( at_i(:,:) > epsi20 ) ; z1_at_i(:,:) = 1._wp / at_i(:,:) 123 ELSEWHERE ; z1_at_i(:,:) = 0._wp 124 END WHERE 129 ALLOCATE( z1_vt_i(jpi,jpj) , z1_vt_s(jpi,jpj) ) 125 130 WHERE( vt_i(:,:) > epsi20 ) ; z1_vt_i(:,:) = 1._wp / vt_i(:,:) 126 131 ELSEWHERE ; z1_vt_i(:,:) = 0._wp … … 135 140 ! 136 141 ! ! mean temperature (K), salinity and age 137 tm_su(:,:) = SUM( t_su(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:)138 142 tm_si(:,:) = SUM( t_si(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:) 139 143 om_i (:,:) = SUM( oa_i(:,:,:) , dim=3 ) * z1_at_i(:,:) … … 153 157 ! ! put rt0 where there is no ice 154 158 WHERE( at_i(:,:)<=epsi20 ) 155 tm_su(:,:) = rt0156 159 tm_si(:,:) = rt0 157 160 tm_i (:,:) = rt0 … … 159 162 END WHERE 160 163 161 DEALLOCATE( z1_ at_i , z1_vt_i , z1_vt_s )164 DEALLOCATE( z1_vt_i , z1_vt_s ) 162 165 ENDIF 166 ! 167 DEALLOCATE( z1_at_i ) 163 168 ! 164 169 END SUBROUTINE ice_var_agg -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/OCE/DIA/diawri.F90
r10425 r11334 26 26 !!---------------------------------------------------------------------- 27 27 USE oce ! ocean dynamics and tracers 28 USE abl ! abl variables in case ln_abl = .true. 28 29 USE dom_oce ! ocean space and time domain 29 30 USE phycst ! physical constants … … 66 67 PUBLIC dia_wri_state 67 68 PUBLIC dia_wri_alloc ! Called by nemogcm module 68 69 #if ! defined key_iomput 70 PUBLIC dia_wri_alloc_abl ! Called by sbcabl module (if ln_abl = .true.) 71 #endif 69 72 INTEGER :: nid_T, nz_T, nh_T, ndim_T, ndim_hT ! grid_T file 70 INTEGER :: nb_T , ndim_bT! grid_T file73 INTEGER :: nb_T, ndim_bT ! grid_T file 71 74 INTEGER :: nid_U, nz_U, nh_U, ndim_U, ndim_hU ! grid_U file 72 75 INTEGER :: nid_V, nz_V, nh_V, ndim_V, ndim_hV ! grid_V file 73 76 INTEGER :: nid_W, nz_W, nh_W ! grid_W file 77 INTEGER :: ndim_A, ndim_hA ! grid_T file 78 INTEGER :: nid_A, nz_A, nh_A ! grid_ABL file 74 79 INTEGER :: ndex(1) ! ??? 75 INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV 76 INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V 80 INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_hT, ndex_hU, ndex_hV, ndex_hA 81 INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_T, ndex_U, ndex_V, ndex_A 77 82 INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) :: ndex_bT 78 83 … … 409 414 & ndex_hU(jpi*jpj) , ndex_U(jpi*jpj*jpk) , & 410 415 & ndex_hV(jpi*jpj) , ndex_V(jpi*jpj*jpk) , STAT=ierr(1) ) 411 ! 412 dia_wri_alloc = MAXVAL(ierr) 416 dia_wri_alloc = MAXVAL(ierr) 413 417 CALL mpp_sum( 'diawri', dia_wri_alloc ) 414 418 ! 415 419 END FUNCTION dia_wri_alloc 416 417 420 421 INTEGER FUNCTION dia_wri_alloc_abl() 422 !!---------------------------------------------------------------------- 423 ALLOCATE( ndex_hA(jpi*jpj), ndex_A (jpi*jpj*jpkam1), STAT=dia_wri_alloc_abl) 424 CALL mpp_sum( 'diawri', dia_wri_alloc_abl ) 425 ! 426 END FUNCTION dia_wri_alloc_abl 427 418 428 SUBROUTINE dia_wri( kt ) 419 429 !!--------------------------------------------------------------------- … … 435 445 INTEGER :: ji, jj, jk ! dummy loop indices 436 446 INTEGER :: ierr ! error code return from allocation 437 INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma ! local integers447 INTEGER :: iimi, iima, ipk, it, itmod, ijmi, ijma, ipka ! local integers 438 448 INTEGER :: jn, ierror ! local integers 439 449 REAL(wp) :: zsto, zout, zmax, zjulian ! local scalars … … 441 451 REAL(wp), DIMENSION(jpi,jpj) :: zw2d ! 2D workspace 442 452 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d ! 3D workspace 453 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d_abl ! 3D workspace 443 454 !!---------------------------------------------------------------------- 444 455 ! … … 449 460 ninist = 0 450 461 ENDIF 462 451 463 ! 452 464 ! 0. Initialisation … … 472 484 ijmi = 1 ; ijma = jpj 473 485 ipk = jpk 486 IF(ln_abl) ipka = jpkam1 474 487 475 488 ! define time axis … … 574 587 & "m", ipk, gdepw_1d, nz_W, "down" ) 575 588 589 IF( ln_abl ) THEN 590 ! Define the ABL grid FILE ( nid_A ) 591 CALL dia_nam( clhstnam, nwrite, 'grid_ABL' ) 592 IF(lwp) WRITE(numout,*) " Name of NETCDF file ", clhstnam ! filename 593 CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit 594 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 595 & nit000-1, zjulian, rdt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set ) 596 CALL histvert( nid_A, "ght_abl", "Vertical T levels", & ! Vertical grid: gdept 597 & "m", ipka, ght_abl(2:jpka), nz_A, "up" ) 598 ! ! Index of ocean points 599 ALLOCATE( zw3d_abl(jpi,jpj,ipka) ) 600 zw3d_abl(:,:,:) = 1._wp 601 CALL wheneq( jpi*jpj*ipka, zw3d_abl, 1, 1., ndex_A , ndim_A ) ! volume 602 CALL wheneq( jpi*jpj , zw3d_abl, 1, 1., ndex_hA, ndim_hA ) ! surface 603 DEALLOCATE(zw3d_abl) 604 ENDIF 576 605 577 606 ! Declare all the output fields as NETCDF variables … … 623 652 CALL histdef( nid_T, "sowindsp", "wind speed at 10m" , "m/s" , & ! wndm 624 653 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 654 ! 655 IF( ln_abl ) THEN 656 ! !!! nid_A : 3D 657 CALL histdef( nid_A, "t_abl", "Potential Temperature" , "K" , & ! t_abl 658 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 659 CALL histdef( nid_A, "q_abl", "Humidity" , "kg/kg" , & ! q_abl 660 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 661 CALL histdef( nid_A, "u_abl", "Atmospheric U-wind " , "m/s" , & ! u_abl 662 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 663 CALL histdef( nid_A, "v_abl", "Atmospheric V-wind " , "m/s" , & ! v_abl 664 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 665 CALL histdef( nid_A, "tke_abl", "Atmospheric TKE " , "m2/s2" , & ! tke_abl 666 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 667 CALL histdef( nid_A, "avm_abl", "Atmospheric turbulent viscosity", "m2/s" , & ! avm_abl 668 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 669 CALL histdef( nid_A, "avt_abl", "Atmospheric turbulent diffusivity", "m2/s2", & ! avt_abl 670 & jpi, jpj, nh_A, ipka, 1, ipka, nz_A, 32, clop, zsto, zout ) 671 CALL histdef( nid_A, "pblh", "Atmospheric boundary layer height " , "m", & ! pblh 672 & jpi, jpj, nh_A, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 673 #if defined key_si3 674 CALL histdef( nid_A, "oce_frac", "Fraction of open ocean" , " ", & ! ato_i 675 & jpi, jpj, nh_A, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 676 #endif 677 CALL histend( nid_A, snc4chunks=snc4set ) 678 ! 679 ENDIF 625 680 ! 626 681 IF( ln_icebergs ) THEN … … 676 731 677 732 clmx ="l_max(only(x))" ! max index on a period 678 !CALL histdef( nid_T, "sobowlin", "Bowl Index" , "W-point", & ! bowl INDEX679 !& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clmx, zsto, zout )733 CALL histdef( nid_T, "sobowlin", "Bowl Index" , "W-point", & ! bowl INDEX 734 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clmx, zsto, zout ) 680 735 #if defined key_diahth 681 736 CALL histdef( nid_T, "sothedep", "Thermocline Depth" , "m" , & ! hth … … 790 845 CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction 791 846 CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed 847 ! 848 IF( ln_abl ) THEN 849 ALLOCATE( zw3d_abl(jpi,jpj,jpka) ) 850 IF( ln_mskland ) THEN 851 DO jk=1,jpka 852 zw3d_abl(:,:,jk) = tmask(:,:,1) 853 END DO 854 ELSE 855 zw3d_abl(:,:,:) = 1._wp 856 ENDIF 857 CALL histwrite( nid_A, "pblh" , it, pblh(:,:) *zw3d_abl(:,:,1 ), ndim_hA, ndex_hA ) ! pblh 858 CALL histwrite( nid_A, "u_abl" , it, u_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! u_abl 859 CALL histwrite( nid_A, "v_abl" , it, v_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! v_abl 860 CALL histwrite( nid_A, "t_abl" , it, tq_abl (:,:,2:jpka,nt_n,1)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! t_abl 861 CALL histwrite( nid_A, "q_abl" , it, tq_abl (:,:,2:jpka,nt_n,2)*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! q_abl 862 CALL histwrite( nid_A, "tke_abl", it, tke_abl (:,:,2:jpka,nt_n )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! tke_abl 863 CALL histwrite( nid_A, "avm_abl", it, avm_abl (:,:,2:jpka )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! avm_abl 864 CALL histwrite( nid_A, "avt_abl", it, avt_abl (:,:,2:jpka )*zw3d_abl(:,:,2:jpka), ndim_A , ndex_A ) ! avt_abl 865 #if defined key_si3 866 CALL histwrite( nid_A, "oce_frac" , it, ato_i(:,:) , ndim_hA, ndex_hA ) ! ato_i 867 #endif 868 DEALLOCATE(zw3d_abl) 869 ENDIF 792 870 ! 793 871 IF( ln_icebergs ) THEN … … 826 904 CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping 827 905 ENDIF 828 !zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1)829 !CALL histwrite( nid_T, "sobowlin", it, zw2d , ndim_hT, ndex_hT ) ! ???906 zw2d(:,:) = FLOAT( nmln(:,:) ) * tmask(:,:,1) 907 CALL histwrite( nid_T, "sobowlin", it, zw2d , ndim_hT, ndex_hT ) ! ??? 830 908 831 909 #if defined key_diahth … … 862 940 CALL histclo( nid_V ) 863 941 CALL histclo( nid_W ) 942 IF(ln_abl) CALL histclo( nid_A ) 864 943 ENDIF 865 944 ! -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/OCE/SBC/sbcabl.F90
r11275 r11334 8 8 !!---------------------------------------------------------------------- 9 9 USE sbc_oce, ONLY : ght_abl, ghw_abl, e3t_abl, e3w_abl 10 USE lib_mpp, ONLY : ctl_stop 10 11 11 12 IMPLICIT NONE … … 31 32 !! 32 33 !!---------------------------------------------------------------------- 33 34 CALL ctl_stop( 'STOP', 'ln_abl = .true. but ABL source directory was not included', & 35 & '(Either switch to ln_abl = .false. or modify your cfg.txt file and recompile)' ) 36 !! 34 37 END SUBROUTINE sbc_abl_init 35 38 -
NEMO/branches/2019/dev_r11265_ASINTER-01_Guillaume_ABL1D/src/OCE/SBC/sbcblk.F90
r11322 r11334 816 816 ! 817 817 INTEGER :: ji, jj ! dummy loop indices 818 REAL(wp) :: zwndi_t , zwndj_t 818 REAL(wp) :: zwndi_t , zwndj_t, zootm_su, zztmp1, zztmp2 ! relative wind components at T-point 819 819 REAL(wp), DIMENSION(jpi,jpj) :: zrhoa ! transfer coefficient for momentum (tau) 820 820 REAL(wp), DIMENSION(jpi,jpj) :: zcd_dui ! transfer coefficient for momentum (tau) … … 855 855 zcd_dui(:,:) = wndm_ice(:,:) * Cd_atm(:,:) 856 856 857 !IF( ln_blk ) THEN857 IF( ln_blk ) THEN 858 858 ! ------------------------------------------------------------ ! 859 859 ! Wind stress relative to the moving ice ( U10m - U_ice ) ! … … 874 874 IF(ln_ctl) CALL prt_ctl( tab2d_1=putaui , clinfo1=' blk_ice: putaui : ' & 875 875 & , tab2d_2=pvtaui , clinfo2=' pvtaui : ' ) 876 ! ELSE 877 ! DO jj = 1, jpj 878 ! DO ji = 1, jpi 879 ! pcd_dui(ji,jj) = zcd_dui(ji,jj) 880 ! pseni (ji,jj) = wndm_ice(ji,jj) * Ch_atm(ji,jj) 881 ! pevpi (ji,jj) = wndm_ice(ji,jj) * Ce_atm(ji,jj) 882 ! pssqi (ji,jj) = 11637800.0_wp * EXP( -5897.8_wp / (ptsui(ji,jj)+rt0) ) / zrhoa(ji,jj) 883 ! END DO 884 ! END DO 885 ! ENDIF 876 ELSE 877 zztmp1 = 11637800.0_wp 878 zztmp2 = -5897.8_wp 879 DO jj = 1, jpj 880 DO ji = 1, jpi 881 pcd_dui(ji,jj) = zcd_dui (ji,jj) 882 pseni (ji,jj) = wndm_ice(ji,jj) * Ch_atm(ji,jj) 883 pevpi (ji,jj) = wndm_ice(ji,jj) * Ce_atm(ji,jj) 884 zootm_su = zztmp2 / ptsui(ji,jj) ! ptsui is in K (it can't be zero ??) 885 pssqi (ji,jj) = zztmp1 * EXP( zootm_su ) / zrhoa(ji,jj) 886 END DO 887 END DO 888 ENDIF 886 889 ! 887 890 IF(ln_ctl) CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice: wndm_ice : ') … … 1248 1251 REAL(wp) :: z0w, z0i, zfmi, zfmw, zfhi, zfhw 1249 1252 REAL(wp) :: zCdn_form_tmp 1250 !!---------------------------------------------------------------------- 1251 1253 !!--------------------------------------------------------------------------- 1254 ! fl: global variable tm_su could be used directly instead of recomputing it 1255 ! 1252 1256 ! mean temperature 1253 1257 WHERE( at_i_b(:,:) > 1.e-20 ) ; ztm_su(:,:) = SUM( t_su(:,:,:) * a_i_b(:,:,:) , dim=3 ) / at_i_b(:,:)
Note: See TracChangeset
for help on using the changeset viewer.
