New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

sbcwave.F90 in branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90 @ 7905

Last change on this file since 7905 was 7905, checked in by jcastill, 7 years ago

Series of small bug fixes and stetic changes:

-Fix possible bug in the calculation of Stokes-Coriolis
-Move all the wave control variables to namelist namsbc_wave
-Use one namelist variable instead of two to set Stokes drift velocity coupling
-Cap the values of the Craig and Banner constant as calculated from wave input fields to take into account small values of the friction velocity
-Add new Phillips parametrization for Stokes drift vertical velocity, using the inverse depth scale as in Breivik 2015, instead of the peak wave number as calculated from wave input fields
-Better control of the wave fields that are read from file depending on the wave parameters

File size: 27.5 KB

Rev	Line
[2990]	1	MODULE sbcwave
	2	!!======================================================================
	3	!! * MODULE sbcwave *
	4	!! Wave module
	5	!!======================================================================
[7471]	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
[7481]	9	!! - ! 2016-12 (G. Madec, E. Clementi) update Stoke drift computation
	10	!! + add sbc_wave_ini routine
[2990]	11	!!----------------------------------------------------------------------
	12
	13	!!----------------------------------------------------------------------
[7481]	14	!! sbc_stokes : calculate 3D Stokes-drift velocities
[7471]	15	!! sbc_wave : wave data from wave model in netcdf files
[7481]	16	!! sbc_wave_init : initialisation fo surface waves
[2990]	17	!!----------------------------------------------------------------------
[7481]	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
[7471]	22	!
	23	USE iom ! I/O manager library
	24	USE in_out_manager ! I/O manager
	25	USE lib_mpp ! distribued memory computing library
[7481]	26	USE fldread ! read input fields
[7471]	27	USE wrk_nemo !
	28	USE phycst ! physical constants
[2990]	29
	30	IMPLICIT NONE
	31	PRIVATE
	32
[7481]	33	PUBLIC sbc_stokes ! routine called in sbccpl
	34	PUBLIC sbc_wave ! routine called in sbcmod
	35	PUBLIC sbc_wave_init ! routine called in sbcmod
[2990]	36
[7471]	37	! Variables checking if the wave parameters are coupled (if not, they are read from file)
[7481]	38	LOGICAL, PUBLIC :: cpl_hsig = .FALSE.
	39	LOGICAL, PUBLIC :: cpl_phioc = .FALSE.
[7905]	40	LOGICAL, PUBLIC :: cpl_sdrft = .FALSE.
[7481]	41	LOGICAL, PUBLIC :: cpl_wper = .FALSE.
[7878]	42	LOGICAL, PUBLIC :: cpl_wfreq = .FALSE.
[7481]	43	LOGICAL, PUBLIC :: cpl_wnum = .FALSE.
[7797]	44	LOGICAL, PUBLIC :: cpl_tauoc = .FALSE.
[7481]	45	LOGICAL, PUBLIC :: cpl_wdrag = .FALSE.
[2990]	46
[7905]	47	INTEGER :: nn_sdrift ! type of parameterization to calculate vertical Stokes drift
	48	INTEGER, PARAMETER :: jp_breivik = 0 ! Breivik 2015: v_z=v_0[exp(2kz)/(1-8k*z)]
	49	INTEGER, PARAMETER :: jp_phillips = 1 ! Phillips: v_z=v_o[exp(2kz)-betasqrt(-2kpiz)erfc(sqrt(-2kz))]
	50	INTEGER, PARAMETER :: jp_peakph = 2 ! Phillips using the peak wave number read from wave model instead of the inverse depth scale
	51
[7481]	52	INTEGER :: jpfld ! number of files to read for stokes drift
	53	INTEGER :: jp_usd ! index of stokes drift (i-component) (m/s) at T-point
	54	INTEGER :: jp_vsd ! index of stokes drift (j-component) (m/s) at T-point
	55	INTEGER :: jp_hsw ! index of significant wave hight (m) at T-point
	56	INTEGER :: jp_wmp ! index of mean wave period (s) at T-point
[7878]	57	INTEGER :: jp_wfr ! index of wave peak frequency (s^-1) at T-point
[7471]	58
[7481]	59	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient
	60	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift
	61	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_wn ! structure of input fields (file informations, fields read) wave number for Qiao
	62	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean
