Context Navigation

surf_ocean_mod.F90 @ 5444

Last change on this file since 5444 was 5285, checked in by abarral, 2 months ago
As discussed internally, remove generic ONLY: ... for new _mod_h modules
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 16.3 KB

Rev	Line
[781]	1	!
[2538]	2	! $Id: surf_ocean_mod.F90 5285 2024-10-28 13:33:29Z fhourdin $
	3	!
[781]	4	MODULE surf_ocean_mod
	5
	6	IMPLICIT NONE
	7
	8	CONTAINS
[3395]	9	!
	10	!******************************************************************************
	11	!
[888]	12	SUBROUTINE surf_ocean(rlon, rlat, swnet, lwnet, alb1, &
[2243]	13	windsp, rmu0, fder, tsurf_in, &
[781]	14	itime, dtime, jour, knon, knindex, &
[4523]	15	p1lay, z1lay, cdragh, cdragm, precip_rain, precip_snow, precip_bs, temp_air, spechum, &
[1067]	16	AcoefH, AcoefQ, BcoefH, BcoefQ, &
	17	AcoefU, AcoefV, BcoefU, BcoefV, &
[2240]	18	ps, u1, v1, gustiness, rugoro, pctsrf, &
[888]	19	snow, qsurf, agesno, &
[2243]	20	z0m, z0h, SFRWL, alb_dir_new, alb_dif_new, evap, fluxsens, fluxlat, &
[1067]	21	tsurf_new, dflux_s, dflux_l, lmt_bils, &
[4370]	22	flux_u1, flux_v1, delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, &
[5022]	23	dt_ds, tkt, tks, taur, sss &
	24	#ifdef ISO
	25	& ,xtprecip_rain, xtprecip_snow,xtspechum,Roce, &
	26	& xtsnow,xtevap,h1 &
	27	#endif
	28	& )
[1067]	29
[3395]	30	use albedo, only: alboc, alboc_cd
[3815]	31	use bulk_flux_m, only: bulk_flux
[3395]	32	USE dimphy, ONLY: klon, zmasq
	33	USE surface_data, ONLY : type_ocean
	34	USE ocean_forced_mod, ONLY : ocean_forced_noice
	35	USE ocean_slab_mod, ONLY : ocean_slab_noice
	36	USE ocean_cpl_mod, ONLY : ocean_cpl_noice
	37	USE indice_sol_mod, ONLY : nbsrf, is_oce
[5022]	38	#ifdef ISO
	39	USE infotrac_phy, ONLY : ntraciso=>ntiso,niso
	40	#ifdef ISOVERIF
	41	USE isotopes_mod, ONLY: iso_eau,ridicule
	42	USE isotopes_verif_mod
	43	#endif
	44	#endif
[5282]	45	USE clesphys_mod_h
[5285]	46	USE yomcst_mod_h
[5274]	47	USE limit_read_mod
[5022]	48	USE config_ocean_skin_m, ONLY: activate_ocean_skin
[3395]	49	!
[5274]	50	! This subroutine will make a call to ocean_XXX_noice according to the ocean mode (force,
	51	! slab or couple). The calculations of albedo and rugosity for the ocean surface are
	52	! done in here because they are identical for the different modes of ocean.
[1785]	53
	54
[781]	55
[5274]	56
[2455]	57	! for cycle_diurne and for iflag_z0_oce==-1 (prescribed z0)
[2178]	58
[3395]	59	! Input variables
	60	!******************************************************************************
[781]	61	INTEGER, INTENT(IN) :: itime, jour, knon
	62	INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
	63	REAL, INTENT(IN) :: dtime
	64	REAL, DIMENSION(klon), INTENT(IN) :: rlon, rlat
[888]	65	REAL, DIMENSION(klon), INTENT(IN) :: swnet ! net shortwave radiation at surface
	66	REAL, DIMENSION(klon), INTENT(IN) :: lwnet ! net longwave radiation at surface
	67	REAL, DIMENSION(klon), INTENT(IN) :: alb1 ! albedo in visible SW interval
