Changeset 14007 for NEMO/trunk/src/OCE/SBC
- Timestamp:
- 2020-12-02T15:30:51+01:00 (4 years ago)
- Location:
- NEMO/trunk/src/OCE/SBC
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/SBC/cpl_oasis3.F90
r13415 r14007 66 66 INTEGER :: nsnd ! total number of fields sent 67 67 INTEGER :: ncplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 68 INTEGER, PUBLIC, PARAMETER :: nmaxfld=6 0! Maximum number of coupling fields68 INTEGER, PUBLIC, PARAMETER :: nmaxfld=62 ! Maximum number of coupling fields 69 69 INTEGER, PUBLIC, PARAMETER :: nmaxcat=5 ! Maximum number of coupling fields 70 70 INTEGER, PUBLIC, PARAMETER :: nmaxcpl=5 ! Maximum number of coupling fields -
NEMO/trunk/src/OCE/SBC/sbc_oce.F90
r13472 r14007 12 12 !! 4.0 ! 2016-06 (L. Brodeau) new unified bulk routine (based on AeroBulk) 13 13 !! 4.0 ! 2019-03 (F. Lemarié, G. Samson) add compatibility with ABL mode 14 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) modified wave parameters in namelist 14 15 !!---------------------------------------------------------------------- 15 16 … … 36 37 LOGICAL , PUBLIC :: ln_blk !: bulk formulation 37 38 LOGICAL , PUBLIC :: ln_abl !: Atmospheric boundary layer model 39 LOGICAL , PUBLIC :: ln_wave !: wave in the system (forced or coupled) 38 40 #if defined key_oasis3 39 41 LOGICAL , PUBLIC :: lk_oasis = .TRUE. !: OASIS used … … 56 58 ! !: = 1 global mean of e-p-r set to zero at each nn_fsbc time step 57 59 ! !: = 2 annual global mean of e-p-r set to zero 58 LOGICAL , PUBLIC :: ln_wave !: true if some coupling with wave model59 LOGICAL , PUBLIC :: ln_cdgw !: true if neutral drag coefficient from wave model60 LOGICAL , PUBLIC :: ln_sdw !: true if 3d stokes drift from wave model61 LOGICAL , PUBLIC :: ln_tauwoc !: true if normalized stress from wave is used62 LOGICAL , PUBLIC :: ln_tauw !: true if ocean stress components from wave is used63 LOGICAL , PUBLIC :: ln_stcor !: true if Stokes-Coriolis term is used64 !65 INTEGER , PUBLIC :: nn_sdrift ! type of parameterization to calculate vertical Stokes drift66 !67 60 LOGICAL , PUBLIC :: ln_icebergs !: Icebergs 68 61 ! … … 71 64 ! !!* namsbc_cpl namelist * 72 65 INTEGER , PUBLIC :: nn_cats_cpl !: Number of sea ice categories over which the coupling is carried out 73 66 ! 67 ! !!* namsbc_wave namelist * 68 LOGICAL , PUBLIC :: ln_sdw !: =T 3d stokes drift from wave model 69 LOGICAL , PUBLIC :: ln_stcor !: =T if Stokes-Coriolis and tracer advection terms are used 70 LOGICAL , PUBLIC :: ln_cdgw !: =T neutral drag coefficient from wave model 71 LOGICAL , PUBLIC :: ln_tauoc !: =T if normalized stress from wave is used 72 LOGICAL , PUBLIC :: ln_wave_test !: =T wave test case (constant Stokes drift) 73 LOGICAL , PUBLIC :: ln_charn !: =T Chranock coefficient from wave model 74 LOGICAL , PUBLIC :: ln_taw !: =T wind stress corrected by wave intake 75 LOGICAL , PUBLIC :: ln_phioc !: =T TKE surface BC from wave model 76 LOGICAL , PUBLIC :: ln_bern_srfc !: Bernoulli head, waves' inuced pressure 77 LOGICAL , PUBLIC :: ln_breivikFV_2016 !: Breivik 2016 profile 78 LOGICAL , PUBLIC :: ln_vortex_force !: vortex force activation 79 LOGICAL , PUBLIC :: ln_stshear !: Stoked Drift shear contribution in zdftke 80 ! 74 81 !!---------------------------------------------------------------------- 75 82 !! switch definition (improve readability) … … 81 88 INTEGER , PUBLIC, PARAMETER :: jp_purecpl = 5 !: Pure ocean-atmosphere Coupled formulation 82 89 INTEGER , PUBLIC, PARAMETER :: jp_none = 6 !: for OPA when doing coupling via SAS module 83 84 !!---------------------------------------------------------------------- 85 !! Stokes drift parametrization definition 86 !!---------------------------------------------------------------------- 87 INTEGER , PUBLIC, PARAMETER :: jp_breivik_2014 = 0 !: Breivik 2014: v_z=v_0*[exp(2*k*z)/(1-8*k*z)] 88 INTEGER , PUBLIC, PARAMETER :: jp_li_2017 = 1 !: Li et al 2017: Stokes drift based on Phillips spectrum (Breivik 2016) 89 ! with depth averaged profile 90 INTEGER , PUBLIC, PARAMETER :: jp_peakfr = 2 !: Li et al 2017: using the peak wave number read from wave model instead 91 ! of the inverse depth scale 92 LOGICAL , PUBLIC :: ll_st_bv2014 = .FALSE. ! logical indicator, .true. if Breivik 2014 parameterisation is active. 93 LOGICAL , PUBLIC :: ll_st_li2017 = .FALSE. ! logical indicator, .true. if Li 2017 parameterisation is active. 94 LOGICAL , PUBLIC :: ll_st_bv_li = .FALSE. ! logical indicator, .true. if either Breivik or Li parameterisation is active. 95 LOGICAL , PUBLIC :: ll_st_peakfr = .FALSE. ! logical indicator, .true. if using Li 2017 with peak wave number 96 90 ! 97 91 !!---------------------------------------------------------------------- 98 92 !! component definition -
NEMO/trunk/src/OCE/SBC/sbcblk.F90
r13501 r14007 314 314 ENDIF 315 315 END DO 316 !317 IF( ln_wave ) THEN318 !Activated wave module but neither drag nor stokes drift activated319 IF( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauwoc .OR. ln_stcor ) ) THEN320 CALL ctl_stop( 'STOP', 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauwoc=F, ln_stcor=F' )321 !drag coefficient read from wave model definable only with mfs bulk formulae and core322 ELSEIF(ln_cdgw .AND. .NOT. ln_NCAR ) THEN323 CALL ctl_stop( 'drag coefficient read from wave model definable only with NCAR and CORE bulk formulae')324 ELSEIF(ln_stcor .AND. .NOT. ln_sdw) THEN325 CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T')326 ENDIF327 ELSE328 IF( ln_cdgw .OR. ln_sdw .OR. ln_tauwoc .OR. ln_stcor ) &329 & CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ', &330 & 'with drag coefficient (ln_cdgw =T) ' , &331 & 'or Stokes Drift (ln_sdw=T) ' , &332 & 'or ocean stress modification due to waves (ln_tauwoc=T) ', &333 & 'or Stokes-Coriolis term (ln_stcori=T)' )334 ENDIF335 316 ! 336 317 IF( ln_abl ) THEN ! ABL: read 3D fields for wind, temperature, humidity and pressure gradient -
NEMO/trunk/src/OCE/SBC/sbcblk_algo_ecmwf.F90
r13460 r14007 17 17 !!---------------------------------------------------------------------- 18 18 !! History : 4.0 ! 2016-02 (L.Brodeau) Original code 19 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) Charnock coeff from wave model 19 20 !!---------------------------------------------------------------------- 20 21 … … 31 32 USE in_out_manager ! I/O manager 32 33 USE prtctl ! Print control 33 USE sbcwave, ONLY : cdn_wave! wave module34 USE sbcwave, ONLY : charn ! wave module 34 35 #if defined key_si3 || defined key_cice 35 36 USE sbc_ice ! Surface boundary condition: ice fields … … 233 234 u_star = 0.035_wp*U_blk*ztmp1/ztmp0 ! (u* = 0.035*Un10) 234 235 235 z0 = charn0*u_star*u_star/grav + 0.11_wp*znu_a/u_star 236 IF (ln_charn) THEN ! Charnock value if wave coupling 237 z0 = charn*u_star*u_star/grav + 0.11_wp*znu_a/u_star 238 ELSE 239 z0 = charn0*u_star*u_star/grav + 0.11_wp*znu_a/u_star 240 ENDIF 241 236 242 z0 = MIN( MAX(ABS(z0), 1.E-9) , 1._wp ) ! (prevents FPE from stupid values from masked region later on) 237 243 … … 302 308 ztmp2 = u_star*u_star 303 309 ztmp1 = znu_a/u_star 304 z0 = MIN( ABS( alpha_M*ztmp1 + charn0*ztmp2/grav ) , 0.001_wp) 310 IF (ln_charn) THEN ! Charnock value if wave coupling 311 z0 = MIN( ABS( alpha_M*ztmp1 + charn*ztmp2/grav ) , 0.001_wp) 312 ELSE 313 z0 = MIN( ABS( alpha_M*ztmp1 + charn0*ztmp2/grav ) , 0.001_wp) 314 ENDIF 305 315 z0t = MIN( ABS( alpha_H*ztmp1 ) , 0.001_wp) ! eq.3.26, Chap.3, p.34, IFS doc - Cy31r1 306 316 z0q = MIN( ABS( alpha_Q*ztmp1 ) , 0.001_wp) -
NEMO/trunk/src/OCE/SBC/sbccpl.F90
r13497 r14007 8 8 !! 3.1 ! 2009_02 (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface 9 9 !! 3.4 ! 2011_11 (C. Harris) more flexibility + multi-category fields 10 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) wave coupling updates 10 11 !!---------------------------------------------------------------------- 11 12 … … 106 107 INTEGER, PARAMETER :: jpr_fraqsr = 42 ! fraction of solar net radiation absorbed in the first ocean level 107 108 INTEGER, PARAMETER :: jpr_mslp = 43 ! mean sea level pressure 108 INTEGER, PARAMETER :: jpr_hsig = 44 ! Hsig 109 INTEGER, PARAMETER :: jpr_phioc = 45 ! Wave=>ocean energy flux 110 INTEGER, PARAMETER :: jpr_sdrftx = 46 ! Stokes drift on grid 1 111 INTEGER, PARAMETER :: jpr_sdrfty = 47 ! Stokes drift on grid 2 109 !** surface wave coupling ** 110 INTEGER, PARAMETER :: jpr_hsig = 44 ! Hsig 111 INTEGER, PARAMETER :: jpr_phioc = 45 ! Wave=>ocean energy flux 112 INTEGER, PARAMETER :: jpr_sdrftx = 46 ! Stokes drift on grid 1 113 INTEGER, PARAMETER :: jpr_sdrfty = 47 ! Stokes drift on grid 2 112 114 INTEGER, PARAMETER :: jpr_wper = 48 ! Mean wave period 113 115 INTEGER, PARAMETER :: jpr_wnum = 49 ! Mean wavenumber 114 INTEGER, PARAMETER :: jpr_ tauwoc= 50 ! Stress fraction adsorbed by waves116 INTEGER, PARAMETER :: jpr_wstrf = 50 ! Stress fraction adsorbed by waves 115 117 INTEGER, PARAMETER :: jpr_wdrag = 51 ! Neutral surface drag coefficient 116 INTEGER, PARAMETER :: jpr_isf = 52 117 INTEGER, PARAMETER :: jpr_icb = 53 118 INTEGER, PARAMETER :: jpr_wfreq = 54 ! Wave peak frequency 119 INTEGER, PARAMETER :: jpr_tauwx = 55 ! x component of the ocean stress from waves 120 INTEGER, PARAMETER :: jpr_tauwy = 56 ! y component of the ocean stress from waves 121 INTEGER, PARAMETER :: jpr_ts_ice = 57 ! Sea ice surface temp 122 123 INTEGER, PARAMETER :: jprcv = 57 ! total number of fields received 118 INTEGER, PARAMETER :: jpr_charn = 52 ! Chranock coefficient 119 INTEGER, PARAMETER :: jpr_twox = 53 ! wave to ocean momentum flux 120 INTEGER, PARAMETER :: jpr_twoy = 54 ! wave to ocean momentum flux 121 INTEGER, PARAMETER :: jpr_tawx = 55 ! net wave-supported stress 122 INTEGER, PARAMETER :: jpr_tawy = 56 ! net wave-supported stress 123 INTEGER, PARAMETER :: jpr_bhd = 57 ! Bernoulli head. waves' induced surface pressure 124 INTEGER, PARAMETER :: jpr_tusd = 58 ! zonal stokes transport 125 INTEGER, PARAMETER :: jpr_tvsd = 59 ! meridional stokes tranmport 126 INTEGER, PARAMETER :: jpr_isf = 60 127 INTEGER, PARAMETER :: jpr_icb = 61 128 INTEGER, PARAMETER :: jpr_ts_ice = 62 ! Sea ice surface temp 129 130 INTEGER, PARAMETER :: jprcv = 62 ! total number of fields received 124 131 125 132 INTEGER, PARAMETER :: jps_fice = 1 ! ice fraction sent to the atmosphere … … 184 191 & sn_snd_thick1, sn_snd_cond, sn_snd_mpnd , sn_snd_sstfrz, sn_snd_ttilyr 185 192 ! ! Received from the atmosphere 186 TYPE(FLD_C) :: sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_ tauw, sn_rcv_dqnsdt, sn_rcv_qsr, &193 TYPE(FLD_C) :: sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, & 187 194 & sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf, sn_rcv_ts_ice 188 195 TYPE(FLD_C) :: sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp, sn_rcv_icb, sn_rcv_isf 189 ! Send to waves196 ! ! Send to waves 190 197 TYPE(FLD_C) :: sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev 191 ! Received from waves192 TYPE(FLD_C) :: sn_rcv_hsig, sn_rcv_phioc, sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper, sn_rcv_wnum, sn_rcv_tauwoc,&193 sn_rcv_wdrag, sn_rcv_wfreq198 ! ! Received from waves 199 TYPE(FLD_C) :: sn_rcv_hsig, sn_rcv_phioc, sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper, sn_rcv_wnum, & 200 & sn_rcv_wstrf, sn_rcv_wdrag, sn_rcv_charn, sn_rcv_taw, sn_rcv_bhd, sn_rcv_tusd, sn_rcv_tvsd 194 201 ! ! Other namelist parameters 195 202 INTEGER :: nn_cplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data … … 274 281 & sn_snd_ifrac , sn_snd_crtw , sn_snd_wlev , sn_rcv_hsig , sn_rcv_phioc , & 275 282 & sn_rcv_w10m , sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr , & 276 & sn_rcv_sdrfx , sn_rcv_sdrfy , sn_rcv_wper , sn_rcv_wnum , sn_rcv_ tauwoc, &277 & sn_rcv_ wdrag , sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal ,&278 & sn_rcv_ iceflx, sn_rcv_co2 , sn_rcv_mslp ,&279 & sn_rcv_ic b , sn_rcv_isf , sn_rcv_wfreq, sn_rcv_tauw , &280 & sn_rcv_ts_ice 283 & sn_rcv_sdrfx , sn_rcv_sdrfy , sn_rcv_wper , sn_rcv_wnum , sn_rcv_wstrf , & 284 & sn_rcv_charn , sn_rcv_taw , sn_rcv_bhd , sn_rcv_tusd , sn_rcv_tvsd, & 285 & sn_rcv_wdrag , sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , & 286 & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_icb , sn_rcv_isf , sn_rcv_ts_ice 287 281 288 !!--------------------------------------------------------------------- 282 289 ! … … 319 326 WRITE(numout,*)' sea ice heat fluxes = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 320 327 WRITE(numout,*)' atm co2 = ', TRIM(sn_rcv_co2%cldes ), ' (', TRIM(sn_rcv_co2%clcat ), ')' 328 WRITE(numout,*)' Sea ice surface skin temperature= ', TRIM(sn_rcv_ts_ice%cldes), ' (', TRIM(sn_rcv_ts_ice%clcat), ')' 329 WRITE(numout,*)' surface waves:' 321 330 WRITE(numout,*)' significant wave heigth = ', TRIM(sn_rcv_hsig%cldes ), ' (', TRIM(sn_rcv_hsig%clcat ), ')' 322 331 WRITE(numout,*)' wave to oce energy flux = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')' … … 325 334 WRITE(numout,*)' Mean wave period = ', TRIM(sn_rcv_wper%cldes ), ' (', TRIM(sn_rcv_wper%clcat ), ')' 326 335 WRITE(numout,*)' Mean wave number = ', TRIM(sn_rcv_wnum%cldes ), ' (', TRIM(sn_rcv_wnum%clcat ), ')' 327 WRITE(numout,*)' Wave peak frequency = ', TRIM(sn_rcv_wfreq%cldes ), ' (', TRIM(sn_rcv_wfreq%clcat ), ')' 328 WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_tauwoc%cldes), ' (', TRIM(sn_rcv_tauwoc%clcat ), ')' 329 WRITE(numout,*)' Stress components by waves = ', TRIM(sn_rcv_tauw%cldes ), ' (', TRIM(sn_rcv_tauw%clcat ), ')' 336 WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')' 330 337 WRITE(numout,*)' Neutral surf drag coefficient = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 331 WRITE(numout,*)' Sea ice surface skin temperature= ', TRIM(sn_rcv_ts_ice%cldes), ' (', TRIM(sn_rcv_ts_ice%clcat), ')'338 WRITE(numout,*)' Charnock coefficient = ', TRIM(sn_rcv_charn%cldes ), ' (', TRIM(sn_rcv_charn%clcat ), ')' 332 339 WRITE(numout,*)' sent fields (multiple ice categories)' 333 340 WRITE(numout,*)' surface temperature = ', TRIM(sn_snd_temp%cldes ), ' (', TRIM(sn_snd_temp%clcat ), ')' … … 351 358 WRITE(numout,*)' - mesh = ', sn_snd_crtw%clvgrd 352 359 ENDIF 353 360 IF( lwp .AND. ln_wave) THEN ! control print 361 WRITE(numout,*)' surface waves:' 362 WRITE(numout,*)' Significant wave heigth = ', TRIM(sn_rcv_hsig%cldes ), ' (', TRIM(sn_rcv_hsig%clcat ), ')' 363 WRITE(numout,*)' Wave to oce energy flux = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')' 364 WRITE(numout,*)' Surface Stokes drift grid u = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')' 365 WRITE(numout,*)' Surface Stokes drift grid v = ', TRIM(sn_rcv_sdrfy%cldes ), ' (', TRIM(sn_rcv_sdrfy%clcat ), ')' 366 WRITE(numout,*)' Mean wave period = ', TRIM(sn_rcv_wper%cldes ), ' (', TRIM(sn_rcv_wper%clcat ), ')' 367 WRITE(numout,*)' Mean wave number = ', TRIM(sn_rcv_wnum%cldes ), ' (', TRIM(sn_rcv_wnum%clcat ), ')' 368 WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')' 369 WRITE(numout,*)' Neutral surf drag coefficient = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 370 WRITE(numout,*)' Charnock coefficient = ', TRIM(sn_rcv_charn%cldes ), ' (', TRIM(sn_rcv_charn%clcat ), ')' 371 WRITE(numout,*)' Transport associated to Stokes drift grid u = ', TRIM(sn_rcv_tusd%cldes ), ' (', TRIM(sn_rcv_tusd%clcat ), ')' 372 WRITE(numout,*)' Transport associated to Stokes drift grid v = ', TRIM(sn_rcv_tvsd%cldes ), ' (', TRIM(sn_rcv_tvsd%clcat ), ')' 373 WRITE(numout,*)' Bernouilli pressure head = ', TRIM(sn_rcv_bhd%cldes ), ' (', TRIM(sn_rcv_bhd%clcat ), ')' 374 WRITE(numout,*)'Wave to ocean momentum flux and Net wave-supported stress = ', TRIM(sn_rcv_taw%cldes ), ' (', TRIM(sn_rcv_taw%clcat ), ')' 375 WRITE(numout,*)' Surface current to waves = ', TRIM(sn_snd_crtw%cldes ), ' (', TRIM(sn_snd_crtw%clcat ), ')' 376 WRITE(numout,*)' - referential = ', sn_snd_crtw%clvref 377 WRITE(numout,*)' - orientation = ', sn_snd_crtw%clvor 378 WRITE(numout,*)' - mesh = ', sn_snd_crtw%clvgrd 379 ENDIF 354 380 ! ! allocate sbccpl arrays 355 381 IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) … … 629 655 cpl_wper = .TRUE. 630 656 ENDIF 631 srcv(jpr_wfreq)%clname = 'O_WFreq' ! wave peak frequency632 IF( TRIM(sn_rcv_wfreq%cldes ) == 'coupled' ) THEN633 srcv(jpr_wfreq)%laction = .TRUE.634 cpl_wfreq = .TRUE.635 ENDIF636 657 srcv(jpr_wnum)%clname = 'O_WNum' ! mean wave number 637 658 IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' ) THEN … … 639 660 cpl_wnum = .TRUE. 640 661 ENDIF 641 srcv(jpr_tauwoc)%clname = 'O_TauOce' ! stress fraction adsorbed by the wave 642 IF( TRIM(sn_rcv_tauwoc%cldes ) == 'coupled' ) THEN 643 srcv(jpr_tauwoc)%laction = .TRUE. 644 cpl_tauwoc = .TRUE. 645 ENDIF 646 srcv(jpr_tauwx)%clname = 'O_Tauwx' ! ocean stress from wave in the x direction 647 srcv(jpr_tauwy)%clname = 'O_Tauwy' ! ocean stress from wave in the y direction 648 IF( TRIM(sn_rcv_tauw%cldes ) == 'coupled' ) THEN 649 srcv(jpr_tauwx)%laction = .TRUE. 650 srcv(jpr_tauwy)%laction = .TRUE. 651 cpl_tauw = .TRUE. 662 srcv(jpr_wstrf)%clname = 'O_WStrf' ! stress fraction adsorbed by the wave 663 IF( TRIM(sn_rcv_wstrf%cldes ) == 'coupled' ) THEN 664 srcv(jpr_wstrf)%laction = .TRUE. 665 cpl_wstrf = .TRUE. 652 666 ENDIF 653 667 srcv(jpr_wdrag)%clname = 'O_WDrag' ! neutral surface drag coefficient … … 656 670 cpl_wdrag = .TRUE. 657 671 ENDIF 658 IF( srcv(jpr_tauwoc)%laction .AND. srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction ) & 659 CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 660 '(sn_rcv_tauwoc=coupled and sn_rcv_tauw=coupled)' ) 672 srcv(jpr_charn)%clname = 'O_Charn' ! Chranock coefficient 673 IF( TRIM(sn_rcv_charn%cldes ) == 'coupled' ) THEN 674 srcv(jpr_charn)%laction = .TRUE. 675 cpl_charn = .TRUE. 676 ENDIF 677 srcv(jpr_bhd)%clname = 'O_Bhd' ! Bernoulli head. waves' induced surface pressure 678 IF( TRIM(sn_rcv_bhd%cldes ) == 'coupled' ) THEN 679 srcv(jpr_bhd)%laction = .TRUE. 680 cpl_bhd = .TRUE. 681 ENDIF 682 srcv(jpr_tusd)%clname = 'O_Tusd' ! zonal stokes transport 683 IF( TRIM(sn_rcv_tusd%cldes ) == 'coupled' ) THEN 684 srcv(jpr_tusd)%laction = .TRUE. 685 cpl_tusd = .TRUE. 686 ENDIF 687 srcv(jpr_tvsd)%clname = 'O_Tvsd' ! meridional stokes tranmport 688 IF( TRIM(sn_rcv_tvsd%cldes ) == 'coupled' ) THEN 689 srcv(jpr_tvsd)%laction = .TRUE. 690 cpl_tvsd = .TRUE. 691 ENDIF 692 693 srcv(jpr_twox)%clname = 'O_Twox' ! wave to ocean momentum flux in the u direction 694 srcv(jpr_twoy)%clname = 'O_Twoy' ! wave to ocean momentum flux in the v direction 695 srcv(jpr_tawx)%clname = 'O_Tawx' ! Net wave-supported stress in the u direction 696 srcv(jpr_tawy)%clname = 'O_Tawy' ! Net wave-supported stress in the v direction 697 IF( TRIM(sn_rcv_taw%cldes ) == 'coupled' ) THEN 698 srcv(jpr_twox)%laction = .TRUE. 699 srcv(jpr_twoy)%laction = .TRUE. 700 srcv(jpr_tawx)%laction = .TRUE. 701 srcv(jpr_tawy)%laction = .TRUE. 702 cpl_taw = .TRUE. 703 ENDIF 661 704 ! 662 705 ! ! ------------------------------- ! … … 1058 1101 ! initialisation of the coupler ! 1059 1102 ! ================================ ! 1060 1061 1103 CALL cpl_define(jprcv, jpsnd, nn_cplmodel) 1062 1104 … … 1071 1113 ENDIF 1072 1114 xcplmask(:,:,0) = 1. - SUM( xcplmask(:,:,1:nn_cplmodel), dim = 3 ) 1115 ! 1073 1116 ! 1074 1117 END SUBROUTINE sbc_cpl_init … … 1146 1189 IF( ncpl_qsr_freq /= 0) ncpl_qsr_freq = 86400 / ncpl_qsr_freq ! used by top 1147 1190 1191 IF ( ln_wave .AND. nn_components == 0 ) THEN 1192 ncpl_qsr_freq = 1; 1193 WRITE(numout,*) 'ncpl_qsr_freq is set to 1 when coupling NEMO with wave (without SAS) ' 1194 ENDIF 1148 1195 ENDIF 1149 1196 ! … … 1320 1367 IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1) 1321 1368 ! 1322 ! ! ========================= !1323 ! ! Wave peak frequency !1324 ! ! ========================= !1325 IF( srcv(jpr_wfreq)%laction ) wfreq(:,:) = frcv(jpr_wfreq)%z3(:,:,1)1326 !1327 1369 ! ! ========================= ! 1328 1370 ! ! Vertical mixing Qiao ! … … 1331 1373 1332 1374 ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode 1333 IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. srcv(jpr_wper)%laction&1334 .OR. srcv(jpr_hsig)%laction .OR. srcv(jpr_wfreq)%laction)THEN1375 IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. & 1376 srcv(jpr_wper)%laction .OR. srcv(jpr_hsig)%laction ) THEN 1335 1377 CALL sbc_stokes( Kmm ) 1336 1378 ENDIF … … 1339 1381 ! ! Stress adsorbed by waves ! 1340 1382 ! ! ========================= ! 1341 IF( srcv(jpr_tauwoc)%laction .AND. ln_tauwoc ) tauoc_wave(:,:) = frcv(jpr_tauwoc)%z3(:,:,1) 1342 1343 ! ! ========================= ! 1344 ! ! Stress component by waves ! 1345 ! ! ========================= ! 1346 IF( srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction .AND. ln_tauw ) THEN 1347 tauw_x(:,:) = frcv(jpr_tauwx)%z3(:,:,1) 1348 tauw_y(:,:) = frcv(jpr_tauwy)%z3(:,:,1) 1349 ENDIF 1350 1383 IF( srcv(jpr_wstrf)%laction .AND. ln_tauoc ) tauoc_wave(:,:) = frcv(jpr_wstrf)%z3(:,:,1) 1384 ! 1351 1385 ! ! ========================= ! 1352 1386 ! ! Wave drag coefficient ! 1353 1387 ! ! ========================= ! 1354 1388 IF( srcv(jpr_wdrag)%laction .AND. ln_cdgw ) cdn_wave(:,:) = frcv(jpr_wdrag)%z3(:,:,1) 1355 1389 ! 1390 ! ! ========================= ! 1391 ! ! Chranock coefficient ! 1392 ! ! ========================= ! 1393 IF( srcv(jpr_charn)%laction .AND. ln_charn ) charn(:,:) = frcv(jpr_charn)%z3(:,:,1) 1394 ! 1395 ! ! ========================= ! 1396 ! ! net wave-supported stress ! 1397 ! ! ========================= ! 1398 IF( srcv(jpr_tawx)%laction .AND. ln_taw ) tawx(:,:) = frcv(jpr_tawx)%z3(:,:,1) 1399 IF( srcv(jpr_tawy)%laction .AND. ln_taw ) tawy(:,:) = frcv(jpr_tawy)%z3(:,:,1) 1400 ! 1401 ! ! ========================= ! 1402 ! !wave to ocean momentum flux! 1403 ! ! ========================= ! 1404 IF( srcv(jpr_twox)%laction .AND. ln_taw ) twox(:,:) = frcv(jpr_twox)%z3(:,:,1) 1405 IF( srcv(jpr_twoy)%laction .AND. ln_taw ) twoy(:,:) = frcv(jpr_twoy)%z3(:,:,1) 1406 ! 1407 ! ! ========================= ! 1408 ! ! wave TKE flux at sfc ! 1409 ! ! ========================= ! 1410 IF( srcv(jpr_phioc)%laction .AND. ln_phioc ) phioc(:,:) = frcv(jpr_phioc)%z3(:,:,1) 1411 ! 1412 ! ! ========================= ! 1413 ! ! Bernoulli head ! 1414 ! ! ========================= ! 1415 IF( srcv(jpr_bhd)%laction .AND. ln_bern_srfc ) bhd_wave(:,:) = frcv(jpr_bhd)%z3(:,:,1) 1416 ! 1417 ! ! ========================= ! 1418 ! ! Stokes transport u dir ! 1419 ! ! ========================= ! 1420 IF( srcv(jpr_tusd)%laction .AND. ln_breivikFV_2016 ) tusd(:,:) = frcv(jpr_tusd)%z3(:,:,1) 1421 ! 1422 ! ! ========================= ! 1423 ! ! Stokes transport v dir ! 1424 ! ! ========================= ! 1425 IF( srcv(jpr_tvsd)%laction .AND. ln_breivikFV_2016 ) tvsd(:,:) = frcv(jpr_tvsd)%z3(:,:,1) 1426 ! 1356 1427 ! Fields received by SAS when OASIS coupling 1357 1428 ! (arrays no more filled at sbcssm stage) -
NEMO/trunk/src/OCE/SBC/sbcmod.F90
r13970 r14007 16 16 !! 4.0 ! 2016-06 (L. Brodeau) new general bulk formulation 17 17 !! 4.0 ! 2019-03 (F. Lemarié & G. Samson) add ABL compatibility (ln_abl=TRUE) 18 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) modified wave forcing and coupling 18 19 !!---------------------------------------------------------------------- 19 20 … … 54 55 USE usrdef_sbc ! user defined: surface boundary condition 55 56 USE closea ! closed sea 57 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 56 58 ! 57 59 USE prtctl ! Print control (prt_ctl routine) … … 70 72 71 73 INTEGER :: nsbc ! type of surface boundary condition (deduced from namsbc informations) 72 74 !! * Substitutions 75 # include "do_loop_substitute.h90" 73 76 !!---------------------------------------------------------------------- 74 77 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 99 102 & nn_ice , ln_ice_embd, & 100 103 & ln_traqsr, ln_dm2dc , & 101 & ln_rnf , nn_fwb , ln_ssr , ln_apr_dyn, & 102 & ln_wave , ln_cdgw , ln_sdw , ln_tauwoc , ln_stcor , & 103 & ln_tauw , nn_lsm, nn_sdrift 104 & ln_rnf , nn_fwb , ln_ssr , ln_apr_dyn, & 105 & ln_wave , nn_lsm 104 106 !!---------------------------------------------------------------------- 105 107 ! … … 133 135 WRITE(numout,*) ' bulk formulation ln_blk = ', ln_blk 134 136 WRITE(numout,*) ' ABL formulation ln_abl = ', ln_abl 137 WRITE(numout,*) ' Surface wave (forced or coupled) ln_wave = ', ln_wave 135 138 WRITE(numout,*) ' Type of coupling (Ocean/Ice/Atmosphere) : ' 136 139 WRITE(numout,*) ' ocean-atmosphere coupled formulation ln_cpl = ', ln_cpl … … 150 153 WRITE(numout,*) ' runoff / runoff mouths ln_rnf = ', ln_rnf 151 154 WRITE(numout,*) ' nb of iterations if land-sea-mask applied nn_lsm = ', nn_lsm 152 WRITE(numout,*) ' surface wave ln_wave = ', ln_wave 153 WRITE(numout,*) ' Stokes drift corr. to vert. velocity ln_sdw = ', ln_sdw 154 WRITE(numout,*) ' vertical parametrization nn_sdrift = ', nn_sdrift 155 WRITE(numout,*) ' wave modified ocean stress ln_tauwoc = ', ln_tauwoc 156 WRITE(numout,*) ' wave modified ocean stress component ln_tauw = ', ln_tauw 157 WRITE(numout,*) ' Stokes coriolis term ln_stcor = ', ln_stcor 158 WRITE(numout,*) ' neutral drag coefficient (CORE,NCAR) ln_cdgw = ', ln_cdgw 159 ENDIF 160 ! 161 IF( .NOT.ln_wave ) THEN 162 ln_sdw = .false. ; ln_cdgw = .false. ; ln_tauwoc = .false. ; ln_tauw = .false. ; ln_stcor = .false. 163 ENDIF 164 IF( ln_sdw ) THEN 165 IF( .NOT.(nn_sdrift==jp_breivik_2014 .OR. nn_sdrift==jp_li_2017 .OR. nn_sdrift==jp_peakfr) ) & 166 CALL ctl_stop( 'The chosen nn_sdrift for Stokes drift vertical velocity must be 0, 1, or 2' ) 167 ENDIF 168 ll_st_bv2014 = ( nn_sdrift==jp_breivik_2014 ) 169 ll_st_li2017 = ( nn_sdrift==jp_li_2017 ) 170 ll_st_bv_li = ( ll_st_bv2014 .OR. ll_st_li2017 ) 171 ll_st_peakfr = ( nn_sdrift==jp_peakfr ) 172 IF( ln_tauwoc .AND. ln_tauw ) & 173 CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 174 '(ln_tauwoc=.true. and ln_tauw=.true.)' ) 175 IF( ln_tauwoc ) & 176 CALL ctl_warn( 'You are subtracting the wave stress to the ocean (ln_tauwoc=.true.)' ) 177 IF( ln_tauw ) & 178 CALL ctl_warn( 'The wave modified ocean stress components are used (ln_tauw=.true.) ', & 179 'This will override any other specification of the ocean stress' ) 155 ENDIF 180 156 ! 181 157 IF( .NOT.ln_usr ) THEN ! the model calendar needs some specificities (except in user defined case) … … 357 333 IF( nn_ice == 3 ) CALL cice_sbc_init( nsbc, Kbb, Kmm ) ! CICE initialization 358 334 ! 359 IF( ln_wave ) CALL sbc_wave_init ! surface wave initialisation 335 IF( ln_wave ) THEN 336 CALL sbc_wave_init ! surface wave initialisation 337 ELSE 338 IF(lwp) WRITE(numout,*) 339 IF(lwp) WRITE(numout,*) ' No surface waves : all wave related logical set to false' 340 ln_sdw = .false. 341 ln_stcor = .false. 342 ln_cdgw = .false. 343 ln_tauoc = .false. 344 ln_wave_test = .false. 345 ln_charn = .false. 346 ln_taw = .false. 347 ln_phioc = .false. 348 ln_bern_srfc = .false. 349 ln_breivikFV_2016 = .false. 350 ln_vortex_force = .false. 351 ln_stshear = .false. 352 ENDIF 360 353 ! 361 354 END SUBROUTINE sbc_init … … 380 373 INTEGER, INTENT(in) :: kt ! ocean time step 381 374 INTEGER, INTENT(in) :: Kbb, Kmm ! ocean time level indices 375 INTEGER :: jj, ji ! dummy loop argument 382 376 ! 383 377 LOGICAL :: ll_sas, ll_opa ! local logical … … 412 406 ! 413 407 IF( .NOT.ll_sas ) CALL sbc_ssm ( kt, Kbb, Kmm ) ! mean ocean sea surface variables (sst_m, sss_m, ssu_m, ssv_m) 414 IF( ln_wave ) CALL sbc_wave( kt, Kmm ) ! surface waves415 416 408 ! 417 409 ! !== sbc formulation ==! 418 410 ! 411 ! 419 412 SELECT CASE( nsbc ) ! Compute ocean surface boundary condition 420 413 ! ! (i.e. utau,vtau, qns, qsr, emp, sfx) … … 423 416 CASE( jp_blk ) 424 417 IF( ll_sas ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OPA-SAS coupling: SAS receiving fields from OPA 418 !!!!!!!!!!! ATTENTION:ln_wave is not only used for oasis coupling !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 419 IF( ln_wave ) THEN 420 IF ( lk_oasis ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OPA-wave coupling 421 CALL sbc_wave ( kt, Kmm ) 422 ENDIF 425 423 CALL sbc_blk ( kt ) ! bulk formulation for the ocean 426 424 ! … … 436 434 IF( ln_mixcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! forced-coupled mixed formulation after forcing 437 435 ! 438 IF ( ln_wave .AND. (ln_tauwoc .OR. ln_tauw) ) CALL sbc_wstress( ) ! Wind stress provided by waves 436 IF( ln_wave .AND. ln_tauoc ) THEN ! Wave stress reduction 437 DO_2D( 0, 0, 0, 0) 438 utau(ji,jj) = utau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji-1,jj) ) * 0.5_wp 439 vtau(ji,jj) = vtau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji,jj-1) ) * 0.5_wp 440 END_2D 441 ! 442 CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) 443 CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) 444 ! 445 taum(:,:) = taum(:,:)*tauoc_wave(:,:) 446 ! 447 IF( kt == nit000 ) CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.', & 448 & 'If not requested select ln_tauoc=.false.' ) 449 ! 450 ELSEIF( ln_wave .AND. ln_taw ) THEN ! Wave stress reduction 451 utau(:,:) = utau(:,:) - tawx(:,:) + twox(:,:) 452 vtau(:,:) = vtau(:,:) - tawy(:,:) + twoy(:,:) 453 CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) 454 CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) 455 ! 456 DO_2D( 0, 0, 0, 0) 457 taum(ji,jj) = sqrt((.5*(utau(ji-1,jj)+utau(ji,jj)))**2 + (.5*(vtau(ji,jj-1)+vtau(ji,jj)))**2) 458 END_2D 459 ! 460 IF( kt == nit000 ) CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.', & 461 & 'If not requested select ln_taw=.false.' ) 462 ! 463 ENDIF 464 CALL lbc_lnk( 'sbcmod', taum(:,:), 'T', 1. ) 439 465 ! 440 466 ! !== Misc. Options ==! -
NEMO/trunk/src/OCE/SBC/sbcwave.F90
r13546 r14007 9 9 !! - ! 2016-12 (G. Madec, E. Clementi) update Stoke drift computation 10 10 !! + add sbc_wave_ini routine 11 !! 4.2 ! 2020-12 (G. Madec, E. Clementi) updates, new Stoke drift computation 12 !! according to Couvelard et al.,2019 11 13 !!---------------------------------------------------------------------- 12 14 13 15 !!---------------------------------------------------------------------- 14 16 !! sbc_stokes : calculate 3D Stokes-drift velocities 15 !! sbc_wave : wave data from wave model in netcdf files17 !! sbc_wave : wave data from wave model: forced (netcdf files) or coupled mode 16 18 !! sbc_wave_init : initialisation fo surface waves 17 19 !!---------------------------------------------------------------------- 18 USE phycst ! physical constants 20 USE phycst ! physical constants 19 21 USE oce ! ocean variables 20 USE sbc_oce ! Surface boundary condition: ocean fields21 USE zdf_oce, ONLY : ln_zdfswm22 USE dom_oce ! ocean domain variables 23 USE sbc_oce ! Surface boundary condition: ocean fields 22 24 USE bdy_oce ! open boundary condition variables 23 25 USE domvvl ! domain: variable volume layers … … 26 28 USE in_out_manager ! I/O manager 27 29 USE lib_mpp ! distribued memory computing library 28 USE fldread 30 USE fldread ! read input fields 29 31 30 32 IMPLICIT NONE … … 32 34 33 35 PUBLIC sbc_stokes ! routine called in sbccpl 34 PUBLIC sbc_wstress ! routine called in sbcmod35 36 PUBLIC sbc_wave ! routine called in sbcmod 36 37 PUBLIC sbc_wave_init ! routine called in sbcmod 37 38 38 39 ! Variables checking if the wave parameters are coupled (if not, they are read from file) 39 LOGICAL, PUBLIC :: cpl_hsig = .FALSE. 40 LOGICAL, PUBLIC :: cpl_phioc = .FALSE. 41 LOGICAL, PUBLIC :: cpl_sdrftx = .FALSE. 42 LOGICAL, PUBLIC :: cpl_sdrfty = .FALSE. 43 LOGICAL, PUBLIC :: cpl_wper = .FALSE. 44 LOGICAL, PUBLIC :: cpl_wfreq = .FALSE. 45 LOGICAL, PUBLIC :: cpl_wnum = .FALSE. 46 LOGICAL, PUBLIC :: cpl_tauwoc = .FALSE. 47 LOGICAL, PUBLIC :: cpl_tauw = .FALSE. 48 LOGICAL, PUBLIC :: cpl_wdrag = .FALSE. 40 LOGICAL, PUBLIC :: cpl_hsig = .FALSE. 41 LOGICAL, PUBLIC :: cpl_phioc = .FALSE. 42 LOGICAL, PUBLIC :: cpl_sdrftx = .FALSE. 43 LOGICAL, PUBLIC :: cpl_sdrfty = .FALSE. 44 LOGICAL, PUBLIC :: cpl_wper = .FALSE. 45 LOGICAL, PUBLIC :: cpl_wnum = .FALSE. 46 LOGICAL, PUBLIC :: cpl_wstrf = .FALSE. 47 LOGICAL, PUBLIC :: cpl_wdrag = .FALSE. 48 LOGICAL, PUBLIC :: cpl_charn = .FALSE. 49 LOGICAL, PUBLIC :: cpl_taw = .FALSE. 50 LOGICAL, PUBLIC :: cpl_bhd = .FALSE. 51 LOGICAL, PUBLIC :: cpl_tusd = .FALSE. 52 LOGICAL, PUBLIC :: cpl_tvsd = .FALSE. 49 53 50 54 INTEGER :: jpfld ! number of files to read for stokes drift … … 53 57 INTEGER :: jp_hsw ! index of significant wave hight (m) at T-point 54 58 INTEGER :: jp_wmp ! index of mean wave period (s) at T-point 55 INTEGER :: jp_wfr ! index of wave peak frequency (1/s) at T-point56 59 57 60 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 58 61 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 59 62 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_wn ! structure of input fields (file informations, fields read) wave number for Qiao 60 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauwoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 61 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauw ! structure of input fields (file informations, fields read) ocean stress components from wave model 62 63 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: cdn_wave !: 64 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: hsw, wmp, wnum !: 65 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wfreq !: 66 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave !: 67 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauw_x, tauw_y !: 68 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d !: 69 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: div_sd !: barotropic stokes drift divergence 70 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ut0sd, vt0sd !: surface Stokes drift velocities at t-point 71 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd , vsd , wsd !: Stokes drift velocities at u-, v- & w-points, resp. 72 63 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 64 65 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: cdn_wave !: Neutral drag coefficient at t-point 66 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: hsw !: Significant Wave Height at t-point 67 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wmp !: Wave Mean Period at t-point 68 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wnum !: Wave Number at t-point 69 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave !: stress reduction factor at t-point 70 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d !: Surface Stokes Drift module at t-point 71 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: div_sd !: barotropic stokes drift divergence 72 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ut0sd, vt0sd !: surface Stokes drift velocities at t-point 73 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd, vsd, wsd !: Stokes drift velocities at u-, v- & w-points, resp.u 74 ! 75 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: charn !: charnock coefficient at t-point 76 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tawx !: Net wave-supported stress, u 77 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tawy !: Net wave-supported stress, v 78 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: twox !: wave-ocean momentum flux, u 79 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: twoy !: wave-ocean momentum flux, v 80 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wavex !: stress reduction factor at, u component 81 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wavey !: stress reduction factor at, v component 82 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: phioc !: tke flux from wave model 83 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: KZN2 !: Kz*N2 84 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: bhd_wave !: Bernoulli head. wave induce pression 85 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tusd, tvsd !: Stokes drift transport 86 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: ZMX !: Kz*N2 73 87 !! * Substitutions 74 88 # include "do_loop_substitute.h90" … … 88 102 !! 2014 (DOI: 10.1175/JPO-D-14-0020.1) 89 103 !! 90 !! ** Method : - Calculate Stokes transport speed 91 !! - Calculate horizontal divergence 92 !! - Integrate the horizontal divergenze from the bottom 93 !! ** action 104 !! ** Method : - Calculate the horizontal Stokes drift velocity (Breivik et al. 2014) 105 !! - Calculate its horizontal divergence 106 !! - Calculate the vertical Stokes drift velocity 107 !! - Calculate the barotropic Stokes drift divergence 108 !! 109 !! ** action : - tsd2d : module of the surface Stokes drift velocity 110 !! - usd, vsd, wsd : 3 components of the Stokes drift velocity 111 !! - div_sd : barotropic Stokes drift divergence 94 112 !!--------------------------------------------------------------------- 95 113 INTEGER, INTENT(in) :: Kmm ! ocean time level index 96 114 INTEGER :: jj, ji, jk ! dummy loop argument 97 115 INTEGER :: ik ! local integer 98 REAL(wp) :: ztransp, zfac, zsp0 99 REAL(wp) :: zdepth, zsqrt_depth, zexp_depth, z_two_thirds, zsqrtpi !sqrt of pi 100 REAL(wp) :: zbot_u, zbot_v, zkb_u, zkb_v, zke3_u, zke3_v, zda_u, zda_v 101 REAL(wp) :: zstokes_psi_u_bot, zstokes_psi_v_bot 102 REAL(wp) :: zdep_u, zdep_v, zkh_u, zkh_v 103 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: zk_t, zk_u, zk_v, zu0_sd, zv0_sd ! 2D workspace 104 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: zstokes_psi_u_top, zstokes_psi_v_top ! 2D workspace 105 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ze3divh ! 3D workspace 106 !!--------------------------------------------------------------------- 107 ! 108 ALLOCATE( ze3divh(jpi,jpj,jpkm1) ) ! jpkm1 -> avoid lbc_lnk on jpk that is not defined 116 REAL(wp) :: ztransp, zfac, ztemp, zsp0, zsqrt, zbreiv16_w 117 REAL(wp) :: zdep_u, zdep_v, zkh_u, zkh_v, zda_u, zda_v, sdtrp 118 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: zk_t, zk_u, zk_v, zu0_sd, zv0_sd ! 2D workspace 119 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ze3divh, zInt_w ! 3D workspace 120 !!--------------------------------------------------------------------- 121 ! 122 ALLOCATE( ze3divh(jpi,jpj,jpkm1) ) ! jpkm1 -> avoid lbc_lnk on jpk that is not defined 123 ALLOCATE( zInt_w(jpi,jpj,jpk) ) 109 124 ALLOCATE( zk_t(jpi,jpj), zk_u(jpi,jpj), zk_v(jpi,jpj), zu0_sd(jpi,jpj), zv0_sd(jpi,jpj) ) 125 zk_t (:,:) = 0._wp 126 zk_u (:,:) = 0._wp 127 zk_v (:,:) = 0._wp 128 zu0_sd (:,:) = 0._wp 129 zv0_sd (:,:) = 0._wp 130 ze3divh (:,:,:) = 0._wp 131 110 132 ! 111 133 ! select parameterization for the calculation of vertical Stokes drift 112 134 ! exp. wave number at t-point 113 IF( ll_st_bv_li ) THEN ! (Eq. (19) in Breivik et al. (2014) ) 135 IF( ln_breivikFV_2016 ) THEN 136 ! Assumptions : ut0sd and vt0sd are surface Stokes drift at T-points 137 ! sdtrp is the norm of Stokes transport 138 ! 139 zfac = 0.166666666667_wp 140 DO_2D( 1, 1, 1, 1 ) ! In the deep-water limit we have ke = ||ust0||/( 6 * ||transport|| ) 141 zsp0 = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj) ) !<-- norm of Surface Stokes drift 142 tsd2d(ji,jj) = zsp0 143 IF( cpl_tusd .AND. cpl_tvsd ) THEN !stokes transport is provided in coupled mode 144 sdtrp = SQRT( tusd(ji,jj)*tusd(ji,jj) + tvsd(ji,jj)*tvsd(ji,jj) ) !<-- norm of Surface Stokes drift transport 145 ELSE 146 ! Stokes drift transport estimated from Hs and Tmean 147 sdtrp = 2.0_wp * rpi / 16.0_wp * & 148 & hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 149 ENDIF 150 zk_t (ji,jj) = zfac * zsp0 / MAX ( sdtrp, 0.0000001_wp ) !<-- ke = ||ust0||/( 6 * ||transport|| ) 151 END_2D 152 !# define zInt_w ze3divh 153 DO_3D( 1, 1, 1, 1, 1, jpk ) ! Compute the primitive of Breivik 2016 function at W-points 154 zfac = - 2._wp * zk_t (ji,jj) * gdepw(ji,jj,jk,Kmm) !<-- zfac should be negative definite 155 ztemp = EXP ( zfac ) 156 zsqrt = SQRT( -zfac ) 157 zbreiv16_w = ztemp - SQRT(rpi)*zsqrt*ERFC(zsqrt) !Eq. 16 Breivik 2016 158 zInt_w(ji,jj,jk) = ztemp - 4._wp * zk_t (ji,jj) * gdepw(ji,jj,jk,Kmm) * zbreiv16_w 159 END_3D 160 ! 161 DO jk = 1, jpkm1 162 zfac = 0.166666666667_wp 163 DO_2D( 1, 1, 1, 1 ) !++ Compute the FV Breivik 2016 function at T-points 164 zsp0 = zfac / MAX(zk_t (ji,jj),0.0000001_wp) 165 ztemp = zInt_w(ji,jj,jk) - zInt_w(ji,jj,jk+1) 166 zu0_sd(ji,jj) = ut0sd(ji,jj) * zsp0 * ztemp * tmask(ji,jj,jk) 167 zv0_sd(ji,jj) = vt0sd(ji,jj) * zsp0 * ztemp * tmask(ji,jj,jk) 168 END_2D 169 DO_2D( 1, 0, 1, 0 ) ! ++ Interpolate at U/V points 170 zfac = 1.0_wp / e3u(ji ,jj,jk,Kmm) 171 usd(ji,jj,jk) = 0.5_wp * zfac * ( zu0_sd(ji,jj)+zu0_sd(ji+1,jj) ) * umask(ji,jj,jk) 172 zfac = 1.0_wp / e3v(ji ,jj,jk,Kmm) 173 vsd(ji,jj,jk) = 0.5_wp * zfac * ( zv0_sd(ji,jj)+zv0_sd(ji,jj+1) ) * vmask(ji,jj,jk) 174 END_2D 175 ENDDO 176 !# undef zInt_w 177 ! 178 ELSE 114 179 zfac = 2.0_wp * rpi / 16.0_wp 115 180 DO_2D( 1, 1, 1, 1 ) … … 128 193 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 129 194 END_2D 130 ELSE IF( ll_st_peakfr ) THEN ! peak wave number calculated from the peak frequency received by the wave model 131 DO_2D( 1, 1, 1, 1 ) 132 zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav 133 END_2D 134 DO_2D( 1, 0, 1, 0 ) 135 zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 136 zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 137 ! 138 zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 139 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 140 END_2D 141 ENDIF 142 ! 195 143 196 ! !== horizontal Stokes Drift 3D velocity ==! 144 IF( ll_st_bv2014 ) THEN 197 145 198 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 146 199 zdep_u = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji+1,jj,jk,Kmm) ) 147 200 zdep_v = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji,jj+1,jk,Kmm) ) 148 ! 201 ! 149 202 zkh_u = zk_u(ji,jj) * zdep_u ! k * depth 150 203 zkh_v = zk_v(ji,jj) * zdep_v … … 156 209 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 157 210 END_3D 158 ELSE IF( ll_st_li2017 .OR. ll_st_peakfr ) THEN159 ALLOCATE( zstokes_psi_u_top(jpi,jpj), zstokes_psi_v_top(jpi,jpj) )160 DO_2D( 1, 0, 1, 0 )161 zstokes_psi_u_top(ji,jj) = 0._wp162 zstokes_psi_v_top(ji,jj) = 0._wp163 END_2D164 zsqrtpi = SQRT(rpi)165 z_two_thirds = 2.0_wp / 3.0_wp166 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! exp. wave number & Stokes drift velocity at u- & v-points167 zbot_u = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji+1,jj,jk+1,Kmm) ) ! 2 * bottom depth168 zbot_v = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji,jj+1,jk+1,Kmm) ) ! 2 * bottom depth169 zkb_u = zk_u(ji,jj) * zbot_u ! 2 * k * bottom depth170 zkb_v = zk_v(ji,jj) * zbot_v ! 2 * k * bottom depth171 !172 zke3_u = MAX(1.e-8_wp, 2.0_wp * zk_u(ji,jj) * e3u(ji,jj,jk,Kmm)) ! 2k * thickness173 zke3_v = MAX(1.e-8_wp, 2.0_wp * zk_v(ji,jj) * e3v(ji,jj,jk,Kmm)) ! 2k * thickness174 175 ! Depth attenuation .... do u component first..176 zdepth = zkb_u177 zsqrt_depth = SQRT(zdepth)178 zexp_depth = EXP(-zdepth)179 zstokes_psi_u_bot = 1.0_wp - zexp_depth &180 & - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &181 & + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )182 zda_u = ( zstokes_psi_u_bot - zstokes_psi_u_top(ji,jj) ) / zke3_u183 zstokes_psi_u_top(ji,jj) = zstokes_psi_u_bot184 185 ! ... and then v component186 zdepth =zkb_v187 zsqrt_depth = SQRT(zdepth)188 zexp_depth = EXP(-zdepth)189 zstokes_psi_v_bot = 1.0_wp - zexp_depth &190 & - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &191 & + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )192 zda_v = ( zstokes_psi_v_bot - zstokes_psi_v_top(ji,jj) ) / zke3_v193 zstokes_psi_v_top(ji,jj) = zstokes_psi_v_bot194 !195 usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)196 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)197 END_3D198 DEALLOCATE( zstokes_psi_u_top, zstokes_psi_v_top )199 211 ENDIF 200 212 … … 235 247 CALL iom_put( "vstokes", vsd ) 236 248 CALL iom_put( "wstokes", wsd ) 237 !238 DEALLOCATE( ze3divh )249 ! ! 250 DEALLOCATE( ze3divh, zInt_w ) 239 251 DEALLOCATE( zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 240 252 ! 241 253 END SUBROUTINE sbc_stokes 242 243 244 SUBROUTINE sbc_wstress( ) 245 !!--------------------------------------------------------------------- 246 !! *** ROUTINE sbc_wstress *** 247 !! 248 !! ** Purpose : Updates the ocean momentum modified by waves 249 !! 250 !! ** Method : - Calculate u,v components of stress depending on stress 251 !! model 252 !! - Calculate the stress module 253 !! - The wind module is not modified by waves 254 !! ** action 255 !!--------------------------------------------------------------------- 256 INTEGER :: jj, ji ! dummy loop argument 257 ! 258 IF( ln_tauwoc ) THEN 259 utau(:,:) = utau(:,:)*tauoc_wave(:,:) 260 vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 261 taum(:,:) = taum(:,:)*tauoc_wave(:,:) 262 ENDIF 263 ! 264 IF( ln_tauw ) THEN 265 DO_2D( 1, 0, 1, 0 ) 266 ! Stress components at u- & v-points 267 utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) ) 268 vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) ) 269 ! 270 ! Stress module at t points 271 taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) ) 272 END_2D 273 CALL lbc_lnk_multi( 'sbcwave', utau(:,:), 'U', -1.0_wp , vtau(:,:), 'V', -1.0_wp , taum(:,:) , 'T', -1.0_wp ) 274 ENDIF 275 ! 276 END SUBROUTINE sbc_wstress 277 278 254 ! 255 ! 279 256 SUBROUTINE sbc_wave( kt, Kmm ) 280 257 !!--------------------------------------------------------------------- 281 258 !! *** ROUTINE sbc_wave *** 282 259 !! 283 !! ** Purpose : read wave parameters from wave model in netcdf files. 284 !! 285 !! ** Method : - Read namelist namsbc_wave 286 !! - Read Cd_n10 fields in netcdf files 287 !! - Read stokes drift 2d in netcdf files 288 !! - Read wave number in netcdf files 289 !! - Compute 3d stokes drift using Breivik et al.,2014 290 !! formulation 291 !! ** action 260 !! ** Purpose : read wave parameters from wave model in netcdf files 261 !! or from a coupled wave mdoel 262 !! 292 263 !!--------------------------------------------------------------------- 293 264 INTEGER, INTENT(in ) :: kt ! ocean time step 294 265 INTEGER, INTENT(in ) :: Kmm ! ocean time index 295 266 !!--------------------------------------------------------------------- 267 ! 268 IF( kt == nit000 .AND. lwp ) THEN 269 WRITE(numout,*) 270 WRITE(numout,*) 'sbc_wave : update the read waves fields' 271 WRITE(numout,*) '~~~~~~~~ ' 272 ENDIF 296 273 ! 297 274 IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN !== Neutral drag coefficient ==! … … 300 277 ENDIF 301 278 302 IF( ln_tauwoc .AND. .NOT. cpl_tauwoc ) THEN !== Wave induced stress ==! 303 CALL fld_read( kt, nn_fsbc, sf_tauwoc ) ! read wave norm stress from external forcing 304 tauoc_wave(:,:) = sf_tauwoc(1)%fnow(:,:,1) * tmask(:,:,1) 305 ENDIF 306 307 IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN !== Wave induced stress ==! 308 CALL fld_read( kt, nn_fsbc, sf_tauw ) ! read ocean stress components from external forcing (T grid) 309 tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1) * tmask(:,:,1) 310 tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1) * tmask(:,:,1) 311 ENDIF 312 313 IF( ln_sdw ) THEN !== Computation of the 3d Stokes Drift ==! 279 IF( ln_tauoc .AND. .NOT. cpl_wstrf ) THEN !== Wave induced stress ==! 280 CALL fld_read( kt, nn_fsbc, sf_tauoc ) ! read stress reduction factor due to wave from external forcing 281 tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) * tmask(:,:,1) 282 ELSEIF ( ln_taw .AND. cpl_taw ) THEN 283 IF (kt < 1) THEN ! The first fields gave by OASIS have very high erroneous values .... 284 twox(:,:)=0._wp 285 twoy(:,:)=0._wp 286 tawx(:,:)=0._wp 287 tawy(:,:)=0._wp 288 tauoc_wavex(:,:) = 1._wp 289 tauoc_wavey(:,:) = 1._wp 290 ELSE 291 tauoc_wavex(:,:) = abs(twox(:,:)/tawx(:,:)) 292 tauoc_wavey(:,:) = abs(twoy(:,:)/tawy(:,:)) 293 ENDIF 294 ENDIF 295 296 IF ( ln_phioc .and. cpl_phioc .and. kt == nit000 ) THEN 297 WRITE(numout,*) 298 WRITE(numout,*) 'sbc_wave : PHIOC from wave model' 299 WRITE(numout,*) '~~~~~~~~ ' 300 ENDIF 301 302 IF( ln_sdw .AND. .NOT. cpl_sdrftx) THEN !== Computation of the 3d Stokes Drift ==! 314 303 ! 315 304 IF( jpfld > 0 ) THEN ! Read from file only if the field is not coupled 316 305 CALL fld_read( kt, nn_fsbc, sf_sd ) ! read wave parameters from external forcing 306 ! ! NB: test case mode, not read as jpfld=0 317 307 IF( jp_hsw > 0 ) hsw (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) * tmask(:,:,1) ! significant wave height 318 308 IF( jp_wmp > 0 ) wmp (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) * tmask(:,:,1) ! wave mean period 319 IF( jp_wfr > 0 ) wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) * tmask(:,:,1) ! Peak wave frequency320 309 IF( jp_usd > 0 ) ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) * tmask(:,:,1) ! 2D zonal Stokes Drift at T point 321 310 IF( jp_vsd > 0 ) vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) * tmask(:,:,1) ! 2D meridional Stokes Drift at T point 322 311 ENDIF 323 312 ! 324 ! Read also wave number if needed, so that it is available in coupling routines 325 IF( ln_zdfswm .AND. .NOT.cpl_wnum ) THEN 326 CALL fld_read( kt, nn_fsbc, sf_wn ) ! read wave parameters from external forcing 327 wnum(:,:) = sf_wn(1)%fnow(:,:,1) * tmask(:,:,1) 328 ENDIF 329 330 ! Calculate only if required fields have been read 331 ! In coupled wave model-NEMO case the call is done after coupling 313 IF( jpfld == 4 .OR. ln_wave_test ) & 314 & CALL sbc_stokes( Kmm ) ! Calculate only if all required fields are read 315 ! ! or in wave test case 316 ! ! ! In coupled case the call is done after (in sbc_cpl) 317 ENDIF 332 318 ! 333 IF( ( ll_st_bv_li .AND. jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0 ) .OR. &334 & ( ll_st_peakfr .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0 ) ) CALL sbc_stokes( Kmm )335 !336 ENDIF337 !338 319 END SUBROUTINE sbc_wave 339 320 … … 343 324 !! *** ROUTINE sbc_wave_init *** 344 325 !! 345 !! ** Purpose : read wave parameters from wave model in netcdf files.326 !! ** Purpose : Initialisation fo surface waves 346 327 !! 347 328 !! ** Method : - Read namelist namsbc_wave 348 !! - Read Cd_n10 fields in netcdf files 349 !! - Read stokes drift 2d in netcdf files 350 !! - Read wave number in netcdf files 351 !! - Compute 3d stokes drift using Breivik et al.,2014 352 !! formulation 329 !! - create the structure used to read required wave fields 330 !! (its size depends on namelist options) 353 331 !! ** action 354 332 !!--------------------------------------------------------------------- … … 357 335 !! 358 336 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 359 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: slf_i , slf_j! array of namelist informations on the fields to read337 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: slf_i ! array of namelist informations on the fields to read 360 338 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, & 361 & sn_hsw, sn_wmp, sn_wfr, sn_wnum, & 362 & sn_tauwoc, sn_tauwx, sn_tauwy ! informations about the fields to be read 363 ! 364 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, & 365 sn_wnum, sn_tauwoc, sn_tauwx, sn_tauwy 366 !!--------------------------------------------------------------------- 339 & sn_hsw, sn_wmp, sn_wnum, sn_tauoc ! informations about the fields to be read 340 ! 341 NAMELIST/namsbc_wave/ cn_dir, sn_cdg, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wnum, sn_tauoc, & 342 & ln_cdgw, ln_sdw, ln_tauoc, ln_stcor, ln_charn, ln_taw, ln_phioc, & 343 & ln_wave_test, ln_bern_srfc, ln_breivikFV_2016, ln_vortex_force, ln_stshear 344 !!--------------------------------------------------------------------- 345 IF(lwp) THEN 346 WRITE(numout,*) 347 WRITE(numout,*) 'sbc_wave_init : surface waves in the system' 348 WRITE(numout,*) '~~~~~~~~~~~~~ ' 349 ENDIF 367 350 ! 368 351 READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 369 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist' 370 352 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist') 353 371 354 READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 372 355 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist' ) 373 356 IF(lwm) WRITE ( numond, namsbc_wave ) 374 357 ! 375 IF( ln_cdgw ) THEN 376 IF( .NOT. cpl_wdrag ) THEN 377 ALLOCATE( sf_cd(1), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 378 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 358 IF(lwp) THEN 359 WRITE(numout,*) ' Namelist namsbc_wave' 360 WRITE(numout,*) ' Stokes drift ln_sdw = ', ln_sdw 361 WRITE(numout,*) ' Breivik 2016 ln_breivikFV_2016 = ', ln_breivikFV_2016 362 WRITE(numout,*) ' Stokes Coriolis & tracer advection terms ln_stcor = ', ln_stcor 363 WRITE(numout,*) ' Vortex Force ln_vortex_force = ', ln_vortex_force 364 WRITE(numout,*) ' Bernouilli Head Pressure ln_bern_srfc = ', ln_bern_srfc 365 WRITE(numout,*) ' wave modified ocean stress ln_tauoc = ', ln_tauoc 366 WRITE(numout,*) ' neutral drag coefficient (CORE bulk only) ln_cdgw = ', ln_cdgw 367 WRITE(numout,*) ' charnock coefficient ln_charn = ', ln_charn 368 WRITE(numout,*) ' Stress modificated by wave ln_taw = ', ln_taw 369 WRITE(numout,*) ' TKE flux from wave ln_phioc = ', ln_phioc 370 WRITE(numout,*) ' Surface shear with Stokes drift ln_stshear = ', ln_stshear 371 WRITE(numout,*) ' Test with constant wave fields ln_wave_test = ', ln_wave_test 372 ENDIF 373 374 ! ! option check 375 IF( .NOT.( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor .OR. ln_charn) ) & 376 & CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_stcor=F') 377 IF( ln_cdgw .AND. ln_blk ) & 378 & CALL ctl_stop( 'drag coefficient read from wave model NOT available yet with aerobulk package') 379 IF( ln_stcor .AND. .NOT.ln_sdw ) & 380 & CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T') 381 382 ! !== Allocate wave arrays ==! 383 ALLOCATE( ut0sd (jpi,jpj) , vt0sd (jpi,jpj) ) 384 ALLOCATE( hsw (jpi,jpj) , wmp (jpi,jpj) ) 385 ALLOCATE( wnum (jpi,jpj) ) 386 ALLOCATE( tsd2d (jpi,jpj) , div_sd(jpi,jpj) , bhd_wave(jpi,jpj) ) 387 ALLOCATE( usd (jpi,jpj,jpk), vsd (jpi,jpj,jpk), wsd (jpi,jpj,jpk) ) 388 ALLOCATE( tusd (jpi,jpj) , tvsd (jpi,jpj) , ZMX (jpi,jpj,jpk) ) 389 usd (:,:,:) = 0._wp 390 vsd (:,:,:) = 0._wp 391 wsd (:,:,:) = 0._wp 392 hsw (:,:) = 0._wp 393 wmp (:,:) = 0._wp 394 ut0sd (:,:) = 0._wp 395 vt0sd (:,:) = 0._wp 396 tusd (:,:) = 0._wp 397 tvsd (:,:) = 0._wp 398 bhd_wave(:,:) = 0._wp 399 ZMX (:,:,:) = 0._wp 400 ! 401 IF( ln_wave_test ) THEN !== Wave TEST case ==! set uniform waves fields 402 jpfld = 0 ! No field read 403 ln_cdgw = .FALSE. ! No neutral wave drag input 404 ln_tauoc = .FALSE. ! No wave induced drag reduction factor 405 ut0sd(:,:) = 0.13_wp * tmask(:,:,1) ! m/s 406 vt0sd(:,:) = 0.00_wp ! m/s 407 hsw (:,:) = 2.80_wp ! meters 408 wmp (:,:) = 8.00_wp ! seconds 409 ! 410 ELSE !== create the structure associated with fields to be read ==! 411 IF( ln_cdgw ) THEN ! wave drag 412 IF( .NOT. cpl_wdrag ) THEN 413 ALLOCATE( sf_cd(1), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 414 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 415 ! 416 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 417 IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 418 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 419 ENDIF 420 ALLOCATE( cdn_wave(jpi,jpj) ) 421 cdn_wave(:,:) = 0._wp 422 ENDIF 423 IF( ln_charn ) THEN ! wave drag 424 IF( .NOT. cpl_charn ) THEN 425 CALL ctl_stop( 'STOP', 'Charnock based wind stress can be used in coupled mode only' ) 426 ENDIF 427 ALLOCATE( charn(jpi,jpj) ) 428 charn(:,:) = 0._wp 429 ENDIF 430 IF( ln_taw ) THEN ! wind stress 431 IF( .NOT. cpl_taw ) THEN 432 CALL ctl_stop( 'STOP', 'wind stress from wave model can be used in coupled mode only, use ln_cdgw instead' ) 433 ENDIF 434 ALLOCATE( tawx(jpi,jpj) ) 435 ALLOCATE( tawy(jpi,jpj) ) 436 ALLOCATE( twox(jpi,jpj) ) 437 ALLOCATE( twoy(jpi,jpj) ) 438 ALLOCATE( tauoc_wavex(jpi,jpj) ) 439 ALLOCATE( tauoc_wavey(jpi,jpj) ) 440 tawx(:,:) = 0._wp 441 tawy(:,:) = 0._wp 442 twox(:,:) = 0._wp 443 twoy(:,:) = 0._wp 444 tauoc_wavex(:,:) = 1._wp 445 tauoc_wavey(:,:) = 1._wp 446 ENDIF 447 448 IF( ln_phioc ) THEN ! TKE flux 449 IF( .NOT. cpl_phioc ) THEN 450 CALL ctl_stop( 'STOP', 'phioc can be used in coupled mode only' ) 451 ENDIF 452 ALLOCATE( phioc(jpi,jpj) ) 453 phioc(:,:) = 0._wp 454 ENDIF 455 456 IF( ln_tauoc ) THEN ! normalized wave stress into the ocean 457 IF( .NOT. cpl_wstrf ) THEN 458 ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc 459 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauoc structure' ) 460 ! 461 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 462 IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 463 CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' ) 464 ENDIF 465 ALLOCATE( tauoc_wave(jpi,jpj) ) 466 tauoc_wave(:,:) = 0._wp 467 ENDIF 468 469 IF( ln_sdw ) THEN ! Stokes drift 470 ! 1. Find out how many fields have to be read from file if not coupled 471 jpfld=0 472 jp_usd=0 ; jp_vsd=0 ; jp_hsw=0 ; jp_wmp=0 473 IF( .NOT. cpl_sdrftx ) THEN 474 jpfld = jpfld + 1 475 jp_usd = jpfld 476 ENDIF 477 IF( .NOT. cpl_sdrfty ) THEN 478 jpfld = jpfld + 1 479 jp_vsd = jpfld 480 ENDIF 481 IF( .NOT. cpl_hsig ) THEN 482 jpfld = jpfld + 1 483 jp_hsw = jpfld 484 ENDIF 485 IF( .NOT. cpl_wper ) THEN 486 jpfld = jpfld + 1 487 jp_wmp = jpfld 488 ENDIF 489 ! 2. Read from file only the non-coupled fields 490 IF( jpfld > 0 ) THEN 491 ALLOCATE( slf_i(jpfld) ) 492 IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd 493 IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd 494 IF( jp_hsw > 0 ) slf_i(jp_hsw) = sn_hsw 495 IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp 496 ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift 497 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 498 ! 499 DO ifpr= 1, jpfld 500 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 501 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 502 END DO 503 ! 504 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 505 ENDIF 379 506 ! 380 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 381 IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 382 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 383 ENDIF 384 ALLOCATE( cdn_wave(jpi,jpj) ) 385 ENDIF 386 387 IF( ln_tauwoc ) THEN 388 IF( .NOT. cpl_tauwoc ) THEN 389 ALLOCATE( sf_tauwoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauwoc 390 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 507 ! 3. Wave number (only needed for Qiao parametrisation, ln_zdfqiao=T) 508 IF( .NOT. cpl_wnum ) THEN 509 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 510 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wn structure' ) 511 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 512 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 513 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 514 ENDIF 391 515 ! 392 ALLOCATE( sf_tauwoc(1)%fnow(jpi,jpj,1) ) 393 IF( sn_tauwoc%ln_tint ) ALLOCATE( sf_tauwoc(1)%fdta(jpi,jpj,1,2) ) 394 CALL fld_fill( sf_tauwoc, (/ sn_tauwoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' ) 395 ENDIF 396 ALLOCATE( tauoc_wave(jpi,jpj) ) 397 ENDIF 398 399 IF( ln_tauw ) THEN 400 IF( .NOT. cpl_tauw ) THEN 401 ALLOCATE( sf_tauw(2), STAT=ierror ) !* allocate and fill sf_wave with sn_tauwx/y 402 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' ) 403 ! 404 ALLOCATE( slf_j(2) ) 405 slf_j(1) = sn_tauwx 406 slf_j(2) = sn_tauwy 407 ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1) ) 408 ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1) ) 409 IF( slf_j(1)%ln_tint ) ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) ) 410 IF( slf_j(2)%ln_tint ) ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) ) 411 CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 412 ENDIF 413 ALLOCATE( tauw_x(jpi,jpj) ) 414 ALLOCATE( tauw_y(jpi,jpj) ) 415 ENDIF 416 417 IF( ln_sdw ) THEN ! Find out how many fields have to be read from file if not coupled 418 jpfld=0 419 jp_usd=0 ; jp_vsd=0 ; jp_hsw=0 ; jp_wmp=0 ; jp_wfr=0 420 IF( .NOT. cpl_sdrftx ) THEN 421 jpfld = jpfld + 1 422 jp_usd = jpfld 423 ENDIF 424 IF( .NOT. cpl_sdrfty ) THEN 425 jpfld = jpfld + 1 426 jp_vsd = jpfld 427 ENDIF 428 IF( .NOT. cpl_hsig .AND. ll_st_bv_li ) THEN 429 jpfld = jpfld + 1 430 jp_hsw = jpfld 431 ENDIF 432 IF( .NOT. cpl_wper .AND. ll_st_bv_li ) THEN 433 jpfld = jpfld + 1 434 jp_wmp = jpfld 435 ENDIF 436 IF( .NOT. cpl_wfreq .AND. ll_st_peakfr ) THEN 437 jpfld = jpfld + 1 438 jp_wfr = jpfld 439 ENDIF 440 441 ! Read from file only the non-coupled fields 442 IF( jpfld > 0 ) THEN 443 ALLOCATE( slf_i(jpfld) ) 444 IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd 445 IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd 446 IF( jp_hsw > 0 ) slf_i(jp_hsw) = sn_hsw 447 IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp 448 IF( jp_wfr > 0 ) slf_i(jp_wfr) = sn_wfr 449 450 ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift 451 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 452 ! 453 DO ifpr= 1, jpfld 454 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 455 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 456 END DO 457 ! 458 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 459 ENDIF 460 ALLOCATE( usd (jpi,jpj,jpk), vsd (jpi,jpj,jpk), wsd(jpi,jpj,jpk) ) 461 ALLOCATE( hsw (jpi,jpj) , wmp (jpi,jpj) ) 462 ALLOCATE( wfreq(jpi,jpj) ) 463 ALLOCATE( ut0sd(jpi,jpj) , vt0sd(jpi,jpj) ) 464 ALLOCATE( div_sd(jpi,jpj) ) 465 ALLOCATE( tsd2d (jpi,jpj) ) 466 467 ut0sd(:,:) = 0._wp 468 vt0sd(:,:) = 0._wp 469 hsw(:,:) = 0._wp 470 wmp(:,:) = 0._wp 471 472 usd(:,:,:) = 0._wp 473 vsd(:,:,:) = 0._wp 474 wsd(:,:,:) = 0._wp 475 ! Wave number needed only if ln_zdfswm=T 476 IF( .NOT. cpl_wnum ) THEN 477 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 478 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wave structure' ) 479 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 480 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 481 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 482 ENDIF 483 ALLOCATE( wnum(jpi,jpj) ) 516 ENDIF 517 ! 484 518 ENDIF 485 519 !
Note: See TracChangeset
for help on using the changeset viewer.