[7809]	63	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_phioc ! structure of input fields (file informations, fields read) wave to ocean energy
[7481]	64	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: cdn_wave !:
	65	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: hsw, wmp, wnum !:
[7878]	66	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wfreq !:
[7809]	67	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: rn_crban !: Craig and Banner constant for surface breaking waves mixing
[7481]	68	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave !:
	69	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d !:
	70	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: div_sd !: barotropic stokes drift divergence
	71	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ut0sd, vt0sd !: surface Stokes drift velocities at t-point
	72	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd , vsd , wsd !: Stokes drift velocities at u-, v- & w-points, resp.
[7471]	73
	74	# include "vectopt_loop_substitute.h90"
[2990]	75	!!----------------------------------------------------------------------
	76	!! NEMO/OPA 4.0 , NEMO Consortium (2011)
[5215]	77	!! $Id$
[2990]	78	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
	79	!!----------------------------------------------------------------------
	80	CONTAINS
	81
[7471]	82	SUBROUTINE sbc_stokes( )
	83	!!---------------------------------------------------------------------
	84	!! * ROUTINE sbc_stokes *
	85	!!
	86	!! ** Purpose : compute the 3d Stokes Drift according to Breivik et al.,
	87	!! 2014 (DOI: 10.1175/JPO-D-14-0020.1)
	88	!!
	89	!! ** Method : - Calculate Stokes transport speed
	90	!! - Calculate horizontal divergence
	91	!! - Integrate the horizontal divergenze from the bottom
	92	!! ** action
	93	!!---------------------------------------------------------------------
[7481]	94	INTEGER :: jj, ji, jk ! dummy loop argument
	95	INTEGER :: ik ! local integer
	96	REAL(wp) :: ztransp, zfac, ztemp
	97	REAL(wp) :: zdep_u, zdep_v, zkh_u, zkh_v, zda_u, zda_v
	98	REAL(wp), DIMENSION(:,:) , POINTER :: zk_t, zk_u, zk_v, zu0_sd, zv0_sd ! 2D workspace
	99	REAL(wp), DIMENSION(:,:,:), POINTER :: ze3divh ! 3D workspace
	100
[7471]	101	!!---------------------------------------------------------------------
	102	!
	103
[7481]	104	CALL wrk_alloc( jpi,jpj,jpk, ze3divh )
	105	CALL wrk_alloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd )
	106	!
[7878]	107	! select parameterization for the calculation of vertical Stokes drift
[7905]	108	! exp. wave number at t-point
	109	IF( nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips ) THEN ! (Eq. (19) in Breivick et al. (2014) )
[7878]	110	zfac = 2.0_wp * rpi / 16.0_wp
[7905]	111	DO jj = 1, jpj
[7878]	112	DO ji = 1, jpi
	113	! Stokes drift velocity estimated from Hs and Tmean
	114	ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp )
	115	! Stokes surface speed
	116	tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)ut0sd(ji,jj) + vt0sd(ji,jj)vt0sd(ji,jj))
	117	! Wavenumber scale
	118	zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp )
	119	END DO
[7471]	120	END DO
[7878]	121	DO jj = 1, jpjm1 ! exp. wave number & Stokes drift velocity at u- & v-points
	122	DO ji = 1, jpim1
	123	zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
	124	zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
	125	!
	126	zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
	127	zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
	128	END DO
[7481]	129	END DO
[7905]	130	ELSE IF( nn_sdrift==jp_peakph ) THEN ! peak wave number calculated from the peak frequency received by the wave model
	131	DO jj = 1, jpjm1
	132	DO ji = 1, jpim1
	133	zk_u(ji,jj) = 0.5_wp * ( wfreq(ji,jj)wfreq(ji,jj) + wfreq(ji+1,jj)wfreq(ji+1,jj) ) / grav
	134	zk_v(ji,jj) = 0.5_wp * ( wfreq(ji,jj)wfreq(ji,jj) + wfreq(ji,jj+1)wfreq(ji,jj+1) ) / grav
	135	!
	136	zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
	137	zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
	138	END DO
	139	END DO
	140	ENDIF
	141	!
	142	! !== horizontal Stokes Drift 3D velocity ==!
	143	IF( nn_sdrift==jp_breivik ) THEN
