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.

sbcflx.F90 in branches/UKMO/AMM15_v3_6_STABLE_package_UKEP/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/AMM15_v3_6_STABLE_package_UKEP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90 @ 8757

Last change on this file since 8757 was 8757, checked in by jcastill, 6 years ago
Changes to the branch to make it merge with branch svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/UKMO/r6232_HZG_WAVE-coupling WARNING: This branch will probably not compile on its own - changes in the original branch will have to be copied here
File size: 14.6 KB

Rev	Line
[888]	1	MODULE sbcflx
	2	!!======================================================================
	3	!! * MODULE sbcflx *
	4	!! Ocean forcing: momentum, heat and freshwater flux formulation
	5	!!=====================================================================
[2528]	6	!! History : 1.0 ! 2006-06 (G. Madec) Original code
	7	!! 3.3 ! 2010-10 (S. Masson) add diurnal cycle
[888]	8	!!----------------------------------------------------------------------
	9
	10	!!----------------------------------------------------------------------
	11	!! namflx : flux formulation namlist
[2528]	12	!! sbc_flx : flux formulation as ocean surface boundary condition (forced mode, fluxes read in NetCDF files)
[888]	13	!!----------------------------------------------------------------------
	14	USE oce ! ocean dynamics and tracers
	15	USE dom_oce ! ocean space and time domain
[2528]	16	USE sbc_oce ! surface boundary condition: ocean fields
	17	USE sbcdcy ! surface boundary condition: diurnal cycle on qsr
[888]	18	USE phycst ! physical constants
	19	USE fldread ! read input fields
	20	USE iom ! IOM library
	21	USE in_out_manager ! I/O manager
	22	USE lib_mpp ! distribued memory computing library
	23	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[8416]	24	USE wrk_nemo ! work arrays
[888]	25
	26	IMPLICIT NONE
	27	PRIVATE
	28
	29	PUBLIC sbc_flx ! routine called by step.F90
	30
[8757]	31	INTEGER :: jpfld = 6 ! maximum number of files to read
	32	INTEGER :: jp_utau ! index of wind stress (i-component) file
	33	INTEGER :: jp_vtau ! index of wind stress (j-component) file
	34	INTEGER :: jp_qtot ! index of total (non solar+solar) heat file
	35	INTEGER :: jp_qsr ! index of solar heat file
	36	INTEGER :: jp_emp ! index of evaporation-precipation file
	37	INTEGER :: jp_press ! index of pressure for UKMO shelf fluxes
[888]	38	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf ! structure of input fields (file informations, fields read)
[8059]	39	LOGICAL , PUBLIC :: ln_shelf_flx = .FALSE. ! UKMO SHELF specific flux flag
[8416]	40	LOGICAL , PUBLIC :: ln_rel_wind = .FALSE. ! UKMO SHELF specific flux flag - relative winds
	41	REAL(wp) :: rn_wfac ! multiplication factor for ice/ocean velocity in the calculation of wind stress (clem)
[8059]	42	INTEGER :: jpfld_local ! maximum number of files to read (locally modified depending on ln_shelf_flx)
[888]	43
	44	!! * Substitutions
	45	# include "domzgr_substitute.h90"
[1029]	46	# include "vectopt_loop_substitute.h90"
[888]	47	!!----------------------------------------------------------------------
[2528]	48	!! NEMO/OPA 3.3 , NEMO-consortium (2010)
[1156]	49	!! $Id$
[2715]	50	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[888]	51	!!----------------------------------------------------------------------
	52	CONTAINS
	53
	54	SUBROUTINE sbc_flx( kt )
	55	!!---------------------------------------------------------------------
	56	!! * ROUTINE sbc_flx *
	57	!!
	58	!! ** Purpose : provide at each time step the surface ocean fluxes
	59	!! (momentum, heat, freshwater and runoff)
	60	!!
	61	!! ** Method : - READ each fluxes in NetCDF files:
	62	!! i-component of the stress utau (N/m2)
	63	!! j-component of the stress vtau (N/m2)
	64	!! net downward heat flux qtot (watt/m2)
	65	!! net downward radiative flux qsr (watt/m2)
	66	!! net upward freshwater (evapo - precip) emp (kg/m2/s)
	67	!!
	68	!! CAUTION : - never mask the surface stress fields
[3625]	69	!! - the stress is assumed to be in the (i,j) mesh referential
[888]	70	!!
	71	!! ** Action : update at each time-step