[3815]	68	REAL, DIMENSION(klon), INTENT(IN) :: windsp ! wind at 10 m, in m s-1
[781]	69	REAL, DIMENSION(klon), INTENT(IN) :: rmu0
	70	REAL, DIMENSION(klon), INTENT(IN) :: fder
[3906]	71	REAL, DIMENSION(klon), INTENT(IN) :: tsurf_in ! defined only for subscripts 1:knon
[2254]	72	REAL, DIMENSION(klon), INTENT(IN) :: p1lay,z1lay ! pression (Pa) et altitude (m) du premier niveau
[1067]	73	REAL, DIMENSION(klon), INTENT(IN) :: cdragh
	74	REAL, DIMENSION(klon), INTENT(IN) :: cdragm
[4523]	75	REAL, DIMENSION(klon), INTENT(IN) :: precip_rain, precip_snow, precip_bs
[781]	76	REAL, DIMENSION(klon), INTENT(IN) :: temp_air, spechum
[1067]	77	REAL, DIMENSION(klon), INTENT(IN) :: AcoefH, AcoefQ, BcoefH, BcoefQ
	78	REAL, DIMENSION(klon), INTENT(IN) :: AcoefU, AcoefV, BcoefU, BcoefV
[781]	79	REAL, DIMENSION(klon), INTENT(IN) :: ps
[2240]	80	REAL, DIMENSION(klon), INTENT(IN) :: u1, v1, gustiness
[781]	81	REAL, DIMENSION(klon), INTENT(IN) :: rugoro
	82	REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf
[5022]	83	#ifdef ISO
	84	REAL, DIMENSION(ntraciso,klon), INTENT(IN) :: xtprecip_rain, xtprecip_snow
	85	REAL, DIMENSION(ntraciso,klon), INTENT(IN) :: xtspechum
	86	#endif
[781]	87
[3395]	88	! In/Output variables
	89	!******************************************************************************
[888]	90	REAL, DIMENSION(klon), INTENT(INOUT) :: snow
	91	REAL, DIMENSION(klon), INTENT(INOUT) :: qsurf
[781]	92	REAL, DIMENSION(klon), INTENT(INOUT) :: agesno
[3906]	93	REAL, DIMENSION(klon), INTENT(inOUT) :: z0h
[5022]	94	#ifdef ISO
	95	REAL, DIMENSION(niso,klon), INTENT(IN) :: xtsnow
	96	REAL, DIMENSION(niso,klon), INTENT(INOUT):: Roce
	97	#endif
[781]	98
[3815]	99	REAL, intent(inout):: delta_sst(:) ! (knon)
	100	! Ocean-air interface temperature minus bulk SST, in K. Defined
	101	! only if activate_ocean_skin >= 1.
	102
	103	real, intent(inout):: delta_sal(:) ! (knon)
	104	! Ocean-air interface salinity minus bulk salinity, in ppt. Defined
	105	! only if activate_ocean_skin >= 1.
	106
	107	REAL, intent(inout):: ds_ns(:) ! (knon)
	108	! "delta salinity near surface". Salinity variation in the
	109	! near-surface turbulent layer. That is subskin salinity minus
	110	! foundation salinity. In ppt.
	111
	112	REAL, intent(inout):: dt_ns(:) ! (knon)
	113	! "delta temperature near surface". Temperature variation in the
	114	! near-surface turbulent layer. That is subskin temperature
	115	! minus foundation temperature. (Can be negative.) In K.
	116
[4370]	117	REAL, intent(inout):: dter(:) ! (knon)
	118	! Temperature variation in the diffusive microlayer, that is
	119	! ocean-air interface temperature minus subskin temperature. In
	120	! K.
	121
	122	REAL, intent(inout):: dser(:) ! (knon)
	123	! Salinity variation in the diffusive microlayer, that is
	124	! ocean-air interface salinity minus subskin salinity. In ppt.
	125
	126	real, intent(inout):: dt_ds(:) ! (knon)
	127	! (tks / tkt) * dTer, in K
	128
