Context Navigation

calcul_fluxs_mod.F90 @ 1074

Last change on this file since 1074 was 1067, checked in by Laurent Fairhead, 16 years ago

Modifications lie au premier niveau du modele pour la diffusion turbulent

du vent.

Preparation pour un couplage des courrant oceaniques.

Property svn:eol-style set to native
Property svn:executable set to *
Property svn:keywords set to Author Date Id Revision

File size: 10.9 KB

Rev	Line
[782]	1	!
	2	MODULE calcul_fluxs_mod
	3
	4
	5	CONTAINS
	6	SUBROUTINE calcul_fluxs( knon, nisurf, dtime, &
	7	tsurf, p1lay, cal, beta, coef1lay, ps, &
	8	precip_rain, precip_snow, snow, qsurf, &
	9	radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, &
	10	petAcoef, peqAcoef, petBcoef, peqBcoef, &
	11	tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
	12
	13	USE dimphy, ONLY : klon
	14
	15	! Cette routine calcule les fluxs en h et q a l'interface et eventuellement
	16	! une temperature de surface (au cas ou ok_veget = false)
	17	!
	18	! L. Fairhead 4/2000
	19	!
	20	! input:
	21	! knon nombre de points a traiter
	22	! nisurf surface a traiter
	23	! tsurf temperature de surface
	24	! p1lay pression 1er niveau (milieu de couche)
	25	! cal capacite calorifique du sol
	26	! beta evap reelle
	27	! coef1lay coefficient d'echange
	28	! ps pression au sol
	29	! precip_rain precipitations liquides
	30	! precip_snow precipitations solides
	31	! snow champs hauteur de neige
	32	! runoff runoff en cas de trop plein
	33	! petAcoef coeff. A de la resolution de la CL pour t
	34	! peqAcoef coeff. A de la resolution de la CL pour q
	35	! petBcoef coeff. B de la resolution de la CL pour t
	36	! peqBcoef coeff. B de la resolution de la CL pour q
	37	! radsol rayonnement net aus sol (LW + SW)
	38	! dif_grnd coeff. diffusion vers le sol profond
	39	!
	40	! output:
	41	! tsurf_new temperature au sol
	42	! qsurf humidite de l'air au dessus du sol
	43	! fluxsens flux de chaleur sensible
	44	! fluxlat flux de chaleur latente
	45	! dflux_s derivee du flux de chaleur sensible / Ts
	46	! dflux_l derivee du flux de chaleur latente / Ts
	47	!
	48
[793]	49	INCLUDE "YOETHF.h"
	50	INCLUDE "FCTTRE.h"
	51	INCLUDE "indicesol.h"
	52	INCLUDE "YOMCST.h"
[782]	53
	54	! Parametres d'entree
	55	!****************************************************************************************
	56	INTEGER, INTENT(IN) :: knon, nisurf
	57	REAL , INTENT(IN) :: dtime
	58	REAL, DIMENSION(klon), INTENT(IN) :: petAcoef, peqAcoef
	59	REAL, DIMENSION(klon), INTENT(IN) :: petBcoef, peqBcoef
	60	REAL, DIMENSION(klon), INTENT(IN) :: ps, q1lay
	61	REAL, DIMENSION(klon), INTENT(IN) :: tsurf, p1lay, cal, beta, coef1lay
	62	REAL, DIMENSION(klon), INTENT(IN) :: precip_rain, precip_snow ! pas utiles
	63	REAL, DIMENSION(klon), INTENT(IN) :: radsol, dif_grnd
	64	REAL, DIMENSION(klon), INTENT(IN) :: t1lay, u1lay, v1lay
	65
	66	! Parametres entree-sorties
	67	!****************************************************************************************
	68	REAL, DIMENSION(klon), INTENT(INOUT) :: snow ! snow pas utile
	69
	70	! Parametres sorties
	71	!****************************************************************************************
	72	REAL, DIMENSION(klon), INTENT(OUT) :: qsurf
	73	REAL, DIMENSION(klon), INTENT(OUT) :: tsurf_new, evap, fluxsens, fluxlat
	74	REAL, DIMENSION(klon), INTENT(OUT) :: dflux_s, dflux_l
	75
	76	! Variables locales
	77	!****************************************************************************************
	78	INTEGER :: i
	79	REAL, DIMENSION(klon) :: zx_mh, zx_nh, zx_oh
	80	REAL, DIMENSION(klon) :: zx_mq, zx_nq, zx_oq
	81	REAL, DIMENSION(klon) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
	82	REAL, DIMENSION(klon) :: zx_sl, zx_k1
	83	REAL, DIMENSION(klon) :: d_ts
	84	REAL :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh
	85	REAL :: qsat_new, q1_new
	86	REAL, PARAMETER :: t_grnd = 271.35, t_coup = 273.15
	87	REAL, PARAMETER :: max_eau_sol = 150.0
	88	CHARACTER (len = 20) :: modname = 'calcul_fluxs'
	89	LOGICAL :: fonte_neige
	90	LOGICAL, SAVE :: check = .FALSE.
	91	!$OMP THREADPRIVATE(check)
	92
	93	! End definition
	94	!****************************************************************************************
	95
	96	IF (check) WRITE(,)'Entree ', modname,' surface = ',nisurf
	97
	98	IF (check) THEN
	99	WRITE(,)' radsol (min, max)', &
	100	MINVAL(radsol(1:knon)), MAXVAL(radsol(1:knon))
	101	CALL flush(6)
	102	ENDIF
	103
	104	! Traitement neige et humidite du sol
	105	!****************************************************************************************
	106	!
	107	!!$ WRITE(,)'test calcul_flux, surface ', nisurf
	108	!!PB test
	109	!!$ if (nisurf == is_oce) then
	110	!!$ snow = 0.
	111	!!$ qsol = max_eau_sol
	112	!!$ else
	113	!!$ where (precip_snow > 0.) snow = snow + (precip_snow * dtime)
	114	!!$ where (snow > epsilon(snow)) snow = max(0.0, snow - (evap * dtime))
	115	!!$! snow = max(0.0, snow + (precip_snow - evap) * dtime)
	116	!!$ where (precip_rain > 0.) qsol = qsol + (precip_rain - evap) * dtime
	117	!!$ endif
	118	!!$ IF (nisurf /= is_ter) qsol = max_eau_sol
	119
	120
	121	!
	122	! Initialisation
	123	!****************************************************************************************
	124	evap = 0.
	125	fluxsens=0.
	126	fluxlat=0.
	127	dflux_s = 0.
	128	dflux_l = 0.
	129	!
	130	! zx_qs = qsat en kg/kg
	131	!****************************************************************************************
	132	DO i = 1, knon
	133	zx_pkh(i) = (ps(i)/ps(i))**RKAPPA
	134	IF (thermcep) THEN
	135	zdelta=MAX(0.,SIGN(1.,rtt-tsurf(i)))
	136	zcvm5 = R5LESRLVTT(1.-zdelta) + R5IESRLSTTzdelta
	137	zcvm5 = zcvm5 / RCPD / (1.0+RVTMP2*q1lay(i))
	138	zx_qs= r2es * FOEEW(tsurf(i),zdelta)/ps(i)
	139	zx_qs=MIN(0.5,zx_qs)
	140	zcor=1./(1.-retv*zx_qs)
	141	zx_qs=zx_qs*zcor
	142	zx_dq_s_dh = FOEDE(tsurf(i),zdelta,zcvm5,zx_qs,zcor) &
	143	/RLVTT / zx_pkh(i)
	144	ELSE
	145	IF (tsurf(i).LT.t_coup) THEN
	146	zx_qs = qsats(tsurf(i)) / ps(i)
	147	zx_dq_s_dh = dqsats(tsurf(i),zx_qs)/RLVTT &
	148	/ zx_pkh(i)
	149	ELSE
	150	zx_qs = qsatl(tsurf(i)) / ps(i)
	151	zx_dq_s_dh = dqsatl(tsurf(i),zx_qs)/RLVTT &
	152	/ zx_pkh(i)
	153	ENDIF
	154	ENDIF
	155	zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh
	156	zx_qsat(i) = zx_qs
	157	zx_coef(i) = coef1lay(i) * &
	158	(1.0+SQRT(u1lay(i)2+v1lay(i)2)) * &
	159	p1lay(i)/(RD*t1lay(i))
	160
	161	ENDDO
	162
	163
	164	! === Calcul de la temperature de surface ===
	165	! zx_sl = chaleur latente d'evaporation ou de sublimation
	166	!****************************************************************************************
	167
	168	DO i = 1, knon
	169	zx_sl(i) = RLVTT
	170	IF (tsurf(i) .LT. RTT) zx_sl(i) = RLSTT
	171	zx_k1(i) = zx_coef(i)
	172	ENDDO
	173
	174
	175	DO i = 1, knon
	176	! Q
	177	zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)
	178	zx_mq(i) = beta(i) * zx_k1(i) * &
	179	(peqAcoef(i) - zx_qsat(i) + &
	180	zx_dq_s_dt(i) * tsurf(i)) &
	181	/ zx_oq(i)
	182	zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &
	183	/ zx_oq(i)
	184
	185	! H
	186	zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)
	187	zx_mh(i) = zx_k1(i) * petAcoef(i) / zx_oh(i)
	188	zx_nh(i) = - (zx_k1(i) * RCPD * zx_pkh(i))/ zx_oh(i)
	189
	190	! Tsurface
	191	tsurf_new(i) = (tsurf(i) + cal(i)/(RCPD * zx_pkh(i)) * dtime * &
	192	(radsol(i) + zx_mh(i) + zx_sl(i) * zx_mq(i)) &
	193	+ dif_grnd(i) * t_grnd * dtime)/ &
	194	( 1. - dtime * cal(i)/(RCPD * zx_pkh(i)) * ( &
	195	zx_nh(i) + zx_sl(i) * zx_nq(i)) &
	196	+ dtime * dif_grnd(i))
	197
	198	!
	199	! Y'a-t-il fonte de neige?
	200	!
	201	! fonte_neige = (nisurf /= is_oce) .AND. &
	202	! & (snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &
	203	! & .AND. (tsurf_new(i) >= RTT)
	204	! if (fonte_neige) tsurf_new(i) = RTT
	205	d_ts(i) = tsurf_new(i) - tsurf(i)
	206	! zx_h_ts(i) = tsurf_new(i) * RCPD * zx_pkh(i)
	207	! zx_q_0(i) = zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
	208
	209	!== flux_q est le flux de vapeur d'eau: kg/(m**2 s) positive vers bas
	210	!== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
	211	evap(i) = - zx_mq(i) - zx_nq(i) * tsurf_new(i)
	212	fluxlat(i) = - evap(i) * zx_sl(i)
	213	fluxsens(i) = zx_mh(i) + zx_nh(i) * tsurf_new(i)
	214
	215	! Derives des flux dF/dTs (W m-2 K-1):
	216	dflux_s(i) = zx_nh(i)
	217	dflux_l(i) = (zx_sl(i) * zx_nq(i))
	218
	219	! Nouvelle valeure de l'humidite au dessus du sol
	220	qsat_new=zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
	221	q1_new = peqAcoef(i) - peqBcoef(i)evap(i)dtime
	222	qsurf(i)=q1_new(1.-beta(i)) + beta(i)qsat_new
	223	!
	224	! en cas de fonte de neige
	225	!
	226	! if (fonte_neige) then
	227	! bilan_f = radsol(i) + fluxsens(i) - (zx_sl(i) * evap (i)) - &
	228	! & dif_grnd(i) * (tsurf_new(i) - t_grnd) - &
	229	! & RCPD * (zx_pkh(i))/cal(i)/dtime * (tsurf_new(i) - tsurf(i))
	230	! bilan_f = max(0., bilan_f)
	231	! fq_fonte = bilan_f / zx_sl(i)
	232	! snow(i) = max(0., snow(i) - fq_fonte * dtime)
	233	! qsol(i) = qsol(i) + (fq_fonte * dtime)
	234	! endif
	235	!!$ if (nisurf == is_ter) &
	236	!!$ & run_off(i) = run_off(i) + max(qsol(i) - max_eau_sol, 0.0)
	237	!!$ qsol(i) = min(qsol(i), max_eau_sol)
	238	ENDDO
	239	!
	240	!****************************************************************************************
	241	!
	242	END SUBROUTINE calcul_fluxs
[1067]	243	!
	244	!****************************************************************************************
	245	!
	246	SUBROUTINE calcul_flux_wind(knon, dtime, &
	247	u0, v0, u1, v1, cdrag_m, &
	248	AcoefU, AcoefV, BcoefU, BcoefV, &
	249	p1lay, t1lay, &
	250	flux_u1, flux_v1)
[782]	251
[1067]	252	USE dimphy
	253	INCLUDE "YOMCST.h"
	254
	255	! Input arguments
	256	!****************************************************************************************
	257	INTEGER, INTENT(IN) :: knon
	258	REAL, INTENT(IN) :: dtime
	259	REAL, DIMENSION(klon), INTENT(IN) :: u0, v0 ! u and v at niveau 0
	260	REAL, DIMENSION(klon), INTENT(IN) :: u1, v1 ! u and v at niveau 1
	261	REAL, DIMENSION(klon), INTENT(IN) :: cdrag_m ! cdrag pour momentum
	262	REAL, DIMENSION(klon), INTENT(IN) :: AcoefU, AcoefV, BcoefU, BcoefV
	263	REAL, DIMENSION(klon), INTENT(IN) :: p1lay ! pression 1er niveau (milieu de couche)
	264	REAL, DIMENSION(klon), INTENT(IN) :: t1lay ! temperature
	265	! Output arguments
	266	!****************************************************************************************
	267	REAL, DIMENSION(klon), INTENT(OUT) :: flux_u1
	268	REAL, DIMENSION(klon), INTENT(OUT) :: flux_v1
	269
	270	! Local variables
	271	!****************************************************************************************
	272	INTEGER :: i
	273	REAL :: mod_wind, buf
	274
	275	!****************************************************************************************
	276	! Calculate the surface flux
	277	!
	278	!****************************************************************************************
	279	DO i=1,knon
	280	mod_wind = 1.0 + SQRT((u1(i) - u0(i))2 + (v1(i)-v0(i))2)
	281	buf = cdrag_m(i) * mod_wind * p1lay(i)/(RD*t1lay(i))
	282	flux_u1(i) = (AcoefU(i) - u0(i)) / (1/buf - BcoefU(i)*dtime )
	283	flux_v1(i) = (AcoefV(i) - v0(i)) / (1/buf - BcoefV(i)*dtime )
	284	END DO
	285
	286	END SUBROUTINE calcul_flux_wind
	287	!
	288	!****************************************************************************************
	289	!
[782]	290	END MODULE calcul_fluxs_mod

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

Download in other formats:

Original Format