Changeset 11277 for branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
- Timestamp:
- 2019-07-17T15:29:15+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
r8058 r11277 4 4 !! Wave module 5 5 !!====================================================================== 6 !! History : 3.3.1 ! 2011-09 (Adani M) Original code: Drag Coefficient 7 !! : 3.4 ! 2012-10 (Adani M) Stokes Drift 6 !! History : 3.3 ! 2011-09 (Adani M) Original code: Drag Coefficient 7 !! : 3.4 ! 2012-10 (Adani M) Stokes Drift 8 !! 3.6 ! 2014-09 (Clementi E, Oddo P)New Stokes Drift Computation 9 !! - ! 2016-12 (G. Madec, E. Clementi) update Stoke drift computation 10 !! + add sbc_wave_ini routine 8 11 !!---------------------------------------------------------------------- 9 USE iom ! I/O manager library 10 USE in_out_manager ! I/O manager 11 USE lib_mpp ! distribued memory computing library 12 USE fldread ! read input fields 13 USE oce 14 USE sbc_oce ! Surface boundary condition: ocean fields 15 USE domvvl 16 17 12 18 13 !!---------------------------------------------------------------------- 19 !! sbc_wave : read drag coefficient from wave model in netcdf files 14 !! sbc_stokes : calculate 3D Stokes-drift velocities 15 !! sbc_wave : wave data from wave model in netcdf files 16 !! sbc_wave_init : initialisation fo surface waves 20 17 !!---------------------------------------------------------------------- 18 USE oce ! ocean variables 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE bdy_oce ! open boundary condition variables 21 USE domvvl ! domain: variable volume layers 22 ! 23 USE iom ! I/O manager library 24 USE in_out_manager ! I/O manager 25 USE lib_mpp ! distribued memory computing library 26 USE fldread ! read input fields 27 USE wrk_nemo ! 28 USE phycst ! physical constants 21 29 22 30 IMPLICIT NONE 23 31 PRIVATE 24 32 25 PUBLIC sbc_wave ! routine called in sbc_blk_core or sbc_blk_mfs 33 PUBLIC sbc_stokes ! routine called in sbccpl 34 PUBLIC sbc_stress ! routine called in sbcmod 35 PUBLIC sbc_wave ! routine called in sbcmod 36 PUBLIC sbc_wave_init ! routine called in sbcmod 26 37 27 INTEGER , PARAMETER :: jpfld = 3 ! maximum number of files to read for srokes drift 28 INTEGER , PARAMETER :: jp_usd = 1 ! index of stokes drift (i-component) (m/s) at T-point 29 INTEGER , PARAMETER :: jp_vsd = 2 ! index of stokes drift (j-component) (m/s) at T-point 30 INTEGER , PARAMETER :: jp_wn = 3 ! index of wave number (1/m) at T-point 31 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 32 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 33 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: cdn_wave 34 REAL(wp),ALLOCATABLE,DIMENSION (:,:) :: usd2d,vsd2d,uwavenum,vwavenum 35 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:,:) :: usd3d,vsd3d,wsd3d 36 37 !! * Substitutions 38 # include "domzgr_substitute.h90" 38 ! 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_sdrft = .FALSE. 42 LOGICAL, PUBLIC :: cpl_wper = .FALSE. 43 LOGICAL, PUBLIC :: cpl_wfreq = .FALSE. 44 LOGICAL, PUBLIC :: cpl_wnum = .FALSE. 45 LOGICAL, PUBLIC :: cpl_tauoc = .FALSE. 46 LOGICAL, PUBLIC :: cpl_tauw = .FALSE. 47 LOGICAL, PUBLIC :: cpl_wdrag = .FALSE. 48 49 INTEGER :: nn_sdrift ! type of parameterization to calculate vertical Stokes drift 50 INTEGER, PARAMETER :: jp_breivik = 0 ! Breivik 2015: v_z=v_0*[exp(2*k*z)/(1-8*k*z)] 51 INTEGER, PARAMETER :: jp_phillips = 1 ! Phillips: v_z=v_o*[exp(2*k*z)-beta*sqrt(-2*k*pi*z)*erfc(sqrt(-2*k*z))] 52 INTEGER, PARAMETER :: jp_peakph = 2 ! Phillips using the peak wave number read from wave model instead of the inverse depth scale 53 54 INTEGER :: jpfld ! number of files to read for stokes drift 55 INTEGER :: jp_usd ! index of stokes drift (i-component) (m/s) at T-point 56 INTEGER :: jp_vsd ! index of stokes drift (j-component) (m/s) at T-point 57 INTEGER :: jp_hsw ! index of significant wave hight (m) at T-point 58 INTEGER :: jp_wmp ! index of mean wave period (s) at T-point 59 INTEGER :: jp_wfr ! index of wave peak frequency (s^-1) at T-point 60 61 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 62 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 63 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_wn ! structure of input fields (file informations, fields read) wave number for Qiao 64 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 65 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauw ! structure of input fields (file informations, fields read) ocean stress components from wave model 66 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_phioc ! structure of input fields (file informations, fields read) wave to ocean energy 67 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: cdn_wave !: 68 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: hsw, wmp, wnum !: 69 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wfreq !: 70 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: rn_crban !: Craig and Banner constant for surface breaking waves mixing 71 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave !: 72 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauw_x !: 73 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauw_y !: 74 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d !: 75 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: div_sd !: barotropic stokes drift divergence 76 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ut0sd, vt0sd !: surface Stokes drift velocities at t-point 77 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd , vsd , wsd !: Stokes drift velocities at u-, v- & w-points, resp. 78 79 # include "vectopt_loop_substitute.h90" 39 80 !!---------------------------------------------------------------------- 40 81 !! NEMO/OPA 4.0 , NEMO Consortium (2011) … … 44 85 CONTAINS 45 86 87 SUBROUTINE sbc_stokes( ) 88 !!--------------------------------------------------------------------- 89 !! *** ROUTINE sbc_stokes *** 90 !! 91 !! ** Purpose : compute the 3d Stokes Drift according to Breivik et al., 92 !! 2014 (DOI: 10.1175/JPO-D-14-0020.1) 93 !! 94 !! ** Method : - Calculate Stokes transport speed 95 !! - Calculate horizontal divergence 96 !! - Integrate the horizontal divergenze from the bottom 97 !! ** action 98 !!--------------------------------------------------------------------- 99 INTEGER :: jj, ji, jk ! dummy loop argument 100 INTEGER :: ik ! local integer 101 REAL(wp) :: ztransp, zfac, ztemp 102 REAL(wp) :: zdep_u, zdep_v, zkh_u, zkh_v, zda_u, zda_v 103 REAL(wp), DIMENSION(:,:) , POINTER :: zk_t, zk_u, zk_v, zu0_sd, zv0_sd ! 2D workspace 104 REAL(wp), DIMENSION(:,:,:), POINTER :: ze3divh ! 3D workspace 105 106 !!--------------------------------------------------------------------- 107 ! 108 109 CALL wrk_alloc( jpi,jpj,jpk, ze3divh ) 110 CALL wrk_alloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 111 ! 112 ! select parameterization for the calculation of vertical Stokes drift 113 ! exp. wave number at t-point 114 IF( nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips ) THEN ! (Eq. (19) in Breivick et al. (2014) ) 115 zfac = 2.0_wp * rpi / 16.0_wp 116 DO jj = 1, jpj 117 DO ji = 1, jpi 118 ! Stokes drift velocity estimated from Hs and Tmean 119 ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 120 ! Stokes surface speed 121 tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj)) 122 ! Wavenumber scale 123 zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp ) 124 END DO 125 END DO 126 DO jj = 1, jpjm1 ! exp. wave number & Stokes drift velocity at u- & v-points 127 DO ji = 1, jpim1 128 zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 129 zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 130 ! 131 zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 132 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 133 END DO 134 END DO 135 ELSE IF( nn_sdrift==jp_peakph ) THEN ! peak wave number calculated from the peak frequency received by the wave model 136 DO jj = 1, jpjm1 137 DO ji = 1, jpim1 138 zk_u(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji+1,jj)*wfreq(ji+1,jj) ) / grav 139 zk_v(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji,jj+1)*wfreq(ji,jj+1) ) / grav 140 ! 141 zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 142 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 143 END DO 144 END DO 145 ENDIF 146 ! 147 ! !== horizontal Stokes Drift 3D velocity ==! 148 IF( nn_sdrift==jp_breivik ) THEN 149 DO jk = 1, jpkm1 150 DO jj = 2, jpjm1 151 DO ji = 2, jpim1 152 zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 153 zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 154 ! 155 zkh_u = zk_u(ji,jj) * zdep_u ! k * depth 156 zkh_v = zk_v(ji,jj) * zdep_v 157 ! ! Depth attenuation 158 zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u ) 159 zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v ) 160 ! 161 usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 162 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 163 END DO 164 END DO 165 END DO 166 ELSE IF( nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph ) THEN 167 DO jk = 1, jpkm1 168 DO jj = 2, jpjm1 169 DO ji = 2, jpim1 170 zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 171 zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 172 ! 173 zkh_u = zk_u(ji,jj) * zdep_u ! k * depth 174 zkh_v = zk_v(ji,jj) * zdep_v 175 ! ! Depth attenuation 176 zda_u = EXP( -2.0_wp*zkh_u ) - SQRT(2.0_wp*rpi*zkh_u) * ERFC(SQRT(2.0_wp*zkh_u)) 177 zda_v = EXP( -2.0_wp*zkh_v ) - SQRT(2.0_wp*rpi*zkh_v) * ERFC(SQRT(2.0_wp*zkh_v)) 178 ! 179 usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 180 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 181 END DO 182 END DO 183 END DO 184 ENDIF 185 186 CALL lbc_lnk( usd(:,:,:), 'U', vsd(:,:,:), 'V', -1. ) 187 ! 188 ! !== vertical Stokes Drift 3D velocity ==! 189 ! 190 DO jk = 1, jpkm1 ! Horizontal e3*divergence 191 DO jj = 2, jpj 192 DO ji = fs_2, jpi 193 ze3divh(ji,jj,jk) = ( e2u(ji ,jj) * e3u_n(ji ,jj,jk) * usd(ji, jj,jk) & 194 & - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk) & 195 & + e1v(ji,jj ) * e3v_n(ji,jj ,jk) * vsd(ji,jj ,jk) & 196 & - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk) ) * r1_e12t(ji,jj) 197 END DO 198 END DO 199 END DO 200 ! 201 IF( .NOT. AGRIF_Root() ) THEN 202 IF( nbondi == 1 .OR. nbondi == 2 ) ze3divh(nlci-1, : ,:) = 0._wp ! east 203 IF( nbondi == -1 .OR. nbondi == 2 ) ze3divh( 2 , : ,:) = 0._wp ! west 204 IF( nbondj == 1 .OR. nbondj == 2 ) ze3divh( : ,nlcj-1,:) = 0._wp ! north 205 IF( nbondj == -1 .OR. nbondj == 2 ) ze3divh( : , 2 ,:) = 0._wp ! south 206 ENDIF 207 ! 208 CALL lbc_lnk( ze3divh, 'T', 1. ) 209 ! 210 IF( .NOT. lk_vvl ) THEN ; ik = 1 ! none zero velocity through the sea surface 211 ELSE ; ik = 2 ! w=0 at the surface (set one for all in sbc_wave_init) 212 ENDIF 213 DO jk = jpkm1, ik, -1 ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero) 214 wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk) 215 END DO 216 #if defined key_bdy 217 IF( lk_bdy ) THEN 218 DO jk = 1, jpkm1 219 wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:) 220 END DO 221 ENDIF 222 #endif 223 ! !== Horizontal divergence of barotropic Stokes transport ==! 224 div_sd(:,:) = 0._wp 225 DO jk = 1, jpkm1 ! 226 div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk) 227 END DO 228 ! 229 CALL iom_put( "ustokes", usd ) 230 CALL iom_put( "vstokes", vsd ) 231 CALL iom_put( "wstokes", wsd ) 232 ! 233 CALL wrk_dealloc( jpi,jpj,jpk, ze3divh ) 234 CALL wrk_dealloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 235 ! 236 END SUBROUTINE sbc_stokes 237 238 239 SUBROUTINE sbc_stress( ) 240 !!--------------------------------------------------------------------- 241 !! *** ROUTINE sbc_stress *** 242 !! 243 !! ** Purpose : Updates the ocean momentum modified by waves 244 !! 245 !! ** Method : - Calculate u,v components of stress depending on stress 246 !! model 247 !! - Calculate the stress module 248 !! - The wind module is not modified by waves 249 !! ** action 250 !!--------------------------------------------------------------------- 251 INTEGER :: jj, ji ! dummy loop argument 252 ! 253 IF( ln_tauoc ) THEN 254 utau(:,:) = utau(:,:)*tauoc_wave(:,:) 255 vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 256 taum(:,:) = taum(:,:)*tauoc_wave(:,:) 257 ENDIF 258 ! 259 IF( ln_tauw ) THEN 260 DO jj = 1, jpjm1 261 DO ji = 1, jpim1 262 ! Stress components at u- & v-points 263 utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) ) 264 vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) ) 265 ! 266 ! Stress module at t points 267 taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) ) 268 END DO 269 END DO 270 CALL lbc_lnk_multi( utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. , taum(:,: ), 'T', -1. ) 271 ENDIF 272 ! 273 END SUBROUTINE sbc_stress 274 275 46 276 SUBROUTINE sbc_wave( kt ) 47 277 !!--------------------------------------------------------------------- 48 !! *** ROUTINE sbc_ apr***49 !! 50 !! ** Purpose : read drag coefficientfrom wave model in netcdf files.278 !! *** ROUTINE sbc_wave *** 279 !! 280 !! ** Purpose : read wave parameters from wave model in netcdf files. 51 281 !! 52 282 !! ** Method : - Read namelist namsbc_wave 53 283 !! - Read Cd_n10 fields in netcdf files 54 284 !! - Read stokes drift 2d in netcdf files 55 !! - Read wave number in netcdf files 56 !! - Compute 3d stokes drift using monochromatic 57 !! ** action : 58 !! 59 !!--------------------------------------------------------------------- 60 USE oce, ONLY : un,vn,hdivn,rotn 61 USE divcur 62 USE wrk_nemo 63 #if defined key_bdy 64 USE bdy_oce, ONLY : bdytmask 65 #endif 66 INTEGER, INTENT( in ) :: kt ! ocean time step 67 INTEGER :: ierror ! return error code 68 INTEGER :: ifpr, jj,ji,jk 69 INTEGER :: ios ! Local integer output status for namelist read 70 REAL(wp),DIMENSION(:,:,:),POINTER :: udummy,vdummy,hdivdummy,rotdummy 71 REAL :: z2dt,z1_2dt 72 TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist informations on the fields to read 285 !! - Read wave number in netcdf files 286 !! - Compute 3d stokes drift using Breivik et al.,2014 287 !! formulation 288 !! ** action 289 !!--------------------------------------------------------------------- 290 INTEGER, INTENT(in ) :: kt ! ocean time step 291 !!--------------------------------------------------------------------- 292 ! 293 IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN !== Neutral drag coefficient ==! 294 CALL fld_read( kt, nn_fsbc, sf_cd ) ! read from external forcing 295 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 296 ! check that the drag coefficient contains proper information even if 297 ! the masks do not match - the momentum stress is not masked! 298 WHERE( cdn_wave < 0.0 ) cdn_wave = 1.5e-3 299 WHERE( cdn_wave > 1.0 ) cdn_wave = 1.5e-3 300 ENDIF 301 302 IF( ln_tauoc .AND. .NOT. cpl_tauoc ) THEN !== Wave induced stress ==! 303 CALL fld_read( kt, nn_fsbc, sf_tauoc ) ! read wave norm stress from external forcing 304 tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 305 WHERE( tauoc_wave < 0.0 ) tauoc_wave = 1.0 306 WHERE( tauoc_wave > 100.0 ) tauoc_wave = 1.0 307 ENDIF 308 309 IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN !== Wave induced stress ==! 310 CALL fld_read( kt, nn_fsbc, sf_tauw ) ! read ocean stress components from external forcing (T grid) 311 tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1) 312 WHERE( tauw_x < -100.0 ) tauw_x = 0.0 313 WHERE( tauw_x > 100.0 ) tauw_x = 0.0 314 315 tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1) 316 WHERE( tauw_y < -100.0 ) tauw_y = 0.0 317 WHERE( tauw_y > 100.0 ) tauw_y = 0.0 318 ENDIF 319 320 IF( ln_phioc .AND. .NOT. cpl_phioc ) THEN !== Wave to ocean energy ==! 321 CALL fld_read( kt, nn_fsbc, sf_phioc ) ! read wave to ocean energy from external forcing 322 rn_crban(:,:) = 29.0 * sf_phioc(1)%fnow(:,:,1) ! ! Alfa is phioc*sqrt(rau0/zrhoa) : rau0=water density, zhroa= air density 323 WHERE( rn_crban > 1.e8 ) rn_crban = 0.0 !remove first mask mistmatch points, then cap values in case of low friction velocity 324 WHERE( rn_crban < 0.0 ) rn_crban = 0.0 325 WHERE( rn_crban > 1000.0 ) rn_crban = 1000.0 326 ENDIF 327 328 IF( ln_sdw .OR. ln_rough ) THEN !== Computation of the 3d Stokes Drift ==! 329 ! 330 IF( jpfld > 0 ) THEN ! Read from file only if the field is not coupled 331 CALL fld_read( kt, nn_fsbc, sf_sd ) ! read wave parameters from external forcing 332 IF( jp_hsw > 0 ) THEN 333 hsw (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) ! significant wave height 334 WHERE( hsw > 100.0 ) hsw = 0.0 335 WHERE( hsw < 0.0 ) hsw = 0.0 336 ENDIF 337 IF( jp_wmp > 0 ) THEN 338 wmp (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) ! wave mean period 339 WHERE( wmp > 100.0 ) wmp = 0.0 340 WHERE( wmp < 0.0 ) wmp = 0.0 341 ENDIF 342 IF( jp_wfr > 0 ) THEN 343 wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) ! Peak wave frequency 344 WHERE( wfreq < 0.0 ) wfreq = 0.0 345 WHERE( wfreq > 100.0 ) wfreq = 0.0 346 ENDIF 347 IF( jp_usd > 0 ) THEN 348 ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) ! 2D zonal Stokes Drift at T point 349 WHERE( ut0sd < -100.0 ) ut0sd = 1.0 350 WHERE( ut0sd > 100.0 ) ut0sd = 1.0 351 ENDIF 352 IF( jp_vsd > 0 ) THEN 353 vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) ! 2D meridional Stokes Drift at T point 354 WHERE( vt0sd < -100.0 ) vt0sd = 1.0 355 WHERE( vt0sd > 100.0 ) vt0sd = 1.0 356 ENDIF 357 ENDIF 358 ENDIF 359 ! 360 IF( ln_sdw ) THEN 361 ! Read also wave number if needed, so that it is available in coupling routines 362 IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 363 CALL fld_read( kt, nn_fsbc, sf_wn ) ! read wave parameters from external forcing 364 wnum(:,:) = sf_wn(1)%fnow(:,:,1) 365 ENDIF 366 367 ! !== Computation of the 3d Stokes Drift ==! 368 ! 369 IF( ((nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) .AND. & 370 jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0) .OR. & 371 (nn_sdrift==jp_peakph .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0) ) & 372 CALL sbc_stokes() ! Calculate only if required fields are read 373 ! ! In coupled wave model-NEMO case the call is done after coupling 374 ! 375 ENDIF 376 ! 377 END SUBROUTINE sbc_wave 378 379 380 SUBROUTINE sbc_wave_init 381 !!--------------------------------------------------------------------- 382 !! *** ROUTINE sbc_wave_init *** 383 !! 384 !! ** Purpose : read wave parameters from wave model in netcdf files. 385 !! 386 !! ** Method : - Read namelist namsbc_wave 387 !! - Read Cd_n10 fields in netcdf files 388 !! - Read stokes drift 2d in netcdf files 389 !! - Read wave number in netcdf files 390 !! - Compute 3d stokes drift using Breivik et al.,2014 391 !! formulation 392 !! ** action 393 !!--------------------------------------------------------------------- 394 INTEGER :: ierror, ios ! local integer 395 INTEGER :: ifpr 396 !! 73 397 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 74 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_wn ! informations about the fields to be read 75 !!--------------------------------------------------------------------- 76 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_wn 77 !!--------------------------------------------------------------------- 78 79 !!---------------------------------------------------------------------- 80 ! 81 ! 82 ! ! -------------------- ! 83 IF( kt == nit000 ) THEN ! First call kt=nit000 ! 84 ! ! -------------------- ! 85 REWIND( numnam_ref ) ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 86 READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 87 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 88 89 REWIND( numnam_cfg ) ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 90 READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 91 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 92 IF(lwm) WRITE ( numond, namsbc_wave ) 93 ! 94 95 IF ( ln_cdgw ) THEN 398 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: slf_i, slf_j ! array of namelist informations on the fields to read 399 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_phioc, & 400 & sn_hsw, sn_wmp, sn_wfr, sn_wnum , & 401 & sn_tauoc, sn_tauwx, sn_tauwy ! information about the fields to be read 402 ! 403 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, sn_wnum, sn_phioc, & 404 sn_tauoc, sn_tauwx, sn_tauwy, & 405 ln_cdgw, ln_sdw, ln_stcor, ln_phioc, ln_tauoc, ln_tauw, ln_zdfqiao, ln_rough, & 406 nn_sdrift, nn_wmix 407 !!--------------------------------------------------------------------- 408 ! 409 REWIND( numnam_ref ) ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 410 READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 411 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 412 413 REWIND( numnam_cfg ) ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 414 READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 415 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 416 IF(lwm) WRITE ( numond, namsbc_wave ) 417 ! 418 IF(lwp) THEN !* Parameter print 419 WRITE(numout,*) 420 WRITE(numout,*) 'sbc_wave_init: wave physics' 421 WRITE(numout,*) '~~~~~~~~' 422 WRITE(numout,*) ' Namelist namsbc_wave : set wave physics parameters' 423 WRITE(numout,*) ' Stokes drift corr. to vert. velocity ln_sdw = ', ln_sdw 424 WRITE(numout,*) ' vertical parametrization nn_sdrift = ', nn_sdrift 425 WRITE(numout,*) ' Stokes coriolis term ln_stcor = ', ln_stcor 426 WRITE(numout,*) ' wave modified ocean stress ln_tauoc = ', ln_tauoc 427 WRITE(numout,*) ' wave modified ocean stress components ln_tauw = ', ln_tauw 428 WRITE(numout,*) ' wave to ocean energy ln_phioc = ', ln_phioc 429 WRITE(numout,*) ' vertical mixing parametrization nn_wmix = ', nn_wmix 430 WRITE(numout,*) ' neutral drag coefficient ln_cdgw = ', ln_cdgw 431 WRITE(numout,*) ' wave roughness length modification ln_rough = ', ln_rough 432 WRITE(numout,*) ' Qiao vertical mixing formulation ln_zdfqiao = ', ln_zdfqiao 433 ENDIF 434 435 IF ( ln_wave ) THEN 436 ! Activated wave physics but no wave physics components activated 437 IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc & 438 .OR. ln_rough .OR. ln_zdfqiao) ) THEN 439 CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_tauw=F, ln_stcor=F ', & 440 'ln_phioc=F, ln_rough=F, ln_zdfqiao=F' ) 441 ELSE 442 IF (ln_stcor .AND. .NOT. ln_sdw) & 443 CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T') 444 IF ( ln_cdgw .AND. .NOT.(nn_drag==jp_ukmo .OR. nn_drag==jp_std .OR. nn_drag==jp_const .OR. nn_drag==jp_mcore) ) & 445 CALL ctl_stop( 'The chosen nn_drag for momentum calculation must be 0, 1, 2, or 3') 446 IF ( ln_cdgw .AND. ln_blk_core .AND. nn_drag==0 ) & 447 CALL ctl_stop( 'The chosen nn_drag for momentum calculation in core forcing must be 1, 2, or 3') 448 IF ( ln_cdgw .AND. ln_flx .AND. nn_drag==3 ) & 449 CALL ctl_stop( 'The chosen nn_drag for momentum calculation in direct forcing must be 0, 1, or 2') 450 IF( ln_phioc .AND. .NOT.(nn_wmix==jp_craigbanner .OR. nn_wmix==jp_janssen) ) & 451 CALL ctl_stop( 'The chosen nn_wmix for wave vertical mixing must be 0, or 1' ) 452 IF( ln_sdw .AND. .NOT.(nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph) ) & 453 CALL ctl_stop( 'The chosen nn_sdrift for Stokes drift vertical velocity must be 0, 1, or 2' ) 454 IF( ln_zdfqiao .AND. .NOT.ln_sdw ) & 455 CALL ctl_stop( 'Qiao vertical mixing can not be used without Stokes drift (ln_sdw)' ) 456 IF( ln_tauoc .AND. ln_tauw ) & 457 CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 458 '(ln_tauoc=.true. and ln_tauw=.true.)' ) 459 IF( ln_tauoc ) & 460 CALL ctl_warn( 'You are subtracting the wave stress to the ocean (ln_tauoc=.true.)' ) 461 IF( ln_tauw ) & 462 CALL ctl_warn( 'The wave modified ocean stress components are used (ln_tauw=.true.) ', & 463 'This will override any other specification of the ocean stress' ) 464 ENDIF 465 ELSE 466 IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao ) & 467 & CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ', & 468 & 'with drag coefficient (ln_cdgw =T) ' , & 469 & 'or Stokes Drift (ln_sdw=T) ' , & 470 & 'or Stokes-Coriolis term (ln_stcor=T)', & 471 & 'or ocean stress modification due to waves (ln_tauoc=T) ', & 472 & 'or ocean stress components from waves (ln_tauw=T) ', & 473 & 'or wave to ocean energy modification (ln_phioc=T) ', & 474 & 'or wave surface roughness (ln_rough=T) ', & 475 & 'or Qiao vertical mixing formulation (ln_zdfqiao=T) ' ) 476 ENDIF 477 ! 478 IF( ln_cdgw ) THEN 479 IF( .NOT. cpl_wdrag ) THEN 96 480 ALLOCATE( sf_cd(1), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 97 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave : unable to allocate sf_wavestructure' )481 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_cd structure' ) 98 482 ! 99 483 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 100 484 IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 101 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 102 ALLOCATE( cdn_wave(jpi,jpj) ) 103 cdn_wave(:,:) = 0.0 104 ENDIF 105 IF ( ln_sdw ) THEN 106 slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 107 ALLOCATE( sf_sd(3), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 108 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 485 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 486 ENDIF 487 ALLOCATE( cdn_wave(jpi,jpj) ) 488 ENDIF 489 490 IF( ln_tauoc ) THEN 491 IF( .NOT. cpl_tauoc ) THEN 492 ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc 493 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauoc structure' ) 109 494 ! 110 DO ifpr= 1, jpfld 111 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 112 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 113 END DO 114 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 115 ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj),uwavenum(jpi,jpj),vwavenum(jpi,jpj) ) 116 ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 117 usd2d(:,:) = 0.0 ; vsd2d(:,:) = 0.0 ; uwavenum(:,:) = 0.0 ; vwavenum(:,:) = 0.0 118 usd3d(:,:,:) = 0.0 ;vsd3d(:,:,:) = 0.0 ; wsd3d(:,:,:) = 0.0 119 ENDIF 120 ENDIF 495 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 496 IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 497 CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 498 ENDIF 499 ALLOCATE( tauoc_wave(jpi,jpj) ) 500 ENDIF 501 502 IF( ln_tauw ) THEN 503 IF( .NOT. cpl_tauw ) THEN 504 ALLOCATE( sf_tauw(2), STAT=ierror ) !* allocate and fill sf_wave with sn_tauwx/y 505 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' ) 506 ! 507 ALLOCATE( slf_j(2) ) 508 slf_j(1) = sn_tauwx 509 slf_j(2) = sn_tauwy 510 ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1) ) 511 ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1) ) 512 IF( slf_j(1)%ln_tint ) ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) ) 513 IF( slf_j(2)%ln_tint ) ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) ) 514 CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 515 ENDIF 516 ALLOCATE( tauw_x(jpi,jpj) ) 517 ALLOCATE( tauw_y(jpi,jpj) ) 518 ENDIF 519 520 IF( ln_phioc ) THEN 521 IF( .NOT. cpl_phioc ) THEN 522 ALLOCATE( sf_phioc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_phioc 523 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_phioc structure' ) 524 ! 525 ALLOCATE( sf_phioc(1)%fnow(jpi,jpj,1) ) 526 IF( sn_phioc%ln_tint ) ALLOCATE( sf_phioc(1)%fdta(jpi,jpj,1,2) ) 527 CALL fld_fill( sf_phioc, (/ sn_phioc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 528 ENDIF 529 ALLOCATE( rn_crban(jpi,jpj) ) 530 ENDIF 531 532 ! Find out how many fields have to be read from file if not coupled 533 jpfld=0 534 jp_usd=0 ; jp_vsd=0 ; jp_hsw=0 ; jp_wmp=0 ; jp_wfr=0 535 IF( ln_sdw ) THEN 536 IF( .NOT. cpl_sdrft ) THEN 537 jpfld = jpfld + 1 538 jp_usd = jpfld 539 jpfld = jpfld + 1 540 jp_vsd = jpfld 541 ENDIF 542 IF( .NOT. cpl_hsig .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN 543 jpfld = jpfld + 1 544 jp_hsw = jpfld 545 ENDIF 546 IF( .NOT. cpl_wper .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN 547 jpfld = jpfld + 1 548 jp_wmp = jpfld 549 ENDIF 550 IF( .NOT. cpl_wfreq .AND. nn_sdrift==jp_peakph ) THEN 551 jpfld = jpfld + 1 552 jp_wfr = jpfld 553 ENDIF 554 ENDIF 555 556 IF( ln_rough .AND. .NOT. cpl_hsig .AND. jp_hsw==0 ) THEN 557 jpfld = jpfld + 1 558 jp_hsw = jpfld 559 ENDIF 560 561 ! Read from file only the non-coupled fields 562 IF( jpfld > 0 ) THEN 563 ALLOCATE( slf_i(jpfld) ) 564 IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd 565 IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd 566 IF( jp_hsw > 0 ) slf_i(jp_hsw) = sn_hsw 567 IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp 568 IF( jp_wfr > 0 ) slf_i(jp_wfr) = sn_wfr 569 ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift 570 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_sd structure' ) 121 571 ! 572 DO ifpr= 1, jpfld 573 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 574 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 575 END DO 122 576 ! 123 IF ( ln_cdgw ) THEN 124 CALL fld_read( kt, nn_fsbc, sf_cd ) !* read drag coefficient from external forcing 125 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 126 ENDIF 127 IF ( ln_sdw ) THEN 128 CALL fld_read( kt, nn_fsbc, sf_sd ) !* read drag coefficient from external forcing 129 130 ! Interpolate wavenumber, stokes drift into the grid_V and grid_V 131 !------------------------------------------------- 132 133 DO jj = 1, jpjm1 134 DO ji = 1, jpim1 135 uwavenum(ji,jj)=0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 136 & + sf_sd(3)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 137 138 vwavenum(ji,jj)=0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 139 & + sf_sd(3)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 140 141 usd2d(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(1)%fnow(ji,jj,1) * tmask(ji,jj,1) & 142 & + sf_sd(1)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 143 144 vsd2d(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(2)%fnow(ji,jj,1) * tmask(ji,jj,1) & 145 & + sf_sd(2)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 146 END DO 147 END DO 148 149 !Computation of the 3d Stokes Drift 150 DO jk = 1, jpk 151 DO jj = 1, jpj-1 152 DO ji = 1, jpi-1 153 usd3d(ji,jj,jk) = usd2d(ji,jj)*exp(2.0*uwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji+1,jj ,jk)))) 154 vsd3d(ji,jj,jk) = vsd2d(ji,jj)*exp(2.0*vwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji ,jj+1,jk)))) 155 END DO 156 END DO 157 usd3d(jpi,:,jk) = usd2d(jpi,:)*exp( 2.0*uwavenum(jpi,:)*(-gdept_0(jpi,:,jk)) ) 158 vsd3d(:,jpj,jk) = vsd2d(:,jpj)*exp( 2.0*vwavenum(:,jpj)*(-gdept_0(:,jpj,jk)) ) 159 END DO 160 161 CALL wrk_alloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 162 163 udummy(:,:,:)=un(:,:,:) 164 vdummy(:,:,:)=vn(:,:,:) 165 hdivdummy(:,:,:)=hdivn(:,:,:) 166 rotdummy(:,:,:)=rotn(:,:,:) 167 un(:,:,:)=usd3d(:,:,:) 168 vn(:,:,:)=vsd3d(:,:,:) 169 CALL div_cur(kt) 170 ! !------------------------------! 171 ! ! Now Vertical Velocity ! 172 ! !------------------------------! 173 z2dt = 2._wp * rdt ! set time step size (Euler/Leapfrog) 174 175 z1_2dt = 1.e0 / z2dt 176 DO jk = jpkm1, 1, -1 ! integrate from the bottom the hor. divergence 177 ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise 178 wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - fse3t_n(:,:,jk) * hdivn(:,:,jk) & 179 & - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) ) & 180 & * tmask(:,:,jk) * z1_2dt 181 #if defined key_bdy 182 wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 183 #endif 184 END DO 185 hdivn(:,:,:)=hdivdummy(:,:,:) 186 rotn(:,:,:)=rotdummy(:,:,:) 187 vn(:,:,:)=vdummy(:,:,:) 188 un(:,:,:)=udummy(:,:,:) 189 CALL wrk_dealloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 190 ENDIF 191 END SUBROUTINE sbc_wave 192 577 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 578 ENDIF 579 580 IF( ln_sdw ) THEN 581 ALLOCATE( usd (jpi,jpj,jpk), vsd (jpi,jpj,jpk), wsd(jpi,jpj,jpk) ) 582 ALLOCATE( wmp (jpi,jpj) ) 583 ALLOCATE( wfreq (jpi,jpj) ) 584 ALLOCATE( ut0sd(jpi,jpj) , vt0sd(jpi,jpj) ) 585 ALLOCATE( div_sd(jpi,jpj) ) 586 ALLOCATE( tsd2d (jpi,jpj) ) 587 usd(:,:,:) = 0._wp 588 vsd(:,:,:) = 0._wp 589 wsd(:,:,:) = 0._wp 590 ! Wave number needed only if ln_zdfqiao=T 591 IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 592 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 593 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wn structure' ) 594 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 595 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 596 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'read wave input', 'namsbc_wave' ) 597 ENDIF 598 ALLOCATE( wnum(jpi,jpj) ) 599 ENDIF 600 601 IF( ln_sdw .OR. ln_rough ) THEN 602 ALLOCATE( hsw (jpi,jpj) ) 603 ENDIF 604 ! 605 END SUBROUTINE sbc_wave_init 606 193 607 !!====================================================================== 194 608 END MODULE sbcwave
Note: See TracChangeset
for help on using the changeset viewer.