[3395]	129	! Output variables
[3815]	130	!**************************************************************************
[3389]	131	REAL, DIMENSION(klon), INTENT(OUT) :: z0m
[3395]	132	!albedo SB >>>
[3906]	133	! REAL, DIMENSION(klon), INTENT(OUT) :: alb1_new ! new albedo in visible SW interval
	134	! REAL, DIMENSION(klon), INTENT(OUT) :: alb2_new ! new albedo in near IR interval
	135	REAL, DIMENSION(6), INTENT(IN) :: SFRWL
	136	REAL, DIMENSION(klon,nsw), INTENT(OUT) :: alb_dir_new,alb_dif_new
[3395]	137	!albedo SB <<<
[781]	138	REAL, DIMENSION(klon), INTENT(OUT) :: evap, fluxsens, fluxlat
[3906]	139	REAL, DIMENSION(klon), INTENT(OUT) :: tsurf_new ! sea surface temperature, in K
[781]	140	REAL, DIMENSION(klon), INTENT(OUT) :: dflux_s, dflux_l
[996]	141	REAL, DIMENSION(klon), INTENT(OUT) :: lmt_bils
[1067]	142	REAL, DIMENSION(klon), INTENT(OUT) :: flux_u1, flux_v1
[781]	143
[3815]	144	REAL, intent(out):: tkt(:) ! (knon)
[5274]	145	! �paisseur (m) de la couche de diffusion thermique (microlayer)
[3815]	146	! cool skin thickness
	147
	148	REAL, intent(out):: tks(:) ! (knon)
[5274]	149	! �paisseur (m) de la couche de diffusion de masse (microlayer)
[3815]	150
	151	REAL, intent(out):: taur(:) ! (knon)
	152	! momentum flux due to rain, in Pa
	153
	154	real, intent(out):: sss(:) ! (klon)
	155	! Bulk salinity of the surface layer of the ocean, in ppt. (Only
	156	! defined for subscripts 1:knon, but we have to declare it with
	157	! size klon because of the coupling machinery.)
	158
[5022]	159	#ifdef ISO
	160	REAL, DIMENSION(ntraciso,klon), INTENT(out) :: xtevap ! isotopes in surface evaporation flux
	161	REAL, DIMENSION(klon), INTENT(out) :: h1 ! just a diagnostic, not useful for the simulation
	162	#endif
	163
[3395]	164	! Local variables
[3815]	165	!*************************************************************************
[2227]	166	INTEGER :: i, k
[1146]	167	REAL :: tmp
	168	REAL, PARAMETER :: cepdu2=(0.1)**2
[3002]	169	REAL, DIMENSION(klon) :: alb_eau, z0_lim
[888]	170	REAL, DIMENSION(klon) :: radsol
[2254]	171	REAL, DIMENSION(klon) :: cdragq ! Cdrag pour l'evaporation
[4523]	172	REAL, DIMENSION(klon) :: precip_totsnow
[2391]	173	CHARACTER(len=20),PARAMETER :: modname="surf_ocean"
[5022]	174	REAL rhoa(knon) ! density of moist air (kg / m3)
[3815]	175	REAL sens_prec_liq(knon)
[781]	176
[3815]	177	REAL t_int(knon) ! ocean-air interface temperature, in K
[5022]	178	REAL s_int(knon) ! ocean-air interface salinity, in ppt
[888]	179
[3815]	180	!**************************************************************************
[888]	181
[5022]	182	#ifdef ISO
	183	#ifdef ISOVERIF
	184	DO i = 1, knon
	185	IF (iso_eau > 0) THEN
	186	CALL iso_verif_egalite_choix(xtspechum(iso_eau,i), &
	187	& spechum(i),'surf_ocean_mod 117', &
	188	& errmax,errmaxrel)
	189	CALL iso_verif_egalite_choix(xtsnow(iso_eau,i), &
	190	& snow(i),'surf_ocean_mod 127', &
	191	& errmax,errmaxrel)
	192	ENDIF !IF (iso_eau > 0) then
	193	ENDDO !DO i=1,klon
	194	#endif
	195	#endif