[1695]	72	!! - utau, vtau i- and j-component of the wind stress
	73	!! - taum wind stress module at T-point
	74	!! - wndm 10m wind module at T-point
[3625]	75	!! - qns non solar heat flux including heat flux due to emp
	76	!! - qsr solar heat flux
	77	!! - emp upward mass flux (evap. - precip.)
	78	!! - sfx salt flux; set to zero at nit000 but possibly non-zero
	79	!! if ice is present (computed in limsbc(_2).F90)
[888]	80	!!----------------------------------------------------------------------
	81	INTEGER, INTENT(in) :: kt ! ocean time step
	82	!!
[1695]	83	INTEGER :: ji, jj, jf ! dummy indices
	84	INTEGER :: ierror ! return error code
[4147]	85	INTEGER :: ios ! Local integer output status for namelist read
[1695]	86	REAL(wp) :: zfact ! temporary scalar
	87	REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3
	88	REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient
	89	REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables
[8059]	90	REAL :: cs ! UKMO SHELF: Friction co-efficient at surface
	91	REAL :: totwindspd ! UKMO SHELF: Magnitude of wind speed vector
[8416]	92	REAL(wp), DIMENSION(:,:), POINTER :: zwnd_i, zwnd_j ! wind speed components at T-point
[8059]	93
	94	REAL(wp) :: rhoa = 1.22 ! Air density kg/m3
	95	REAL(wp) :: cdrag = 1.5e-3 ! drag coefficient
[888]	96	!!
	97	CHARACTER(len=100) :: cn_dir ! Root directory for location of flx files
	98	TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist information structures
[8416]	99	TYPE(FLD_N) :: sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp, sn_press ! informations about the fields to be read
	100	LOGICAL :: ln_foam_flx = .FALSE. ! UKMO FOAM specific flux flag
	101	NAMELIST/namsbc_flx/ cn_dir, sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp, &
	102	& ln_foam_flx, sn_press, ln_shelf_flx, ln_rel_wind, &
	103	& rn_wfac
[888]	104	!!---------------------------------------------------------------------
[2528]	105	!
[8757]	106	ln_readtau = .NOT. (ln_wave .AND. ln_tauw )
	107
	108	! prepare the index of the fields that have to be read
	109	jpfld = 0
	110	IF( ln_readtau ) THEN
	111	jp_utau = jpfld+1
	112	jp_vtau = jpfld+2
	113	jpfld = jpfld+2
	114	ELSE
	115	jp_utau = 0 ; jp_vtau = 0
	116	ENDIF
	117	jp_qtot = jpfld+1
	118	jp_qsr = jpfld+2
	119	jp_emp = jpfld+3
	120	jp_press = jpfld+4
	121	jpfld = jpfld+4
	122
[2528]	123	IF( kt == nit000 ) THEN ! First call kt=nit000
[888]	124	! set file information
[4147]	125	REWIND( numnam_ref ) ! Namelist namsbc_flx in reference namelist : Files for fluxes
	126	READ ( numnam_ref, namsbc_flx, IOSTAT = ios, ERR = 901)
	127	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_flx in reference namelist', lwp )
	128
	129	REWIND( numnam_cfg ) ! Namelist namsbc_flx in configuration namelist : Files for fluxes
	130	READ ( numnam_cfg, namsbc_flx, IOSTAT = ios, ERR = 902 )
	131	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_flx in configuration namelist', lwp )
[4624]	132	IF(lwm) WRITE ( numond, namsbc_flx )
[2528]	133	!
	134	! ! check: do we plan to use ln_dm2dc with non-daily forcing?
	135	IF( ln_dm2dc .AND. sn_qsr%nfreqh /= 24 ) &
	136	& CALL ctl_stop( 'sbc_blk_core: ln_dm2dc can be activated only with daily short-wave forcing' )
	137	!
	138	! ! store namelist information in an array
[888]	139	slf_i(jp_utau) = sn_utau ; slf_i(jp_vtau) = sn_vtau
	140	slf_i(jp_qtot) = sn_qtot ; slf_i(jp_qsr ) = sn_qsr
	141	slf_i(jp_emp ) = sn_emp
[2528]	142	!
[8059]	143	ALLOCATE( sf(jpfld), STAT=ierror ) ! set sf structure
	144	IF( ln_shelf_flx ) slf_i(jp_press) = sn_press
	145
	146	! define local jpfld depending on shelf_flx logical
	147	IF( ln_shelf_flx ) THEN
	148	jpfld_local = jpfld
	149	ELSE
	150	jpfld_local = jpfld-1
	151	ENDIF
	152	!