[7878]	144	DO jk = 1, jpkm1
	145	DO jj = 2, jpjm1
	146	DO ji = 2, jpim1
	147	zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) )
	148	zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) )
	149	!
	150	zkh_u = zk_u(ji,jj) * zdep_u ! k * depth
	151	zkh_v = zk_v(ji,jj) * zdep_v
	152	! ! Depth attenuation
	153	zda_u = EXP( -2.0_wpzkh_u ) / ( 1.0_wp + 8.0_wpzkh_u )
	154	zda_v = EXP( -2.0_wpzkh_v ) / ( 1.0_wp + 8.0_wpzkh_v )
	155	!
	156	usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
	157	vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
	158	END DO
	159	END DO
	160	END DO
[7905]	161	ELSE IF( nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph ) THEN
[7878]	162	DO jk = 1, jpkm1
	163	DO jj = 2, jpjm1
	164	DO ji = 2, jpim1
	165	zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) )
	166	zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) )
	167	!
	168	zkh_u = zk_u(ji,jj) * zdep_u ! k * depth
	169	zkh_v = zk_v(ji,jj) * zdep_v
[7905]	170	! ! Depth attenuation
	171	zda_u = EXP( -2.0_wpzkh_u ) - SQRT(2.0_wprpizkh_u) ERFC(SQRT(2.0_wp*zkh_u))
	172	zda_v = EXP( -2.0_wpzkh_v ) - SQRT(2.0_wprpizkh_v) ERFC(SQRT(2.0_wp*zkh_v))
[7878]	173	!
	174	usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
	175	vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
	176	END DO
	177	END DO
	178	END DO
[7905]	179	ENDIF
[7878]	180
[7481]	181	CALL lbc_lnk( usd(:,:,:), 'U', vsd(:,:,:), 'V', -1. )
[7471]	182	!
[7481]	183	! !== vertical Stokes Drift 3D velocity ==!
[7471]	184	!
[7481]	185	DO jk = 1, jpkm1 ! Horizontal e3*divergence
[7471]	186	DO jj = 2, jpj
	187	DO ji = fs_2, jpi
[7481]	188	ze3divh(ji,jj,jk) = ( e2u(ji ,jj) * e3u_n(ji ,jj,jk) * usd(ji, jj,jk) &
	189	& - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk) &
	190	& + e1v(ji,jj ) * e3v_n(ji,jj ,jk) * vsd(ji,jj ,jk) &
	191	& - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk) ) * r1_e12t(ji,jj)
[7471]	192	END DO
	193	END DO
	194	END DO
	195	!
	196	IF( .NOT. AGRIF_Root() ) THEN
[7481]	197	IF( nbondi == 1 .OR. nbondi == 2 ) ze3divh(nlci-1, : ,:) = 0._wp ! east
	198	IF( nbondi == -1 .OR. nbondi == 2 ) ze3divh( 2 , : ,:) = 0._wp ! west
	199	IF( nbondj == 1 .OR. nbondj == 2 ) ze3divh( : ,nlcj-1,:) = 0._wp ! north
	200	IF( nbondj == -1 .OR. nbondj == 2 ) ze3divh( : , 2 ,:) = 0._wp ! south
[7471]	201	ENDIF
	202	!
[7481]	203	CALL lbc_lnk( ze3divh, 'T', 1. )
[7471]	204	!
[7481]	205	IF( .NOT. lk_vvl ) THEN ; ik = 1 ! none zero velocity through the sea surface
	206	ELSE ; ik = 2 ! w=0 at the surface (set one for all in sbc_wave_init)
	207	ENDIF
	208	DO jk = jpkm1, ik, -1 ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero)
	209	wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk)
[7471]	210	END DO
	211	#if defined key_bdy
	212	IF( lk_bdy ) THEN
	213	DO jk = 1, jpkm1
[7481]	214	wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:)
[7471]	215	END DO
	216	ENDIF
	217	#endif
[7481]	218	! !== Horizontal divergence of barotropic Stokes transport ==!
	219	div_sd(:,:) = 0._wp
	220	DO jk = 1, jpkm1 !
	221	div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk)
	222	END DO