[3815]	196
[3395]	197	!******************************************************************************
	198	! Calculate total net radiance at surface
	199	!
	200	!******************************************************************************
[2719]	201	radsol(1:klon) = 0.0 ! initialisation a priori inutile
[888]	202	radsol(1:knon) = swnet(1:knon) + lwnet(1:knon)
	203
[4523]	204
	205	!****************************************************************************************
	206	!Total solid precip
	207
	208	IF (ok_bs) THEN
[4526]	209	precip_totsnow(:)=precip_snow(:)+precip_bs(:)
[4523]	210	ELSE
[4526]	211	precip_totsnow(:)=precip_snow(:)
[4523]	212	ENDIF
	213
	214
[3395]	215	!******************************************************************************
	216	! Cdragq computed from cdrag
	217	! The difference comes only from a factor (f_z0qh_oce) on z0, so that
	218	! it can be computed inside surf_ocean
	219	! More complicated appraches may require the propagation through
	220	! pbl_surface of an independant cdragq variable.
	221	!******************************************************************************
[2254]	222
	223	IF ( f_z0qh_oce .ne. 1.) THEN
[3395]	224	! Si on suit les formulations par exemple de Tessel, on
	225	! a z0h=0.4nu/u, z0q=0.62nu/u, d'ou f_z0qh_oce=0.62/0.4=1.55
[2718]	226	cdragq(1:knon)=cdragh(1:knon)* &
[3395]	227	log(z1lay(1:knon)/z0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*z0h(1:knon)))
[2254]	228	ELSE
[2718]	229	cdragq(1:knon)=cdragh(1:knon)
[2254]	230	ENDIF
	231
[3815]	232	rhoa = PS(:KNON) / (Rd * temp_air(:knon) * (1. + retv * spechum(:knon)))
[3395]	233	!******************************************************************************
	234	! Switch according to type of ocean (couple, slab or forced)
	235	!******************************************************************************
[996]	236	SELECT CASE(type_ocean)
[781]	237	CASE('couple')
[888]	238	CALL ocean_cpl_noice( &
	239	swnet, lwnet, alb1, &
[1067]	240	windsp, fder, &
[781]	241	itime, dtime, knon, knindex, &
[4523]	242	p1lay, cdragh, cdragq, cdragm, precip_rain, precip_totsnow,temp_air,spechum,&
[1067]	243	AcoefH, AcoefQ, BcoefH, BcoefQ, &
	244	AcoefU, AcoefV, BcoefU, BcoefV, &
[3815]	245	ps, u1, v1, gustiness, tsurf_in, &
[888]	246	radsol, snow, agesno, &
[1067]	247	qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
[3815]	248	tsurf_new, dflux_s, dflux_l, sens_prec_liq, sss, delta_sal, rhoa, &
[4370]	249	delta_sst, dTer, dSer, dt_ds)
[781]	250
	251	CASE('slab')