[1133]	153	IF( ierror > 0 ) THEN
	154	CALL ctl_stop( 'sbc_flx: unable to allocate sf structure' ) ; RETURN
[888]	155	ENDIF
[8059]	156	DO ji= 1, jpfld_local
[2528]	157	ALLOCATE( sf(ji)%fnow(jpi,jpj,1) )
	158	IF( slf_i(ji)%ln_tint ) ALLOCATE( sf(ji)%fdta(jpi,jpj,1,2) )
[1200]	159	END DO
[2528]	160	! ! fill sf with slf_i and control print
[1133]	161	CALL fld_fill( sf, slf_i, cn_dir, 'sbc_flx', 'flux formulation for ocean surface boundary condition', 'namsbc_flx' )
[888]	162	!
[3625]	163	sfx(:,:) = 0.0_wp ! salt flux due to freezing/melting (non-zero only if ice is present; set in limsbc(_2).F90)
	164	!
[888]	165	ENDIF
	166
[2528]	167	CALL fld_read( kt, nn_fsbc, sf ) ! input fields provided at the current time-step
	168
	169	IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN ! update ocean fluxes at each SBC frequency
[888]	170
[8416]	171	!!UKMO SHELF wind speed relative to surface currents - put here to allow merging with coupling branch
[8757]	172	IF( ln_shelf_flx .AND. ln_readtau ) THEN
[8416]	173	CALL wrk_alloc( jpi,jpj, zwnd_i, zwnd_j )
	174
	175	IF( ln_rel_wind ) THEN
	176	DO jj = 2, jpjm1
	177	DO ji = fs_2, fs_jpim1 ! vect. opt.
	178	zwnd_i(ji,jj) = ( sf(jp_utau)%fnow(ji,jj,1) - rn_wfac * 0.5 * ( ssu_m(ji-1,jj ) + ssu_m(ji,jj) ))
	179	zwnd_j(ji,jj) = ( sf(jp_vtau)%fnow(ji,jj,1) - rn_wfac * 0.5 * ( ssv_m(ji ,jj-1) + ssv_m(ji,jj) ))
	180	END DO
	181	END DO
	182	CALL lbc_lnk( zwnd_i(:,:) , 'T', -1. )
	183	CALL lbc_lnk( zwnd_j(:,:) , 'T', -1. )
	184	ELSE
	185	zwnd_i(:,:) = sf(jp_utau)%fnow(:,:,1)
	186	zwnd_j(:,:) = sf(jp_vtau)%fnow(:,:,1)
	187	ENDIF
	188	ENDIF
	189
[2528]	190	IF( ln_dm2dc ) THEN ; qsr(:,:) = sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) ! modify now Qsr to include the diurnal cycle
	191	ELSE ; qsr(:,:) = sf(jp_qsr)%fnow(:,:,1)
	192	ENDIF
[888]	193	!CDIR COLLAPSE
[8059]	194	!!UKMO SHELF effect of atmospheric pressure on SSH
	195	! If using ln_apr_dyn, this is done there so don't repeat here.
	196	IF( ln_shelf_flx .AND. .NOT. ln_apr_dyn) THEN
	197	DO jj = 1, jpjm1
	198	DO ji = 1, jpim1
	199	apgu(ji,jj) = (-1.0/rau0)*(sf(jp_press)%fnow(ji+1,jj,1)-sf(jp_press)%fnow(ji,jj,1))/e1u(ji,jj)
	200	apgv(ji,jj) = (-1.0/rau0)*(sf(jp_press)%fnow(ji,jj+1,1)-sf(jp_press)%fnow(ji,jj,1))/e2v(ji,jj)
	201	END DO
	202	END DO
	203	ENDIF ! ln_shelf_flx
[8416]	204
[2528]	205	DO jj = 1, jpj ! set the ocean fluxes from read fields
[1274]	206	DO ji = 1, jpi
[8059]	207	IF( ln_shelf_flx ) THEN
	208	!! UKMO SHELF - need atmospheric pressure to calculate Haney forcing
	209	pressnow(ji,jj) = sf(jp_press)%fnow(ji,jj,1)
	210	!! UKMO SHELF flux files contain wind speed not wind stress
[8757]	211	IF( ln_readtau) THEN
[8416]	212	totwindspd = sqrt(zwnd_i(ji,jj)zwnd_i(ji,jj) + zwnd_j(ji,jj)zwnd_j(ji,jj))
[8059]	213	cs = 0.63 + (0.066 * totwindspd)
[8416]	214	utau(ji,jj) = cs * (rhoa/rau0) * zwnd_i(ji,jj) * totwindspd
	215	vtau(ji,jj) = cs * (rhoa/rau0) * zwnd_j(ji,jj) * totwindspd