[7471]	223	!
[7600]	224	CALL iom_put( "ustokes", usd )
	225	CALL iom_put( "vstokes", vsd )
	226	CALL iom_put( "wstokes", wsd )
	227	!
[7481]	228	CALL wrk_dealloc( jpi,jpj,jpk, ze3divh )
	229	CALL wrk_dealloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd )
	230	!
[7471]	231	END SUBROUTINE sbc_stokes
	232
[7481]	233
	234	SUBROUTINE sbc_wave( kt )
[7471]	235	!!---------------------------------------------------------------------
[7481]	236	!! * ROUTINE sbc_wave *
[7471]	237	!!
[7481]	238	!! ** Purpose : read wave parameters from wave model in netcdf files.
[7471]	239	!!
[7481]	240	!! ** Method : - Read namelist namsbc_wave
	241	!! - Read Cd_n10 fields in netcdf files
	242	!! - Read stokes drift 2d in netcdf files
	243	!! - Read wave number in netcdf files
	244	!! - Compute 3d stokes drift using Breivik et al.,2014
	245	!! formulation
[7471]	246	!! ** action
	247	!!---------------------------------------------------------------------
[7481]	248	INTEGER, INTENT(in ) :: kt ! ocean time step
	249	!!---------------------------------------------------------------------
	250	!
	251	IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN !== Neutral drag coefficient ==!
	252	CALL fld_read( kt, nn_fsbc, sf_cd ) ! read from external forcing
	253	cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1)
[7792]	254	! check that the drag coefficient contains proper information even if
	255	! the masks do not match - the momentum stress is not masked!
	256	WHERE( cdn_wave < 0.0 ) cdn_wave = 1.5e-3
	257	WHERE( cdn_wave > 1.0 ) cdn_wave = 1.5e-3
[7481]	258	ENDIF
[7471]	259
[7797]	260	IF( ln_tauoc .AND. .NOT. cpl_tauoc ) THEN !== Wave induced stress ==!
[7481]	261	CALL fld_read( kt, nn_fsbc, sf_tauoc ) ! read wave norm stress from external forcing
	262	tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1)
[7792]	263	WHERE( tauoc_wave < -100.0 ) tauoc_wave = 1.0
	264	WHERE( tauoc_wave > 100.0 ) tauoc_wave = 1.0
[7481]	265	ENDIF
	266
[7809]	267	IF( ln_phioc .AND. .NOT. cpl_phioc ) THEN !== Wave to ocean energy ==!
	268	CALL fld_read( kt, nn_fsbc, sf_phioc ) ! read wave to ocean energy from external forcing
[7853]	269	rn_crban(:,:) = 29.0 * sf_phioc(1)%fnow(:,:,1) ! ! Alfa is phioc*sqrt(rau0/zrhoa) : rau0=water density, zhroa= air density
[7905]	270	WHERE( rn_crban > 1.e8 ) rn_crban = 0.0 !remove first mask mistmatch points, then cap values in case of low friction velocity
	271	WHERE( rn_crban < 0.0 ) rn_crban = 0.0
	272	WHERE( rn_crban > 1000.0 ) rn_crban = 1000.0
[7809]	273	ENDIF
	274
[7905]	275	IF( ln_sdw .OR. ln_rough ) THEN !== Computation of the 3d Stokes Drift ==!
[7481]	276	!
	277	IF( jpfld > 0 ) THEN ! Read from file only if the field is not coupled
	278	CALL fld_read( kt, nn_fsbc, sf_sd ) ! read wave parameters from external forcing
[7853]	279	IF( jp_hsw > 0 ) THEN
	280	hsw (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) ! significant wave height
	281	WHERE( hsw > 100.0 ) hsw = 0.0
	282	WHERE( hsw < 0.0 ) hsw = 0.0
	283	ENDIF
	284	IF( jp_wmp > 0 ) THEN
	285	wmp (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) ! wave mean period
	286	WHERE( wmp > 100.0 ) wmp = 0.0
	287	WHERE( wmp < 0.0 ) wmp = 0.0
	288	ENDIF
[7878]	289	IF( jp_wfr > 0 ) THEN
	290	wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) ! Peak wave frequency