[888]	252	CALL ocean_slab_noice( &
[996]	253	itime, dtime, jour, knon, knindex, &
[4523]	254	p1lay, cdragh, cdragq, cdragm, precip_rain, precip_totsnow, temp_air, spechum,&
[1067]	255	AcoefH, AcoefQ, BcoefH, BcoefQ, &
	256	AcoefU, AcoefV, BcoefU, BcoefV, &
[2240]	257	ps, u1, v1, gustiness, tsurf_in, &
[2209]	258	radsol, snow, &
[1067]	259	qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
[996]	260	tsurf_new, dflux_s, dflux_l, lmt_bils)
[3395]	261
[781]	262	CASE('force')
[888]	263	CALL ocean_forced_noice( &
	264	itime, dtime, jour, knon, knindex, &
[4523]	265	p1lay, cdragh, cdragq, cdragm, precip_rain, precip_totsnow, &
[781]	266	temp_air, spechum, &
[1067]	267	AcoefH, AcoefQ, BcoefH, BcoefQ, &
	268	AcoefU, AcoefV, BcoefU, BcoefV, &
[3815]	269	ps, u1, v1, gustiness, tsurf_in, &
[888]	270	radsol, snow, agesno, &
[1067]	271	qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
[5022]	272	tsurf_new, dflux_s, dflux_l, sens_prec_liq, rhoa &
	273	#ifdef ISO
	274	,xtprecip_rain, xtprecip_snow, xtspechum,Roce,rlat, &
	275	xtsnow,xtevap,h1 &
	276	#endif
	277	)
[781]	278	END SELECT
	279
[3395]	280	!******************************************************************************
	281	! fcodron: compute lmt_bils forced case (same as wfbils_oce / 1.-contfracatm)
	282	!******************************************************************************
[2057]	283	IF (type_ocean.NE.'slab') THEN
[3395]	284	lmt_bils(1:klon)=0.
	285	DO i=1,knon
	286	lmt_bils(knindex(i))=(swnet(i)+lwnet(i)+fluxsens(i)+fluxlat(i)) &
	287	*pctsrf(knindex(i),is_oce)/(1.-zmasq(knindex(i)))
	288	END DO
[2057]	289	END IF
	290
[3395]	291	!******************************************************************************
	292	! Calculate ocean surface albedo
	293	!******************************************************************************
	294	!albedo SB >>>
	295	IF (iflag_albedo==0) THEN
	296	!--old parametrizations of ocean surface albedo
	297	!
	298	IF (iflag_cycle_diurne.GE.1) THEN
	299	!
	300	CALL alboc_cd(rmu0,alb_eau)
	301	!
	302	!--ad-hoc correction for model radiative balance tuning
	303	!--now outside alboc_cd routine
	304	alb_eau(1:klon) = fmagic*alb_eau(1:klon) + pmagic
	305	alb_eau(1:klon)=MIN(MAX(alb_eau(1:klon),0.0),1.0)
	306	!
	307	ELSE
	308	!
	309	CALL alboc(REAL(jour),rlat,alb_eau)
	310	!--ad-hoc correction for model radiative balance tuning
	311	!--now outside alboc routine
	312	alb_eau(1:klon) = fmagic*alb_eau(1:klon) + pmagic
	313	alb_eau(1:klon)=MIN(MAX(alb_eau(1:klon),0.04),0.60)
	314	!
	315	ENDIF
	316	!
	317	DO i =1, knon
	318	DO k=1,nsw
	319	alb_dir_new(i,k) = alb_eau(knindex(i))
	320	ENDDO
	321	ENDDO
	322	!IM 09122015 next line corresponds to the old way of doing in LMDZ5A/IPSLCM5A versions
	323	!albedo for diffuse radiation is taken the same as for direct radiation
	324	alb_dif_new(1:knon,:)=alb_dir_new(1:knon,:)
	325	!IM 09122015 end
	326	!
	327	ELSE IF (iflag_albedo==1) THEN
	328	!--new parametrization of ocean surface albedo by Sunghye Baek
	329	!--albedo for direct and diffuse radiation are different
	330	!
	331	CALL ocean_albedo(knon,rmu0,knindex,windsp,SFRWL,alb_dir_new,alb_dif_new)
	332	!
	333	!--ad-hoc correction for model radiative balance tuning
	334	alb_dir_new(1:knon,:) = fmagic*alb_dir_new(1:knon,:) + pmagic
	335	alb_dif_new(1:knon,:) = fmagic*alb_dif_new(1:knon,:) + pmagic
	336	alb_dir_new(1:knon,:)=MIN(MAX(alb_dir_new(1:knon,:),0.0),1.0)
	337	alb_dif_new(1:knon,:)=MIN(MAX(alb_dif_new(1:knon,:),0.0),1.0)
	338	!
	339	ELSE IF (iflag_albedo==2) THEN
	340	! F. Codron albedo read from limit.nc
	341	CALL limit_read_rug_alb(itime, dtime, jour,&
	342	knon, knindex, z0_lim, alb_eau)
	343	DO i =1, knon
	344	DO k=1,nsw
	345	alb_dir_new(i,k) = alb_eau(i)
	346	ENDDO
	347	ENDDO
	348	alb_dif_new=alb_dir_new