[8757]	216	ENDIF
[8059]	217	ELSE
[8757]	218	IF( ln_readtau ) THEN
[8059]	219	utau(ji,jj) = sf(jp_utau)%fnow(ji,jj,1)
	220	vtau(ji,jj) = sf(jp_vtau)%fnow(ji,jj,1)
[8757]	221	ENDIF
[8059]	222	ENDIF
	223	qsr (ji,jj) = sf(jp_qsr )%fnow(ji,jj,1)
	224	IF( ln_foam_flx .OR. ln_shelf_flx ) THEN
	225	!! UKMO FOAM flux files contain non-solar heat flux (qns) rather than total heat flux (qtot)
	226	qns (ji,jj) = sf(jp_qtot)%fnow(ji,jj,1)
	227	!! UKMO FOAM flux files contain the net DOWNWARD freshwater flux P-E rather then E-P
	228	emp (ji,jj) = -1. * sf(jp_emp )%fnow(ji,jj,1)
	229	ELSE
	230	qns (ji,jj) = sf(jp_qtot)%fnow(ji,jj,1) - sf(jp_qsr)%fnow(ji,jj,1)
	231	emp (ji,jj) = sf(jp_emp )%fnow(ji,jj,1)
	232	ENDIF
[1274]	233	END DO
[888]	234	END DO
[3625]	235	! ! add to qns the heat due to e-p
	236	qns(:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp ! mass flux is at SST
	237	!
[8059]	238
	239	!! UKMO FOAM wind fluxes need lbc_lnk calls owing to a bug in interp.exe
	240	IF( ln_foam_flx ) THEN
	241	CALL lbc_lnk( utau(:,:), 'U', -1. )
	242	CALL lbc_lnk( vtau(:,:), 'V', -1. )
	243	ENDIF
	244
[2528]	245	! ! module of wind stress and wind speed at T-point
	246	zcoef = 1. / ( zrhoa * zcdrag )
[1695]	247	!CDIR NOVERRCHK
	248	DO jj = 2, jpjm1
	249	!CDIR NOVERRCHK
	250	DO ji = fs_2, fs_jpim1 ! vect. opt.
	251	ztx = utau(ji-1,jj ) + utau(ji,jj)
	252	zty = vtau(ji ,jj-1) + vtau(ji,jj)
	253	zmod = 0.5 * SQRT( ztx * ztx + zty * zty )
	254	taum(ji,jj) = zmod
	255	wndm(ji,jj) = SQRT( zmod * zcoef )
	256	END DO
	257	END DO
[4990]	258	taum(:,:) = taum(:,:) * tmask(:,:,1) ; wndm(:,:) = wndm(:,:) * tmask(:,:,1)
[1695]	259	CALL lbc_lnk( taum(:,:), 'T', 1. ) ; CALL lbc_lnk( wndm(:,:), 'T', 1. )
	260
[2528]	261	IF( nitend-nit000 <= 100 .AND. lwp ) THEN ! control print (if less than 100 time-step asked)
[888]	262	WRITE(numout,*)
[1274]	263	WRITE(numout,*) ' read daily momentum, heat and freshwater fluxes OK'
[8059]	264	DO jf = 1, jpfld_local
[1274]	265	IF( jf == jp_utau .OR. jf == jp_vtau ) zfact = 1.
	266	IF( jf == jp_qtot .OR. jf == jp_qsr ) zfact = 0.1
	267	IF( jf == jp_emp ) zfact = 86400.
	268	WRITE(numout,*)
	269	WRITE(numout,) ' day: ', ndastp , TRIM(sf(jf)%clvar), ' ', zfact
	270	CALL prihre( sf(jf)%fnow, jpi, jpj, 1, jpi, 20, 1, jpj, 10, zfact, numout )
	271	END DO
	272	CALL FLUSH(numout)
	273	ENDIF
	274	!
[8757]	275	IF( ln_shelf_flx .AND. ln_readtau ) THEN
[8416]	276	CALL wrk_dealloc( jpi,jpj, zwnd_i, zwnd_j )
	277	ENDIF
	278	!
[888]	279	ENDIF
	280	!
	281	END SUBROUTINE sbc_flx
	282
	283	!!======================================================================
	284	END MODULE sbcflx

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

Download in other formats: