- Timestamp:
- 2016-11-25T16:41:40+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
r7194 r7340 4 4 !! Wave module 5 5 !!====================================================================== 6 !! History : 3.3 ! 2011-09 ( Adani M) Original code: Drag Coefficient7 !! : 3.4 ! 2012-10 ( Adani M)Stokes Drift8 !! 3.6 ! 2014-09 ( Clementi E, Oddo P)New Stokes Drift Computation6 !! History : 3.3 ! 2011-09 (M. Adani) Original code: Drag Coefficient 7 !! : 3.4 ! 2012-10 (M. Adani) Stokes Drift 8 !! 3.6 ! 2014-09 (E. Clementi,P. Oddo) New Stokes Drift Computation 9 9 !!---------------------------------------------------------------------- 10 10 … … 13 13 !!---------------------------------------------------------------------- 14 14 USE oce ! 15 USE sbc_oce ! Surface boundary condition: ocean fields15 USE sbc_oce ! Surface boundary condition: ocean fields 16 16 USE bdy_oce ! 17 17 USE domvvl ! 18 !19 18 USE iom ! I/O manager library 20 19 USE in_out_manager ! I/O manager 21 20 USE lib_mpp ! distribued memory computing library 22 USE fldread ! read input fields21 USE fldread ! read input fields 23 22 USE wrk_nemo ! 24 23 USE phycst ! physical constants … … 29 28 PUBLIC sbc_wave ! routine called in sbcmod 30 29 31 INTEGER 32 INTEGER 33 INTEGER 34 INTEGER 35 INTEGER 30 INTEGER, PARAMETER :: jpfld = 4 ! number of files to read for stokes drift 31 INTEGER, PARAMETER :: jp_usd = 1 ! index of stokes drift (i-component) (m/s) at T-point 32 INTEGER, PARAMETER :: jp_vsd = 2 ! index of stokes drift (j-component) (m/s) at T-point 33 INTEGER, PARAMETER :: jp_swh = 3 ! index of significant wave hight (m) at T-point 34 INTEGER, PARAMETER :: jp_wmp = 4 ! index of mean wave period (s) at T-point 36 35 37 36 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient … … 43 42 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave 44 43 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d 45 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: usd2d, vsd2d 46 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: usd3d, vsd3d, wsd3d 47 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: usd3dt, vsd3dt 48 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: zusd2dt, zvsd2dt 44 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd3d, vsd3d, wsd3d 45 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd3dt , vsd3dt 49 46 50 47 !! * Substitutions … … 74 71 USE zdf_oce, ONLY : ln_zdfqiao 75 72 76 INTEGER, INTENT( in ) :: kt ! ocean time step 77 ! 78 INTEGER :: ierror ! return error code 79 INTEGER :: ifpr, jj,ji,jk 80 INTEGER :: ios ! Local integer output status for namelist read 81 ! 82 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 83 REAL(wp) :: ztransp, zsp0, zk, zus, zvs 84 REAL(wp), DIMENSION(jpi,jpj) :: zfac 85 REAL(wp), DIMENSION(:,:,:), POINTER :: ze3hdiv ! 3D workspace 86 TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist informations on the fields to read 87 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, & 88 & sn_swh, sn_wmp, sn_wnum, sn_tauoc ! informations about the fields to be read 89 !! 73 INTEGER, INTENT( in ) :: kt ! ocean time step 74 ! 75 INTEGER :: ierror ! return error code 76 INTEGER :: ifpr, jj,ji,jk ! dummy loop indice 77 INTEGER :: ios ! Local integer output status for namelist read 78 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 79 REAL(wp) :: ztransp, zfac, zsp0, zk, zus, zvs 80 REAL(wp), DIMENSION(jpi,jpj) :: zusd2dt, zvsd2dt 81 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3hdiv ! 3D workspace 82 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd ! informations about the fields to be read 83 TYPE(FLD_N) :: sn_swh, sn_wmp, sn_wnum, sn_tauoc ! " " " " " " " 84 TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist informations on the fields to read 90 85 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_swh, sn_wmp, sn_wnum, sn_tauoc 91 86 !!--------------------------------------------------------------------- … … 103 98 IF(lwm) WRITE ( numond, namsbc_wave ) 104 99 ! 105 IF 106 ALLOCATE( sf_cd(1), STAT=ierror ) ! *allocate and fill sf_wave with sn_cdg100 IF( ln_cdgw ) THEN 101 ALLOCATE( sf_cd(1), STAT=ierror ) ! allocate and fill sf_wave with sn_cdg 107 102 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 108 103 ! 109 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1))110 IF( sn_cdg%ln_tint ) 104 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 105 IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 111 106 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 112 107 ALLOCATE( cdn_wave(jpi,jpj) ) 113 cdn_wave(:,:) = 0.0 114 ENDIF 115 116 IF ( ln_tauoc ) THEN 117 ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc 108 ENDIF 109 110 IF( ln_tauoc ) THEN 111 ALLOCATE( sf_tauoc(1), STAT=ierror ) ! allocate and fill sf_wave with sn_tauoc 118 112 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 119 113 ! 120 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1))121 IF( sn_tauoc%ln_tint ) 114 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 115 IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 122 116 CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 123 117 ALLOCATE( tauoc_wave(jpi,jpj) ) 124 118 tauoc_wave(:,:) = 0.0 125 ENDIF126 127 IF 119 ENDIF 120 121 IF( ln_sdw ) THEN 128 122 slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; 129 123 slf_i(jp_swh) = sn_swh ; slf_i(jp_wmp) = sn_wmp; 130 ALLOCATE( sf_sd(jpfld), STAT=ierror ) ! *allocate and fill sf_sd with stokes drift124 ALLOCATE( sf_sd(jpfld), STAT=ierror ) ! allocate and fill sf_sd with stokes drift 131 125 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 132 126 ! 133 DO ifpr = 1, jpfld127 DO ifpr = 1, jpfld 134 128 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 135 129 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) … … 137 131 138 132 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 139 ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj) ) 140 ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 141 ALLOCATE( usd3dt(jpi,jpj,jpk),vsd3dt(jpi,jpj,jpk) ) 133 ALLOCATE( usd3dt(jpi,jpj,jpk), vsd3dt(jpi,jpj,jpk), wsd3d(jpi,jpj,jpk) ) 134 ALLOCATE( usd3d (jpi,jpj,jpk), vsd3d (jpi,jpj,jpk) ) 142 135 ALLOCATE( swh(jpi,jpj), wmp(jpi,jpj) ) 143 ALLOCATE( zusd2dt(jpi,jpj), zvsd2dt(jpi,jpj) ) 144 usd3d(:,:,:) = 0._wp ; usd2d(:,:) = 0._wp ; 145 vsd3d(:,:,:) = 0._wp ; vsd2d(:,:) = 0._wp ; 146 wsd3d(:,:,:) = 0._wp ; 147 usd3dt(:,:,:) = 0._wp ; vsd3dt(:,:,:) = 0._wp ; 148 swh (:,:) = 0._wp ; wmp (:,:) = 0._wp ; 149 IF ( ln_zdfqiao ) THEN !== Vertical mixing enhancement using Qiao,2010 ==! 136 usd3d(:,:,:) = 0._wp 137 vsd3d(:,:,:) = 0._wp 138 wsd3d(:,:,:) = 0._wp 139 IF( ln_zdfqiao ) THEN !== Vertical mixing enhancement using Qiao,2010 ==! 150 140 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 151 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' )152 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1))141 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 142 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 153 143 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 154 144 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 155 145 ALLOCATE( wnum(jpi,jpj),tsd2d(jpi,jpj) ) 156 wnum(:,:) = 0._wp ; tsd2d(:,:) = 0._wp157 146 ENDIF 158 147 ENDIF 159 148 ENDIF 160 149 ! 161 IF 162 CALL fld_read( kt, nn_fsbc, sf_cd ) ! read from external forcing150 IF( ln_cdgw ) THEN !== Neutral drag coefficient ==! 151 CALL fld_read( kt, nn_fsbc, sf_cd ) ! read from external forcing 163 152 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 164 153 ENDIF 165 154 166 IF 167 CALL fld_read( kt, nn_fsbc, sf_tauoc ) !*read wave norm stress from external forcing155 IF( ln_tauoc ) THEN !== Wave induced stress ==! 156 CALL fld_read( kt, nn_fsbc, sf_tauoc ) ! read wave norm stress from external forcing 168 157 tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 169 158 ENDIF 170 159 171 IF ( ln_sdw ) THEN !== Computation of the 3d Stokes Drift ==! 172 ! 173 CALL fld_read( kt, nn_fsbc, sf_sd ) !* read wave parameters from external forcing 160 IF( ln_sdw ) THEN !== Computation of the 3d Stokes Drift ==! 161 CALL fld_read( kt, nn_fsbc, sf_sd ) ! read wave parameters from external forcing 174 162 swh(:,:) = sf_sd(jp_swh)%fnow(:,:,1) ! significant wave height 175 163 wmp(:,:) = sf_sd(jp_wmp)%fnow(:,:,1) ! wave mean period … … 180 168 !(DOI: 10.1175/JPO-D-14-0020.1)==! 181 169 ! 182 CALL wrk_alloc( jpi,jpj,jpk, ze3hdiv )183 170 DO jk = 1, jpk 184 171 DO jj = 1, jpj 185 172 DO ji = 1, jpi 186 ! On T grid 187 ! Stokes transport speed estimated from Hs and Tmean 188 ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 189 ! Stokes surface speed 190 zsp0 = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2) 191 ! Wavenumber scale 192 zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 193 ! Depth attenuation 194 zfac(ji,jj) = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 195 END DO 196 END DO 197 ! 198 DO jj = 1, jpj 199 DO ji = 1, jpi 200 usd3dt(ji,jj,jk) = zfac(ji,jj) * zusd2dt (ji,jj) * tmask(ji,jj,jk) 201 vsd3dt(ji,jj,jk) = zfac(ji,jj) * zvsd2dt (ji,jj) * tmask(ji,jj,jk) 173 ! On T grid 174 ! Stokes transport speed estimated from Hs and Tmean 175 ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 176 ! Stokes surface speed 177 zsp0 = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2 ) 178 ! Wavenumber scale 179 zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 180 ! Depth attenuation 181 zfac = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 182 ! 183 usd3dt(ji,jj,jk) = zfac * zusd2dt(ji,jj) * tmask(ji,jj,jk) 184 vsd3dt(ji,jj,jk) = zfac * zvsd2dt(ji,jj) * tmask(ji,jj,jk) 202 185 END DO 203 186 END DO … … 207 190 DO jj = 1, jpjm1 208 191 DO ji = 1, jpim1 209 usd3d(ji,jj,jk) = 0.5 * umask(ji,jj,jk) *&210 & (usd3dt(ji,jj,jk) + usd3dt(ji+1,jj,jk))211 vsd3d(ji,jj,jk) = 0.5 * vmask(ji,jj,jk) *&212 & (vsd3dt(ji,jj,jk) + vsd3dt(ji,jj+1,jk))192 usd3d(ji,jj,jk) = 0.5 * umask(ji,jj,jk) * & 193 & ( usd3dt(ji,jj,jk) + usd3dt(ji+1,jj,jk) ) 194 vsd3d(ji,jj,jk) = 0.5 * vmask(ji,jj,jk) * & 195 & ( vsd3dt(ji,jj,jk) + vsd3dt(ji,jj+1,jk) ) 213 196 END DO 214 197 END DO … … 218 201 CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) 219 202 ! 220 DO jk = 1, jpkm1 ! e3t *Horizontal divergence221 DO jj = 2, jpj m1222 DO ji = fs_2, fs_jpim1 ! vector opt.223 ze3hdiv(ji,jj,jk) = ( e2u(ji ,jj) * fse3u_n(ji ,jj,jk) *usd3d(ji ,jj,jk) &224 & - e2u(ji-1,jj) * fse3u_n(ji-1,jj,jk) *usd3d(ji-1,jj,jk) &225 & + e1v(ji,jj ) * fse3v_n(ji,jj ,jk) *vsd3d(ji,jj ,jk) &226 & - e1v(ji,jj-1) * fse3v_n(ji,jj-1,jk) *vsd3d(ji,jj-1,jk) ) * r1_e1e2t(ji,jj)203 DO jk = 1, jpkm1 ! Horizontal divergence 204 DO jj = 2, jpj 205 DO ji = fs_2, jpi 206 ze3hdiv(ji,jj,jk) = ( e2u(ji ,jj) * usd3d(ji ,jj,jk) & 207 & - e2u(ji-1,jj) * usd3d(ji-1,jj,jk) & 208 & + e1v(ji,jj ) * vsd3d(ji,jj ,jk) & 209 & - e1v(ji,jj-1) * vsd3d(ji,jj-1,jk) ) * r1_e1e2t(ji,jj) 227 210 END DO 228 211 END DO 229 IF( .NOT. AGRIF_Root() ) THEN 230 IF( nbondi == 1 .OR. nbondi == 2 ) ze3hdiv(nlci-1, : ,jk) = 0._wp ! east 231 IF( nbondi == -1 .OR. nbondi == 2 ) ze3hdiv( 2 , : ,jk) = 0._wp ! west 232 IF( nbondj == 1 .OR. nbondj == 2 ) ze3hdiv( : ,nlcj-1,jk) = 0._wp ! north 233 IF( nbondj == -1 .OR. nbondj == 2 ) ze3hdiv( : , 2 ,jk) = 0._wp ! south 234 ENDIF 235 END DO 212 END DO 213 ! 214 IF( .NOT. AGRIF_Root() ) THEN 215 IF( nbondi == 1 .OR. nbondi == 2 ) ze3hdiv(nlci-1, : ,:) = 0._wp ! east 216 IF( nbondi == -1 .OR. nbondi == 2 ) ze3hdiv( 2 , : ,:) = 0._wp ! west 217 IF( nbondj == 1 .OR. nbondj == 2 ) ze3hdiv( : ,nlcj-1,:) = 0._wp ! north 218 IF( nbondj == -1 .OR. nbondj == 2 ) ze3hdiv( : , 2 ,:) = 0._wp ! south 219 ENDIF 220 ! 236 221 CALL lbc_lnk( ze3hdiv, 'T', 1. ) 237 222 ! 238 DO jk = jpkm1, 1, -1 ! *integrate from the bottom the e3t * hor. divergence239 wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - ze3hdiv(:,:,jk)223 DO jk = jpkm1, 1, -1 ! integrate from the bottom the e3t * hor. divergence 224 wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - fse3t_n(:,:,jk) * ze3hdiv(:,:,jk) 240 225 END DO 241 226 #if defined key_bdy … … 246 231 ENDIF 247 232 #endif 248 CALL wrk_dealloc( jpi,jpj,jpk, ze3hdiv )249 233 250 234 IF ( ln_zdfqiao ) THEN 251 CALL fld_read( kt, nn_fsbc, sf_wn ) ! *read wave parameters from external forcing235 CALL fld_read( kt, nn_fsbc, sf_wn ) ! read wave parameters from external forcing 252 236 wnum(:,:) = sf_wn(1)%fnow(:,:,1) 253 237 … … 256 240 DO jj = 1, jpj 257 241 DO ji = 1, jpi 258 tsd2d(ji,jj) = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2)242 tsd2d(ji,jj) = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2 ) 259 243 END DO 260 244 END DO
Note: See TracChangeset
for help on using the changeset viewer.