[781]	349	ENDIF
[3395]	350	!albedo SB <<<
[2227]	351
[3395]	352	!******************************************************************************
	353	! Calculate the rugosity
	354	!******************************************************************************
	355	IF (iflag_z0_oce==0) THEN
	356	DO i = 1, knon
	357	tmp = MAX(cepdu2,gustiness(i)+u1(i)2+v1(i)2)
	358	z0m(i) = 0.018cdragm(i) (gustiness(i)+u1(i)2+v1(i)2)/RG &
	359	+ 0.1114e-6 / SQRT(cdragm(i) tmp)
	360	z0m(i) = MAX(1.5e-05,z0m(i))
	361	ENDDO
	362	z0h(1:knon)=z0m(1:knon) ! En attendant mieux
[2243]	363
[3395]	364	ELSE IF (iflag_z0_oce==1) THEN
	365	DO i = 1, knon
	366	tmp = MAX(cepdu2,gustiness(i)+u1(i)2+v1(i)2)
	367	z0m(i) = 0.018cdragm(i) (gustiness(i)+u1(i)2+v1(i)2)/RG &
	368	+ 0.1114e-6 / SQRT(cdragm(i) tmp)
	369	z0m(i) = MAX(1.5e-05,z0m(i))
	370	z0h(i)=0.414e-6 / SQRT(cdragm(i) tmp)
	371	ENDDO
	372	ELSE IF (iflag_z0_oce==-1) THEN
	373	DO i = 1, knon
	374	z0m(i) = z0min
	375	z0h(i) = z0min
	376	ENDDO
	377	ELSE
[2391]	378	CALL abort_physic(modname,'version non prevue',1)
[3395]	379	ENDIF
[3815]	380
	381	if (activate_ocean_skin >= 1) then
	382	if (type_ocean /= 'couple') sss(:knon) = 35.
	383	call bulk_flux(tkt, tks, taur, dter, dser, t_int, s_int, ds_ns, dt_ns, &
	384	u = windsp(:knon), t_ocean_1 = tsurf_new(:knon), s1 = sss(:knon), &
[4523]	385	rain = precip_rain(:knon) + precip_totsnow(:knon), &
[3815]	386	hf = - fluxsens(:knon), hlb = - fluxlat(:knon), &
	387	rnl = - lwnet(:knon), &
	388	tau = sqrt(flux_u1(:knon)2 + flux_v1(:knon)2), rhoa = rhoa, &
	389	xlv = [(rlvtt, i = 1, knon)], rf = - sens_prec_liq, dtime = dtime, &
	390	rns = swnet(:knon))
	391	delta_sst = t_int - tsurf_new(:knon)
	392	delta_sal = s_int - sss(:knon)
[4370]	393
	394	if (activate_ocean_skin == 2) then
	395	tsurf_new(:knon) = t_int
	396	if (type_ocean == 'couple') dt_ds = (tks / tkt) * dter
	397	end if
[3815]	398	end if
	399
[781]	400	END SUBROUTINE surf_ocean
[3815]	401	!****************************************************************************
[3395]	402	!
[781]	403	END MODULE surf_ocean_mod

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

Download in other formats:

Original Format