[7905]	291	WHERE( wfreq < 0.0 ) wfreq = 0.0
	292	WHERE( wfreq > 100.0 ) wfreq = 0.0
[7878]	293	ENDIF
[7853]	294	IF( jp_usd > 0 ) THEN
	295	ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) ! 2D zonal Stokes Drift at T point
	296	WHERE( ut0sd < -100.0 ) ut0sd = 1.0
	297	WHERE( ut0sd > 100.0 ) ut0sd = 1.0
	298	ENDIF
	299	IF( jp_vsd > 0 ) THEN
	300	vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) ! 2D meridional Stokes Drift at T point
	301	WHERE( vt0sd < -100.0 ) vt0sd = 1.0
	302	WHERE( vt0sd > 100.0 ) vt0sd = 1.0
	303	ENDIF
[7481]	304	ENDIF
[7905]	305	ENDIF
	306	!
	307	IF( ln_sdw ) THEN
[7481]	308	! Read also wave number if needed, so that it is available in coupling routines
	309	IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN
	310	CALL fld_read( kt, nn_fsbc, sf_wn ) ! read wave parameters from external forcing
	311	wnum(:,:) = sf_wn(1)%fnow(:,:,1)
	312	ENDIF
	313
	314	! !== Computation of the 3d Stokes Drift ==!
	315	!
[7905]	316	IF( ((nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) .AND. &
	317	jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0) .OR. &
	318	(nn_sdrift==jp_peakph .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0) ) &
[7878]	319	CALL sbc_stokes() ! Calculate only if required fields are read
	320	! ! In coupled wave model-NEMO case the call is done after coupling
[7481]	321	!
	322	ENDIF
[7471]	323	!
[7481]	324	END SUBROUTINE sbc_wave
[7471]	325
[7481]	326
	327	SUBROUTINE sbc_wave_init
[2990]	328	!!---------------------------------------------------------------------
[7481]	329	!! * ROUTINE sbc_wave_init *
[2990]	330	!!
[7471]	331	!! ** Purpose : read wave parameters from wave model in netcdf files.
[2990]	332	!!
	333	!! ** Method : - Read namelist namsbc_wave
	334	!! - Read Cd_n10 fields in netcdf files
[3680]	335	!! - Read stokes drift 2d in netcdf files
[7471]	336	!! - Read wave number in netcdf files
	337	!! - Compute 3d stokes drift using Breivik et al.,2014
	338	!! formulation
	339	!! ** action
[2990]	340	!!---------------------------------------------------------------------
[7481]	341	INTEGER :: ierror, ios ! local integer
	342	INTEGER :: ifpr
	343	!!
[3680]	344	CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files
[7471]	345	TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: slf_i ! array of namelist informations on the fields to read
[7878]	346	TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_phioc, &
	347	& sn_hsw, sn_wmp, sn_wfr, sn_wnum , &
	348	& sn_tauoc ! informations about the fields to be read
[7481]	349	!
[7905]	350	NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, sn_wnum, sn_tauoc, sn_phioc, &
	351	ln_cdgw, ln_sdw, ln_stcor, ln_phioc, ln_tauoc, ln_zdfqiao, ln_rough, &
	352	nn_drag, nn_sdrift, nn_wmix
[2990]	353	!!---------------------------------------------------------------------
	354	!
[7481]	355	REWIND( numnam_ref ) ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model
	356	READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901)
	357	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp )
	358
	359	REWIND( numnam_cfg ) ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model
	360	READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 )
	361	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp )
	362	IF(lwm) WRITE ( numond, namsbc_wave )
	363	!
