- Timestamp:
- 2015-12-02T17:12:45+01:00 (8 years ago)
- Location:
- branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM
- Files:
-
- 2 added
- 11 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/CONFIG/SHARED/namelist_ref
r5930 r5983 253 253 ln_cdgw = .false. ! Neutral drag coefficient read from wave model (T => fill namsbc_wave) 254 254 ln_sdw = .false. ! Computation of 3D stokes drift (T => fill namsbc_wave) 255 ln_tauoc= .false. ! Activate ocean stress modified by external wave induced stress (T => fill namsbc_wave) 256 ln_stcor= .false. ! Activate Stokes Coriolis term (T => fill namsbc_wave) 255 257 nn_lsm = 0 ! =0 land/sea mask for input fields is not applied (keep empty land/sea mask filename field) , 256 258 ! =1:n number of iterations of land/sea mask application for input fields (fill land/sea mask filename field) … … 915 917 ln_zdfexp = .false. ! time-stepping: split-explicit (T) or implicit (F) time stepping 916 918 nn_zdfexp = 3 ! number of sub-timestep for ln_zdfexp=T 919 ln_zdfqiao = .false. ! Enhanced wave vertical mixing Qiao (2010) 917 920 / 918 921 !----------------------------------------------------------------------- … … 1240 1243 sn_usd = 'sdw_wave' , 1 , 'u_sd2d' , .true. , .false. , 'daily' , '' , '' , '' 1241 1244 sn_vsd = 'sdw_wave' , 1 , 'v_sd2d' , .true. , .false. , 'daily' , '' , '' , '' 1242 sn_wn = 'sdw_wave' , 1 , 'wave_num' , .true. , .false. , 'daily' , '' , '' , '' 1245 sn_swh = 'sdw_wave' , 1 , 'hs' , .true. , .false. , 'daily' , '' , '' , '' 1246 sn_wmp = 'sdw_wave' , 1 , 'wmp' , .true. , .false. , 'daily' , '' , '' , '' 1247 sn_wnum = 'sdw_wave' , 1 , 'wave_num' , .true. , .false. , 'daily' , '' , '' , '' 1248 sn_tauoc - 'sdw_wave' , 1 , 'wave_stress', .true. , .false. , 'daily' , '' , '' , '' 1243 1249 ! 1244 1250 cn_dir_cdg = './' ! root directory for the location of drag coefficient files -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r5930 r5983 41 41 USE sbcapr ! surface boundary condition: atmospheric pressure 42 42 USE dynadv, ONLY: ln_dynadv_vec 43 USE sbcwave, ONLY: usd2d, vsd2d 43 44 #if defined key_agrif 44 45 USE agrif_opa_interp ! agrif … … 422 423 ENDIF 423 424 ! 425 ! Add Stokes Coriolis if defined 426 IF ( ln_stcor ) THEN 427 DO jj = 1, jpjm1 428 DO ji = 1, fs_jpim1 ! vector opt. 429 430 zy1 = ff(ji ,jj-1) * ( vsd2d(ji ,jj-1) + vsd2d(ji+1,jj-1) ) 431 zy2 = ff(ji ,jj ) * ( vsd2d(ji ,jj ) + vsd2d(ji+1,jj ) ) 432 zx1 = ff(ji-1,jj ) * ( usd2d(ji-1,jj ) + usd2d(ji-1,jj+1) ) 433 zx2 = ff(ji ,jj ) * ( usd2d(ji ,jj ) + usd2d(ji ,jj+1) ) 434 435 zu_frc(ji,jj) = zu_frc(ji,jj) + 0.25 * ( zy1 + zy2 ) * hur(ji,jj) 436 zv_frc(ji,jj) = zv_frc(ji,jj) - 0.25 * ( zx1 + zx2 ) * hvr(ji,jj) 437 END DO 438 END DO 439 ENDIF 440 ! 424 441 IF ( ln_apr_dyn ) THEN ! Add atm pressure forcing 425 442 IF (ln_bt_fw) THEN -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90
r5930 r5983 9 9 !! - ! 2010-09 (D.Storkey and E.O'Dea) bug fixes for BDY module 10 10 !! 3.3 ! 2011-10 (M. Leclair) split former ssh_wzv routine and remove all vvl related work 11 !! 3.6 ! 2014-10 (E. Clementi, P. Oddo) add wave contribution to surface vertical velocity 11 12 !!---------------------------------------------------------------------- 12 13 … … 38 39 USE wrk_nemo ! Memory Allocation 39 40 USE timing ! Timing 41 USE sbcwave, ONLY: usd2dt, vsd2dt,wsd3d 40 42 41 43 IMPLICIT NONE … … 160 162 REAL(wp), POINTER, DIMENSION(:,: ) :: z2d 161 163 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d, zhdiv 164 REAL(wp), ALLOCATABLE, DIMENSION(:,: ) :: sshnu, sshnv, dsshnu, dsshnv 162 165 !!---------------------------------------------------------------------- 163 166 ! … … 206 209 END DO 207 210 ENDIF 208 211 ! 212 ! In case ln_wave and ln_sdw the surface vertical velocity is modified 213 ! accounting for Sokes Drift velocity 214 IF ( ln_wave .AND. ln_sdw ) THEN 215 ALLOCATE(sshnu(jpi,jpj),sshnv(jpi,jpj),dsshnu(jpi,jpj),dsshnv(jpi,jpj)) 216 sshnu (:,:) = 0._wp 217 sshnv (:,:) = 0._wp 218 dsshnu(:,:) = 0._wp 219 dsshnv(:,:) = 0._wp 220 ! sshn interpolated on U-V grid 221 !--------------------------------- 222 DO jj = 1 , jpjm1 223 DO ji = 1 , jpim1 224 sshnu(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * & 225 & ( sshn(ji ,jj) * tmask(ji ,jj,1) & 226 & + sshn(ji+1,jj) * tmask(ji+1,jj,1) ) 227 sshnv(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * & 228 & ( sshn(ji,jj ) * tmask(ji,jj ,1) & 229 & + sshn(ji,jj+1) * tmask(ji,jj+1,1) ) 230 ENDDO 231 ENDDO 232 ! Compute d(ssh)/dx and d(ssh)/dy 233 !--------------------------------- 234 DO jj = 1 , jpjm1 235 DO ji = 1 , jpim1 236 dsshnu(ji,jj) = ( sshnu(ji+1,jj) - sshnu(ji,jj) ) / e1u(ji,jj) 237 dsshnv(ji,jj) = ( sshnv(ji,jj+1) - sshnv(ji,jj) ) / e2v(ji,jj) 238 ENDDO 239 ENDDO 240 ! Compute the surface vertical velocity accounting for the Stokes Drift 241 !--------------------------------------------------------------------- 242 wn(:,:,1) = wn(:,:,1) + usd2dt(:,:) * dsshnu(:,:) & 243 & + vsd2dt(:,:) * dsshnv(:,:) & 244 & - ( wsd3d (:,:,1) ) * tmask(:,:,1) 245 ENDIF 246 ! 209 247 #if defined key_bdy 210 248 IF( lk_bdy ) THEN -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90
r5836 r5983 65 65 LOGICAL , PUBLIC :: ln_cdgw !: true if neutral drag coefficient from wave model 66 66 LOGICAL , PUBLIC :: ln_sdw !: true if 3d stokes drift from wave model 67 LOGICAL , PUBLIC :: ln_tauoc !: true if normalized stress from wave is used 68 LOGICAL , PUBLIC :: ln_stcor !: true if Stokes-Coriolis term is used 67 69 ! 68 70 LOGICAL , PUBLIC :: ln_icebergs !: Icebergs -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r5836 r5983 88 88 & ln_blk_mfs, ln_apr_dyn, nn_ice, nn_ice_embd, ln_dm2dc , ln_rnf , & 89 89 & ln_ssr , nn_isf , nn_fwb, ln_cdgw , ln_wave , ln_sdw , & 90 & nn_lsm , nn_limflx , nn_components, ln_cpl90 & ln_tauoc , ln_stcor , nn_lsm, nn_limflx , nn_components, ln_cpl 91 91 INTEGER :: ios 92 92 INTEGER :: ierr, ierr0, ierr1, ierr2, ierr3, jpm … … 360 360 CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice ) ! OPA-SAS coupling: OPA receiving fields from SAS 361 361 END SELECT 362 362 IF (ln_wave .AND. ln_tauoc) THEN 363 utau(:,:) = utau(:,:)*tauoc_wave(:,:) 364 vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 365 taum(:,:) = taum(:,:)*tauoc_wave(:,:) 366 ! 367 SELECT CASE( nsbc ) 368 CASE( 0,1,2,3,5,-1 ) ; 369 IF(lwp) WRITE(numout,*) 'WARNING: You are subtracting the wave stress to the ocean. & 370 & If not requested select ln_tauoc=.false' 371 END SELECT 372 ! 373 END IF 363 374 IF( ln_mixcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice ) ! forced-coupled mixed formulation after forcing 364 375 -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
r5860 r5983 6 6 !! History : 3.3 ! 2011-09 (Adani M) Original code: Drag Coefficient 7 7 !! : 3.4 ! 2012-10 (Adani M) Stokes Drift 8 !! History : 3.6 !2014-09 (Clementi E, Oddo P)New Stokes Drift Computation 8 9 !!---------------------------------------------------------------------- 9 10 … … 12 13 !!---------------------------------------------------------------------- 13 14 USE oce ! 14 USE sbc_oce 15 USE sbc_oce ! Surface boundary condition: ocean fields 15 16 USE bdy_oce ! 16 17 USE domvvl ! … … 19 20 USE in_out_manager ! I/O manager 20 21 USE lib_mpp ! distribued memory computing library 21 USE fldread 22 USE fldread ! read input fields 22 23 USE wrk_nemo ! 24 USE phycst ! physical constants 23 25 24 26 IMPLICIT NONE … … 27 29 PUBLIC sbc_wave ! routine called in sbc_blk_core or sbc_blk_mfs 28 30 29 INTEGER , PARAMETER :: jpfld = 3 ! maximum number of files to read for srokes drift 30 INTEGER , PARAMETER :: jp_usd = 1 ! index of stokes drift (i-component) (m/s) at T-point 31 INTEGER , PARAMETER :: jp_vsd = 2 ! index of stokes drift (j-component) (m/s) at T-point 32 INTEGER , PARAMETER :: jp_wn = 3 ! index of wave number (1/m) at T-point 33 34 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 35 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 36 37 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:) :: cdn_wave 38 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:,:) :: usd3d, vsd3d, wsd3d 39 REAL(wp), ALLOCATABLE, DIMENSION (:,:) :: usd2d, vsd2d, uwavenum, vwavenum 31 INTEGER , PARAMETER :: jpfld = 4 ! number of files to read for stokes drift 32 INTEGER , PARAMETER :: jp_usd = 1 ! index of stokes drift (i-component) (m/s) at T-point 33 INTEGER , PARAMETER :: jp_vsd = 2 ! index of stokes drift (j-component) (m/s) at T-point 34 INTEGER , PARAMETER :: jp_swh = 3 ! index of significant wave hight (m) at T-point 35 INTEGER , PARAMETER :: jp_wmp = 4 ! index of mean wave period (s) at T-point 36 ! 37 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 38 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 39 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_wn ! structure of input fields (file informations, fields read) wave number for Qiao 40 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 41 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: cdn_wave 42 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: usd2d,vsd2d 43 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: swh,wmp,wnum 44 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: usd2dt,vsd2dt,tsd2d 45 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:,:) :: usd3d,vsd3d,wsd3d 46 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: tauoc_wave 40 47 41 48 !! * Substitutions … … 58 65 !! - Read Cd_n10 fields in netcdf files 59 66 !! - Read stokes drift 2d in netcdf files 60 !! - Read wave number in netcdf files 61 !! - Compute 3d stokes drift using monochromatic 67 !! - Read wave number in netcdf files 68 !! - Compute 3d stokes drift using Breivik et al.,2014 69 !! formulation 62 70 !! ** action : 63 71 !!--------------------------------------------------------------------- 72 USE zdf_oce, ONLY : ln_zdfqiao 73 64 74 INTEGER, INTENT( in ) :: kt ! ocean time step 65 75 ! 66 76 INTEGER :: ierror ! return error code 67 77 INTEGER :: ifpr, jj,ji,jk 68 INTEGER :: ios ! Local integer output status for namelist read 78 INTEGER :: ios ! Local integer output status for namelist read 79 80 REAL(wp) :: ztransp,zsp0, zk, zus,zvs 81 REAL(wp), DIMENSION(jpi,jpj) :: zfac 82 69 83 TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist informations on the fields to read 70 84 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 71 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_wn ! informations about the fields to be read 85 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, & 86 & sn_swh, sn_wmp, sn_wnum, sn_tauoc ! informations about the fields to be read 72 87 REAL(wp), DIMENSION(:,:,:), POINTER :: zusd_t, zvsd_t, ze3hdiv ! 3D workspace 73 88 !! 74 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_ wn89 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_swh, sn_wmp, sn_wnum, sn_tauoc 75 90 !!--------------------------------------------------------------------- 76 91 ! … … 97 112 cdn_wave(:,:) = 0.0 98 113 ENDIF 114 ! 115 IF ( ln_tauoc ) THEN 116 ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc 117 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 118 ! 119 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 120 IF( sn_cdg%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 121 CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 122 ALLOCATE( tauoc_wave(jpi,jpj) ) 123 tauoc_wave(:,:) = 0.0 124 ENDIF 125 ! 99 126 IF ( ln_sdw ) THEN 100 slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 101 ALLOCATE( sf_sd(3), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 127 slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; 128 slf_i(jp_swh) = sn_swh ; slf_i(jp_wmp) = sn_wmp; 129 ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift 102 130 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 103 131 ! … … 106 134 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 107 135 END DO 136 ! 108 137 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 109 ALLOCATE( usd2d(jpi,jpj) , vsd2d(jpi,jpj) , uwavenum(jpi,jpj) , vwavenum(jpi,jpj))138 ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj),usd2dt(jpi,jpj),vsd2dt(jpi,jpj)) 110 139 ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 111 usd3d(:,:,:) = 0._wp ; usd2d(:,:) = 0._wp ; uwavenum(:,:) = 0._wp 112 vsd3d(:,:,:) = 0._wp ; vsd2d(:,:) = 0._wp ; vwavenum(:,:) = 0._wp 113 wsd3d(:,:,:) = 0._wp 140 ALLOCATE( swh(jpi,jpj), wmp(jpi,jpj) ) 141 usd3d(:,:,:) = 0._wp ; usd2d(:,:) = 0._wp ; usd2dt(:,:) = 0._wp ; 142 vsd3d(:,:,:) = 0._wp ; vsd2d(:,:) = 0._wp ; vsd2dt(:,:) = 0._wp ; 143 wsd3d(:,:,:) = 0._wp ; 144 swh (:,:) = 0._wp ; wmp (:,:) = 0._wp ; 145 IF ( ln_zdfqiao ) THEN !== Vertical mixing enhancement using Qiao,2010 ==! 146 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 147 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable toallocate sf_wave structure' ) 148 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 149 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 150 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 151 ALLOCATE( wnum(jpi,jpj),tsd2d(jpi,jpj) ) 152 wnum(:,:) = 0._wp ; tsd2d(:,:) = 0._wp 153 ENDIF 114 154 ENDIF 115 155 ENDIF … … 119 159 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 120 160 ENDIF 121 ! 122 IF ( ln_sdw ) THEN !== Computation of the 3d Stokes Drift ==! 123 ! 124 CALL fld_read( kt, nn_fsbc, sf_sd ) !* read drag coefficient from external forcing 125 ! 126 ! 127 CALL wrk_alloc( jpi,jpj,jpk, zusd_t, zvsd_t, ze3hdiv ) 128 ! !* distribute it on the vertical 129 DO jk = 1, jpkm1 130 zusd_t(:,:,jk) = sf_sd(jp_usd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * fsdept_n(:,:,jk) ) 131 zvsd_t(:,:,jk) = sf_sd(jp_vsd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * fsdept_n(:,:,jk) ) 132 END DO 133 ! !* interpolate the stokes drift from t-point to u- and v-points 134 DO jk = 1, jpkm1 135 DO jj = 1, jpjm1 136 DO ji = 1, jpim1 137 usd3d(ji,jj,jk) = 0.5_wp * ( zusd_t(ji ,jj,jk) + zusd_t(ji+1,jj,jk) ) * umask(ji,jj,jk) 138 vsd3d(ji,jj,jk) = 0.5_wp * ( zvsd_t(ji ,jj,jk) + zvsd_t(ji,jj+1,jk) ) * vmask(ji,jj,jk) 139 END DO 140 END DO 141 END DO 161 ! 162 IF ( ln_tauoc ) THEN !== Wave induced stress ==! 163 CALL fld_read( kt, nn_fsbc, sf_tauoc ) !* read wave norm stress from external forcing 164 tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 165 ENDIF 166 ! 167 IF ( ln_sdw ) THEN !== Computation of the 3d Stokes Drift ==! 168 ! 169 CALL fld_read( kt, nn_fsbc, sf_sd ) !* read wave parameters from external forcing 170 swh(:,:) = sf_sd(jp_swh)%fnow(:,:,1) ! significant wave height 171 wmp(:,:) = sf_sd(jp_wmp)%fnow(:,:,1) ! wave mean period 172 usd2dt(:,:) = sf_sd(jp_usd)%fnow(:,:,1) ! 2D zonal Stokes Drift 173 vsd2dt(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) ! 2D meridional Stokes Drift 174 ! 175 !== Computation of the 3d Stokes Drift according to Breivik et al.,2014 176 !(DOI: 10.1175/JPO-D-14-0020.1)==! 177 ! 178 CALL wrk_alloc( jpi,jpj,jpk, ze3hdiv ) 179 DO jk = 1, jpk 180 DO jj = 1, jpj 181 DO ji = 1, jpi 182 ! On T grid 183 ! Stokes transport speed estimated from Hs and Tmean 184 ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 185 ! Stokes surface speed 186 zsp0 = SQRT( sf_sd(jp_usd)%fnow(ji,jj,1)**2 + sf_sd(jp_vsd)%fnow(ji,jj,1)**2) 187 ! Wavenumber scale 188 zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 189 ! Depth attenuation 190 zfac(ji,jj) = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 191 END DO 192 END DO 193 ! 194 DO jj = 1, jpjm1 195 DO ji = 1, jpim1 196 ! Into the U and V Grid 197 zus = 0.5 * ( 2. - umask(ji,jj,1) ) * ( zfac(ji,jj) * tmask(ji,jj,1) & 198 & + zfac(ji+1,jj) * tmask(ji+1,jj,1) ) 199 ! 200 zvs = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( zfac(ji,jj) * tmask(ji,jj,1) & 201 & + zfac(ji,jj+1) * tmask(ji,jj+1,1) ) 202 ! 203 usd2d(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( usd2dt(ji,jj) * tmask(ji,jj,1) & 204 & + usd2dt(ji+1,jj) * tmask(ji+1,jj,1) ) 205 ! 206 vsd2d(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( vsd2dt(ji,jj) * tmask(ji,jj,1) & 207 & + vsd2dt(ji,jj+1) * tmask(ji,jj+1,1) ) 208 ! 209 usd3d(ji,jj,jk) = usd2d(ji,jj)*zus 210 vsd3d(ji,jj,jk) = vsd2d(ji,jj)*zvs 211 END DO 212 END DO 213 END DO 214 ! 142 215 CALL lbc_lnk( usd3d(:,:,:), 'U', -1. ) 143 216 CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) … … 171 244 ENDIF 172 245 #endif 173 CALL wrk_dealloc( jpi,jpj,jpk, zusd_t, zvsd_t, ze3hdiv ) 174 ! 246 ! CALL wrk_dealloc( jpi,jpj,jpk, zusd_t, zvsd_t, ze3hdiv ) 247 CALL wrk_dealloc( jpi,jpj,jpk, ze3hdiv ) 248 ! 249 IF ( ln_zdfqiao ) THEN 250 wnum(:,:) = sf_wn(1)%fnow(:,:,1) 251 ! Calculate the module of the stokes drift on T grid 252 !------------------------------------------------- 253 DO jj = 1, jpj 254 DO ji = 1, jpi 255 tsd2d(ji,jj) = ((sf_sd(jp_usd)%fnow(ji,jj,1) * tmask(ji,jj,1))**2.0 + & 256 & (sf_sd(jp_vsd)%fnow(ji,jj,1) * tmask(ji,jj,1))**2.0)**0.5 257 END DO 258 END DO 259 ENDIF 260 ! 175 261 ENDIF 176 262 ! -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90
r5930 r5983 9 9 !! 3.7 ! 2014-05 (G. Madec) Add 2nd/4th order cases for CEN and FCT schemes 10 10 !! - ! 2014-12 (G. Madec) suppression of cross land advection option 11 !! 3.6 ! 2015-06 (E. Clementi) Addition of Stokes drift in case of wave coupling 11 12 !!---------------------------------------------------------------------- 12 13 … … 34 35 USE wrk_nemo ! Memory Allocation 35 36 USE timing ! Timing 36 37 USE diaptr ! Poleward heat transport 37 USE sbcwave ! wave module 38 USE sbc_oce ! surface boundary condition: ocean 39 40 USE diaptr ! Poleward heat transport 41 USE sbcwave ! wave module 42 USE sbc_oce ! surface boundary condition: ocean 38 43 39 44 IMPLICIT NONE … … 95 100 ! 96 101 ! ! set time step 102 zun(:,:,:) = 0.0 103 zvn(:,:,:) = 0.0 104 zwn(:,:,:) = 0.0 105 ! 97 106 IF( neuler == 0 .AND. kt == nit000 ) THEN ! at nit000 98 107 r2dtra(:) = rdttra(:) ! = rdtra (restarting with Euler time stepping) … … 102 111 ! 103 112 ! !== effective transport ==! 113 IF (ln_wave .AND. ln_sdw) THEN 114 DO jk = 1, jpkm1 115 zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) * & 116 & ( un(:,:,jk) + usd3d(:,:,jk) ) !eulerian transport + Stokes Drift 117 zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) * & 118 & ( vn(:,:,jk) + vsd3d(:,:,jk) ) 119 zwn(:,:,jk) = e1e2t(:,:) * & 120 & ( wn(:,:,jk) + wsd3d(:,:,jk) ) 121 END DO 122 ELSE 104 123 DO jk = 1, jpkm1 105 124 zun(:,:,jk) = e2u (:,:) * fse3u(:,:,jk) * un(:,:,jk) ! eulerian transport only … … 107 126 zwn(:,:,jk) = e1e2t(:,:) * wn(:,:,jk) 108 127 END DO 128 ENDIF 109 129 ! 110 130 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! add z-tilde and/or vvl corrections -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90
r5836 r5983 35 35 INTEGER , PUBLIC :: nn_npc !: non penetrative convective scheme call frequency 36 36 INTEGER , PUBLIC :: nn_npcp !: non penetrative convective scheme print frequency 37 LOGICAL , PUBLIC :: ln_zdfqiao !: Enhanced wave vertical mixing Qiao(2010) formulation flag 37 38 38 39 -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90
r5836 r5983 51 51 INTEGER :: ioptio, ios ! local integers 52 52 !! 53 NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp, & 54 & ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp 53 NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp, & 54 & ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp, & 55 & ln_zdfqiao 55 56 !!---------------------------------------------------------------------- 56 57 … … 81 82 WRITE(numout,*) ' npc call frequency nn_npc = ', nn_npc 82 83 WRITE(numout,*) ' npc print frequency nn_npcp = ', nn_npcp 84 WRITE(numout,*) ' Qiao formulation flag ln_zdfqiao=', ln_zdfqiao 83 85 ENDIF 84 86 -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/step.F90
r5930 r5983 26 26 !! 3.6 ! 2012-07 (J. Simeon, G. Madec. C. Ethe) Online coarsening of outputs 27 27 !! 3.6 ! 2014-04 (F. Roquet, G. Madec) New equations of state 28 !! 3.6 ! 2014-10 (E. Clementi, P. Oddo) Add Qiao vertical mixing in case of waves 28 29 !! 3.7 ! 2014-10 (G. Madec) LDF simplication 29 30 !! - ! 2014-12 (G. Madec) remove KPP scheme … … 75 76 !! -8- Outputs and diagnostics 76 77 !!---------------------------------------------------------------------- 77 INTEGER :: j k! dummy loop indice78 INTEGER :: ji,jj,jk ! dummy loop indice 78 79 INTEGER :: indic ! error indicator if < 0 79 80 INTEGER :: kcall ! optional integer argument (dom_vvl_sf_nxt) … … 132 133 IF( lk_zdftke ) CALL zdf_tke( kstp ) ! TKE closure scheme for Kz 133 134 IF( lk_zdfgls ) CALL zdf_gls( kstp ) ! GLS closure scheme for Kz 135 IF( ln_zdfqiao ) THEN 136 CALL zdf_qiao(kstp ) ! Qiao vertical mixing 137 DO jk = 1, jpkm1 138 DO jj = 1, jpj 139 DO ji = 1, jpi 140 avmu(ji,jj,jk) = (avmu(ji,jj,jk) + QBvu(ji,jj,jk)) * umask(ji,jj,jk) 141 avmv(ji,jj,jk) = (avmv(ji,jj,jk) + QBvv(ji,jj,jk)) * vmask(ji,jj,jk) 142 avt( ji,jj,jk) = (avt( ji,jj,jk) + QBv(ji,jj,jk)) * tmask(ji,jj,jk) 143 END DO 144 END DO 145 END DO 146 ENDIF 147 ! 134 148 IF( lk_zdfcst ) THEN ! Constant Kz (reset avt, avm[uv] to the background value) 135 149 avt (:,:,:) = rn_avt0 * wmask (:,:,:) … … 218 232 CALL dyn_vor ( kstp ) ! vorticity term including Coriolis 219 233 CALL dyn_ldf ( kstp ) ! lateral mixing 234 IF( ln_stcor ) CALL dyn_stcor ( kstp ) ! Stokes-Coriolis forcing 220 235 CALL dyn_hpg ( kstp ) ! horizontal gradient of Hydrostatic pressure 221 236 CALL dyn_spg ( kstp ) ! surface pressure gradient -
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/step_oce.F90
r5930 r5983 19 19 USE sbcapr ! surface boundary condition: atmospheric pressure 20 20 USE sbctide ! Tide initialisation 21 USE sbcwave ! Wave intialisation 21 22 22 23 USE traqsr ! solar radiation penetration (tra_qsr routine) … … 41 42 USE dynzdf ! vertical diffusion (dyn_zdf routine) 42 43 USE dynspg ! surface pressure gradient (dyn_spg routine) 44 USE dynstcor ! simp. form of Stokes-Coriolis 43 45 44 46 USE dynnxt ! time-stepping (dyn_nxt routine) … … 71 73 USE zdfric ! Richardson vertical mixing (zdf_ric routine) 72 74 USE zdfmxl ! Mixed-layer depth (zdf_mxl routine) 75 USE zdfqiao !Qiao module wave induced mixing (zdf_qiao routine) 73 76 74 77 USE zpshde ! partial step: hor. derivative (zps_hde routine)
Note: See TracChangeset
for help on using the changeset viewer.