[7905]	364	IF(lwp) THEN !* Parameter print
	365	WRITE(numout,*)
	366	WRITE(numout,*) 'sbc_wave_init: wave physics'
	367	WRITE(numout,*) '~~~~~~~~'
	368	WRITE(numout,*) ' Namelist namsbc_wave : set wave physics parameters'
	369	WRITE(numout,*) ' Stokes drift corr. to vert. velocity ln_sdw = ', ln_sdw
	370	WRITE(numout,*) ' vertical parametrization nn_sdrift = ', nn_sdrift
	371	WRITE(numout,*) ' Stokes coriolis term ln_stcor = ', ln_stcor
	372	WRITE(numout,*) ' wave modified ocean stress ln_tauoc = ', ln_tauoc
	373	WRITE(numout,*) ' wave to ocean energy ln_phioc = ', ln_phioc
	374	WRITE(numout,*) ' vertical mixing parametrization nn_wmix = ', nn_wmix
	375	WRITE(numout,*) ' neutral drag coefficient ln_cdgw = ', ln_cdgw
	376	WRITE(numout,*) ' momentum formulation nn_drag = ', nn_drag
	377	WRITE(numout,*) ' wave roughness length modification ln_rough = ', ln_rough
	378	WRITE(numout,*) ' Qiao vertical mixing formulation ln_zdfqiao = ', ln_zdfqiao
	379	ENDIF
	380
	381	IF ( ln_wave ) THEN
	382	! Activated wave physics but no wave physics components activated
	383	IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao) ) THEN
	384	CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_stcor=F, ln_phioc=F ', &
	385	'ln_rough=F, ln_zdfqiao=F' )
	386	ELSEIF (ln_stcor .AND. .NOT. ln_sdw) THEN
	387	CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T')
	388	ENDIF
	389	IF ( ln_cdgw .AND. .NOT.(nn_drag==jp_ukmo .OR. nn_drag==jp_std .OR. nn_drag==jp_const .OR. nn_drag==jp_mcore) ) &
	390	CALL ctl_stop( 'The chosen nn_drag for momentum calculation must be 0, 1, 2, or 3')
	391	IF ( ln_cdgw .AND. ln_blk_core .AND. nn_drag==0 ) &
	392	CALL ctl_stop( 'The chosen nn_drag for momentum calculation in core forcing must be 1, 2, or 3')
	393	IF ( ln_cdgw .AND. ln_flx .AND. nn_drag==3 ) &
	394	CALL ctl_stop( 'The chosen nn_drag for momentum calculation in direct forcing must be 0, 1, or 2')
	395	IF( ln_phioc .AND. .NOT.(nn_wmix==jp_craigbanner .OR. nn_wmix==jp_janssen) ) &
	396	CALL ctl_stop( 'The chosen nn_wmix for wave vertical mixing must be 0, or 1' )
	397	IF( ln_sdw .AND. .NOT.(nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph) ) &
	398	CALL ctl_stop( 'The chosen nn_sdrift for Stokes drift vertical velocity must be 0, 1, or 2' )
	399	IF( ln_zdfqiao .AND. .NOT.ln_sdw ) &
	400	CALL ctl_stop( 'Qiao vertical mixing can not be used without Stokes drift (ln_sdw)' )
	401	ELSE
	402	IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao ) &
	403	& CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ', &
	404	& 'with drag coefficient (ln_cdgw =T) ' , &
	405	& 'or Stokes Drift (ln_sdw=T) ' , &
	406	& 'or Stokes-Coriolis term (ln_stcor=T)', &
	407	& 'or ocean stress modification due to waves (ln_tauoc=T) ', &
	408	& 'or wave to ocean energy modification (ln_phioc=T) ', &
	409	& 'or wave surface roughness (ln_rough=T) ', &
	410	& 'or Qiao vertical mixing formulation (ln_zdfqiao=T) ' )
	411	ENDIF
	412	!
[7481]	413	IF( ln_cdgw ) THEN
	414	IF( .NOT. cpl_wdrag ) THEN
	415	ALLOCATE( sf_cd(1), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg
[7809]	416	IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_cd structure' )
[7481]	417	!
	418	ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) )
	419	IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) )
[7905]	420	CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
[7471]	421	ENDIF
[7481]	422	ALLOCATE( cdn_wave(jpi,jpj) )
	423	ENDIF
[7471]	424
[7481]	425	IF( ln_tauoc ) THEN
[7797]	426	IF( .NOT. cpl_tauoc ) THEN
[7481]	427	ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc
[7809]	428	IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauoc structure' )
[7481]	429	!
	430	ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) )
	431	IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) )
[7905]	432	CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
[7471]	433	ENDIF
[7481]	434	ALLOCATE( tauoc_wave(jpi,jpj) )
	435	ENDIF
[7471]	436
[7809]	437	IF( ln_phioc ) THEN
	438	IF( .NOT. cpl_phioc ) THEN
	439	ALLOCATE( sf_phioc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_phioc
	440	IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_phioc structure' )
	441	!
	442	ALLOCATE( sf_phioc(1)%fnow(jpi,jpj,1) )
	443	IF( sn_phioc%ln_tint ) ALLOCATE( sf_phioc(1)%fdta(jpi,jpj,1,2) )
[7905]	444	CALL fld_fill( sf_phioc, (/ sn_phioc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
[7809]	445	ENDIF
	446	ALLOCATE( rn_crban(jpi,jpj) )
	447	ENDIF
	448
[7905]	449	! Find out how many fields have to be read from file if not coupled
	450	jpfld=0
	451	jp_usd=0 ; jp_vsd=0 ; jp_hsw=0 ; jp_wmp=0 ; jp_wfr=0
	452	IF( ln_sdw ) THEN
	453	IF( .NOT. cpl_sdrft ) THEN
[7481]	454	jpfld = jpfld + 1
	455	jp_usd = jpfld
	456	jpfld = jpfld + 1
	457	jp_vsd = jpfld
	458	ENDIF
[7905]	459	IF( .NOT. cpl_hsig .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN
[7481]	460	jpfld = jpfld + 1
	461	jp_hsw = jpfld
	462	ENDIF
[7905]	463	IF( .NOT. cpl_wper .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN
[7481]	464	jpfld = jpfld + 1
	465	jp_wmp = jpfld
	466	ENDIF
[7905]	467	IF( .NOT. cpl_wfreq .AND. nn_sdrift==jp_peakph ) THEN
[7878]	468	jpfld = jpfld + 1
	469	jp_wfr = jpfld
	470	ENDIF
[7905]	471	ENDIF
[2990]	472
[7905]	473	IF( ln_rough .AND. .NOT. cpl_hsig .AND. jp_hsw==0 ) THEN
	474	jpfld = jpfld + 1
	475	jp_hsw = jpfld
	476	ENDIF
	477
	478	! Read from file only the non-coupled fields
	479	IF( jpfld > 0 ) THEN
	480	ALLOCATE( slf_i(jpfld) )
	481	IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd
	482	IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd
	483	IF( jp_hsw > 0 ) slf_i(jp_hsw) = sn_hsw
	484	IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp
	485	IF( jp_wfr > 0 ) slf_i(jp_wfr) = sn_wfr
	486	ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift
	487	IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_sd structure' )
	488	!
	489	DO ifpr= 1, jpfld
	490	ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) )
	491	IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) )
	492	END DO
	493	!
	494	CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
	495	ENDIF
	496
	497	IF( ln_sdw ) THEN
[7481]	498	ALLOCATE( usd (jpi,jpj,jpk), vsd (jpi,jpj,jpk), wsd(jpi,jpj,jpk) )
[7905]	499	ALLOCATE( wmp (jpi,jpj) )
[7878]	500	ALLOCATE( wfreq (jpi,jpj) )
[7481]	501	ALLOCATE( ut0sd(jpi,jpj) , vt0sd(jpi,jpj) )
	502	ALLOCATE( div_sd(jpi,jpj) )
	503	ALLOCATE( tsd2d (jpi,jpj) )
	504	usd(:,:,:) = 0._wp
	505	vsd(:,:,:) = 0._wp
	506	wsd(:,:,:) = 0._wp
	507	! Wave number needed only if ln_zdfqiao=T
[7606]	508	IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN
[7481]	509	ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum
[7809]	510	IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wn structure' )
[7481]	511	ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) )
	512	IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) )
[7905]	513	CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'read wave input', 'namsbc_wave' )
[7471]	514	ENDIF
[7481]	515	ALLOCATE( wnum(jpi,jpj) )
[3680]	516	ENDIF
[7905]	517
	518	IF( ln_sdw .OR. ln_rough ) THEN
	519	ALLOCATE( hsw (jpi,jpj) )
	520	ENDIF
[7471]	521	!
[7481]	522	END SUBROUTINE sbc_wave_init
	523
[2990]	524	!!======================================================================
	525	END MODULE sbcwave

Note: See TracBrowser for help on using the repository browser.

Download in other formats: