Context Navigation

thermcell_main.F90 @ 2276

Last change on this file since 2276 was 2232, checked in by aboissinot, 5 years ago
The thermal plume model is able to manage several plumes in the same column and work without the convective adjustment.
File size: 18.2 KB

Rev	Line
[2060]	1	!
	2	!
	3	!
[2127]	4	SUBROUTINE thermcell_main(ngrid,nlay,nq,ptimestep,firstcall, &
	5	pplay,pplev,pphi,zpopsk, &
	6	pu,pv,pt,pq, &
	7	pduadj,pdvadj,pdtadj,pdqadj, &
[2232]	8	fm_tot,entr_tot,detr_tot,zw2_tot,fraca)
[2060]	9
[2101]	10
[2127]	11	!===============================================================================
[2060]	12	! Auteurs: Frederic Hourdin, Catherine Rio, Anne Mathieu
	13	! Version du 09.02.07
	14	! Calcul du transport vertical dans la couche limite en presence
	15	! de "thermiques" explicitement representes avec processus nuageux
	16	!
	17	! Reecriture a partir d'un listing papier a Habas, le 14/02/00
	18	!
	19	! le thermique est suppose homogene et dissipe par melange avec
	20	! son environnement. la longueur l_mix controle l'efficacite du
	21	! melange
	22	!
	23	! Le calcul du transport des differentes especes se fait en prenant
	24	! en compte:
	25	! 1. un flux de masse montant
	26	! 2. un flux de masse descendant
	27	! 3. un entrainement
	28	! 4. un detrainement
	29	!
	30	! Modif 2013/01/04 (FH hourdin@lmd.jussieu.fr)
	31	! Introduction of an implicit computation of vertical advection in
	32	! the environment of thermal plumes in thermcell_dq
[2127]	33	! impl = 0 : explicit ; impl = 1 : implicit ; impl =-1 : old version
[2060]	34	! controled by iflag_thermals =
	35	! 15, 16 run with impl=-1 : numerical convergence with NPv3
	36	! 17, 18 run with impl=1 : more stable
	37	! 15 and 17 correspond to the activation of the stratocumulus "bidouille"
[2127]	38	!
	39	! Major changes 2018-19 (AB alexandre.boissinot@lmd.jussieu.fr)
	40	! New detr and entre formulae (no longer alimentation)
	41	! lmin can be greater than 1
[2177]	42	! Mix every tracer
[2232]	43	! Can stack verticaly multiple plumes (it makes thermcell_dv2 unusable for the moment)
[2127]	44	!
	45	!===============================================================================
[2060]	46
	47	USE thermcell_mod
[2143]	48	USE print_control_mod, ONLY: prt_level
[2060]	49
	50	IMPLICIT NONE
	51
	52
[2127]	53	!===============================================================================
[2060]	54	! Declaration
[2127]	55	!===============================================================================
[2060]	56
[2127]	57	! Inputs:
	58	! -------
[2060]	59
[2177]	60	INTEGER, INTENT(in) :: ngrid
	61	INTEGER, INTENT(in) :: nlay
	62	INTEGER, INTENT(in) :: nq
[2060]	63
[2177]	64	REAL, INTENT(in) :: ptimestep
	65	REAL, INTENT(in) :: pplay(ngrid,nlay) ! Layer pressure
	66	REAL, INTENT(in) :: pplev(ngrid,nlay+1) ! Level pressure
	67	REAL, INTENT(in) :: pphi(ngrid,nlay) ! Geopotential
[2232]	68	REAL, INTENT(in) :: zpopsk(ngrid,nlay) ! Exner function
[2060]	69
[2177]	70	REAL, INTENT(in) :: pu(ngrid,nlay) ! Zonal wind
	71	REAL, INTENT(in) :: pv(ngrid,nlay) ! Meridional wind
	72	REAL, INTENT(in) :: pt(ngrid,nlay) ! Temperature
	73	REAL, INTENT(in) :: pq(ngrid,nlay,nq) ! Tracers mass mixing ratio
[2060]	74
[2177]	75	LOGICAL, INTENT(in) :: firstcall
[2060]	76
[2127]	77	! Outputs:
	78	! --------
[2060]	79
[2177]	80	REAL, INTENT(out) :: pduadj(ngrid,nlay) ! u convective variations
	81	REAL, INTENT(out) :: pdvadj(ngrid,nlay) ! v convective variations
	82	REAL, INTENT(out) :: pdtadj(ngrid,nlay) ! t convective variations
	83	REAL, INTENT(out) :: pdqadj(ngrid,nlay,nq) ! q convective variations
[2060]	84
[2232]	85	REAL, INTENT(inout) :: fm_tot(ngrid,nlay+1) ! Total mass flux
	86	REAL, INTENT(inout) :: entr_tot(ngrid,nlay) ! Total entrainment
	87	REAL, INTENT(inout) :: detr_tot(ngrid,nlay) ! Total detrainment
[2177]	88
[2232]	89	REAL, INTENT(out) :: fraca(ngrid,nlay+1) ! Updraft fraction
	90	REAL, INTENT(out) :: zw2_tot(ngrid,nlay+1) ! Total plume vertical speed
	91
[2127]	92	! Local:
	93	! ------
[2060]	94
[2127]	95	INTEGER ig, k, l, iq
[2232]	96	INTEGER lmax(ngrid) ! Highest layer reached by the plume
	97	INTEGER lmix(ngrid) ! Layer in which plume vertical speed is maximal
	98	INTEGER lmin(ngrid) ! First unstable layer
[2060]	99
[2177]	100	REAL zmix(ngrid) ! Altitude of maximal vertical speed
	101	REAL zmax(ngrid) ! Maximal altitudes where plumes are active
	102	REAL zmin(ngrid) ! Minimal altitudes where plumes are active
[2143]	103
[2177]	104	REAL zlay(ngrid,nlay) ! Layers altitudes
	105	REAL zlev(ngrid,nlay+1) ! Levels altitudes
	106	REAL rho(ngrid,nlay) ! Layers densities
	107	REAL rhobarz(ngrid,nlay) ! Levels densities
	108	REAL masse(ngrid,nlay) ! Layers masses
[2060]	109
[2177]	110	REAL zu(ngrid,nlay) ! u environment
	111	REAL zv(ngrid,nlay) ! v environment
	112	REAL zt(ngrid,nlay) ! TR environment
	113	REAL zqt(ngrid,nlay) ! qt environment
	114	REAL zql(ngrid,nlay) ! ql environment
	115	REAL zhl(ngrid,nlay) ! TP environment
	116	REAL ztv(ngrid,nlay) ! TRPV environment
	117	REAL zqs(ngrid,nlay) ! qsat environment
[2060]	118
[2177]	119	REAL zua(ngrid,nlay) ! u plume
	120	REAL zva(ngrid,nlay) ! v plume
	121	REAL zqla(ngrid,nlay) ! qv plume
	122	REAL zqta(ngrid,nlay) ! qt plume
	123	REAL zhla(ngrid,nlay) ! TP plume
	124	REAL ztva(ngrid,nlay) ! TRPV plume
	125	REAL zqsa(ngrid,nlay) ! qsat plume
[2060]	126
[2177]	127	REAL zqa(ngrid,nlay,nq) ! q plume (ql=0, qv=qt)
[2060]	128
[2177]	129	REAL f_star(ngrid,nlay+1) ! Normalized mass flux
[2232]	130	REAL entr_star(ngrid,nlay) ! Normalized entrainment
	131	REAL detr_star(ngrid,nlay) ! Normalized detrainment
[2060]	132
[2232]	133	REAL f(ngrid) ! Mass flux norm
[2177]	134	REAL fm(ngrid,nlay+1) ! Mass flux
[2232]	135	REAL entr(ngrid,nlay) ! Entrainment
	136	REAL detr(ngrid,nlay) ! Detrainment
[2060]	137
[2232]	138	REAL zw2(ngrid,nlay+1) ! Plume vertical speed
[2177]	139	REAL wmax(ngrid) ! Maximal vertical speed
	140	REAL zdthladj(ngrid,nlay) ! Potential temperature variations
	141	REAL dummy(ngrid,nlay) ! Dummy argument for thermcell_dq()
[2060]	142
[2232]	143	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	144	INTEGER lbot(ngrid)
	145	LOGICAL re_tpm
	146	INTEGER while_loop_counter
	147	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	148
[2127]	149	!===============================================================================
[2060]	150	! Initialization
[2127]	151	!===============================================================================
[2060]	152
[2232]	153	fm_tot(:,:) = 0.
	154	entr_tot(:,:) = 0.
	155	detr_tot(:,:) = 0.
	156	zw2_tot(:,:) = 0.
[2060]	157
[2127]	158	pduadj(:,:) = 0.0
	159	pdvadj(:,:) = 0.0
	160	pdtadj(:,:) = 0.0
	161	pdqadj(:,:,:) = 0.0
	162
[2177]	163	zdthladj(:,:) = 0.0
[2060]	164
[2232]	165	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	166	re_tpm = .true.
	167	lbot(:) = linf
	168	while_loop_counter = 0.
	169	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	170
[2127]	171	!===============================================================================
	172	! Environment settings
	173	!===============================================================================
[2060]	174
[2127]	175	!-------------------------------------------------------------------------------
[2060]	176	! Calcul de T,q,ql a partir de Tl et qt dans l environnement
[2127]	177	!-------------------------------------------------------------------------------
[2060]	178
[2127]	179	CALL thermcell_env(ngrid,nlay,nq,pq,pt,pu,pv,pplay,pplev, &
	180	& zqt,zql,zt,ztv,zhl,zu,zv,zpopsk,zqs)
[2060]	181
[2127]	182	!-------------------------------------------------------------------------------
	183	! Levels and layers altitudes
	184	!-------------------------------------------------------------------------------
[2060]	185
	186	DO l=2,nlay
	187	zlev(:,l) = 0.5 * (pphi(:,l) + pphi(:,l-1)) / RG
	188	ENDDO
	189
	190	zlev(:,1) = 0.
[2101]	191	zlev(:,nlay+1) = (2. * pphi(:,nlay) - pphi(:,nlay-1)) / RG
[2060]	192
	193	DO l=1,nlay
[2143]	194	zlay(:,l) = pphi(:,l) / RG
[2060]	195	ENDDO
	196
[2127]	197	!-------------------------------------------------------------------------------
	198	! Levels and layers densities
	199	!-------------------------------------------------------------------------------
[2060]	200
[2127]	201	rho(:,:) = pplay(:,:) / (zpopsk(:,:) * RD * ztv(:,:))
[2060]	202
[2177]	203	rhobarz(:,1) = rho(:,1)
[2060]	204	IF (prt_level.ge.10) THEN
[2143]	205	print *, 'WARNING: density in the first layer is equal to density at the first level!'
[2060]	206	ENDIF
	207
	208	DO l=2,nlay
	209	rhobarz(:,l) = 0.5 * (rho(:,l) + rho(:,l-1))
	210	ENDDO
	211
[2127]	212	!-------------------------------------------------------------------------------
	213	! Layers masses
	214	!-------------------------------------------------------------------------------
[2060]	215
	216	DO l=1,nlay
	217	masse(:,l) = (pplev(:,l) - pplev(:,l+1)) / RG
	218	ENDDO
	219
[2127]	220	!===============================================================================
	221	! Explicative schemes
	222	!===============================================================================
	223
	224	!-------------------------------------------------------------------------------
	225	! Thermal plume variables
	226	!-------------------------------------------------------------------------------
[2060]	227
[2127]	228	! top of the model
	229	! ===========================
	230	!
	231	! ---------------------------
[2177]	232	! _
	233	! ----- F_lmax+1=0 ------zmax \
	234	! lmax \|
	235	! ------F_lmax>0------------- \|
	236	! \|
	237	! --------------------------- \|
	238	! \|
	239	! --------------------------- \|
	240	! \|
	241	! ------------------wmax,zmix \|
	242	! lmix \|
	243	! --------------------------- \|
	244	! \|
	245	! --------------------------- \|
	246	! \| E, D
	247	! --------------------------- \|
	248	! \|
[2232]	249	! --------------------------- rhobarz, f_star, fm, zw2, fraca
[2177]	250	! zt, zu, zv, zo, rho \|
	251	! --------------------------- \|
	252	! \|
	253	! --------------------------- \|
	254	! \|
	255	! --------------------------- \|
	256	! \|
	257	! ------F_lmin+1>0----------- \|
	258	! lmin \|
	259	! ----- F_lmin=0 ------------ _/
[2127]	260	!
	261	! ---------------------------
	262	!
	263	! ===========================
	264	! bottom of the model
[2060]	265
[2127]	266	!-------------------------------------------------------------------------------
	267	! Zoom on layers k and k-1
	268	!-------------------------------------------------------------------------------
[2060]	269
[2127]	270	! \| /\|\ \| \|
	271	! \|---- \| F_k+1 -----------\|--------------------------\| level k+1
	272	! \| \| w_k+1 \| \|
	273	! \| --\|--> D_k \|
	274	! \| \| \| layer k
	275	! \| <--\|-- E_k \|
	276	! \| /\|\ \| \|
	277	! \|---- \| F_k ----------\|-----------------------------\| level k
	278	! \| \| w_k \| \|
	279	! \| --\|--> D_k-1 \|
	280	! \| \| \| layer k-1
	281	! \| <--\|-- E_k-1 \|
	282	! \| /\|\ \| \|
	283	! \|---- \| F_k-1 -----\|--------------------------------\| level k-1
	284	! \| w_k-1
	285	! 0 fraca 1
	286	! \__________________/ \______________________________/
	287	! plume (fraca) environment (1-fraca)
[2060]	288
[2127]	289	!===============================================================================
	290	! Thermal plumes computation
	291	!===============================================================================
[2060]	292
[2232]	293	DO WHILE (re_tpm.and.(while_loop_counter<nlay))
	294	while_loop_counter = while_loop_counter + 1
	295
[2127]	296	!-------------------------------------------------------------------------------
[2232]	297	! Thermal plumes speeds, normalized fluxes, tracers and temperatures
[2127]	298	!-------------------------------------------------------------------------------
[2232]	299
	300	CALL thermcell_plume(ngrid,nlay,nq,ptimestep, &
	301	& ztv,zhl,zqt,zql,zlev,pplev,zpopsk, &
	302	& detr_star,entr_star,f_star, &
	303	& ztva,zhla,zqta,zqla,zqsa, &
	304	& zw2,lbot,lmin)
	305
[2127]	306	!-------------------------------------------------------------------------------
[2101]	307	! Thermal plumes characteristics: zmax, zmix, wmax
[2127]	308	!-------------------------------------------------------------------------------
[2232]	309
[2127]	310	! AB: Careful, zw2 became its square root in thermcell_height!
[2232]	311	CALL thermcell_height(ngrid,nlay,lmin,lmix,lmax, &
	312	& zlev,zmin,zmix,zmax,zw2,wmax,f_star)
	313
[2127]	314	!-------------------------------------------------------------------------------
	315	! Closure
	316	!-------------------------------------------------------------------------------
[2232]	317
	318	CALL thermcell_closure(ngrid,nlay,ptimestep,rho,zlev, &
	319	& lmax,entr_star,zmin,zmax,wmax,f)
	320
[2177]	321	! FH: Test valable seulement en 1D mais pas genant
[2232]	322	IF (.not. (f(1).ge.0.) ) THEN
	323	print *, 'ERROR: mass flux norm is not positive!'
	324	print *, 'f =', f(1)
	325	CALL abort
	326	ENDIF
	327
[2127]	328	!-------------------------------------------------------------------------------
[2101]	329	! Mass fluxes
[2127]	330	!-------------------------------------------------------------------------------
[2232]	331
	332	CALL thermcell_flux(ngrid,nlay,ptimestep,masse, &
	333	& lmin,lmax,entr_star,detr_star, &
	334	& f,rhobarz,zlev,zw2,fm,entr,detr)
	335
[2127]	336	!-------------------------------------------------------------------------------
[2232]	337	!
[2127]	338	!-------------------------------------------------------------------------------
[2232]	339
	340	re_tpm = .false.
	341	DO ig=1,ngrid
	342	IF(lmax(ig) > lmin(ig)) THEN
	343	lbot(ig) = lmax(ig)
	344	re_tpm = .true.
	345	ELSE
	346	lbot(ig) = nlay
	347	ENDIF
	348	ENDDO
	349
	350	!-------------------------------------------------------------------------------
	351	! Thermal plumes stacking
	352	!-------------------------------------------------------------------------------
	353
	354	zw2_tot(:,:) = zw2_tot(:,:) + zw2(:,:)
	355	entr_tot(:,:) = entr_tot(:,:) + entr(:,:)
	356	detr_tot(:,:) = detr_tot(:,:) + detr(:,:)
	357	fm_tot(:,:) = fm_tot(:,:) + fm(:,:)
	358
	359	ENDDO
[2101]	360
[2232]	361	!===============================================================================
	362	! Updraft fraction
	363	!===============================================================================
	364
[2101]	365	DO ig=1,ngrid
	366	fraca(ig,1) = 0.
	367	fraca(ig,nlay+1) = 0.
	368	ENDDO
	369
	370	DO l=2,nlay
	371	DO ig=1,ngrid
[2232]	372	IF (zw2_tot(ig,l) > 0.) THEN
	373	fraca(ig,l) = fm_tot(ig,l) / (rhobarz(ig,l) * zw2_tot(ig,l))
[2101]	374	ELSE
	375	fraca(ig,l) = 0.
	376	ENDIF
	377	ENDDO
	378	ENDDO
	379
[2127]	380	!===============================================================================
[2232]	381	! Vertical transport
[2127]	382	!===============================================================================
[2060]	383
[2127]	384	!-------------------------------------------------------------------------------
	385	! Calcul du transport vertical de la temperature potentielle
	386	!-------------------------------------------------------------------------------
[2101]	387
[2232]	388	CALL thermcell_dq(ngrid,nlay,ptimestep,fm_tot,entr_tot,detr_tot, &
	389	& masse,zhl,zdthladj,dummy)
[2060]	390
	391	DO l=1,nlay
	392	DO ig=1,ngrid
[2127]	393	pdtadj(ig,l) = zdthladj(ig,l) * zpopsk(ig,l)
[2060]	394	ENDDO
	395	ENDDO
	396
[2127]	397	!-------------------------------------------------------------------------------
	398	! Calcul du transport vertical des traceurs
	399	!-------------------------------------------------------------------------------
[2060]	400
[2127]	401	DO iq=1,nq
[2232]	402	CALL thermcell_dq(ngrid,nlay,ptimestep,fm_tot,entr_tot,detr_tot, &
	403	& masse,pq(:,:,iq),pdqadj(:,:,iq),zqa(:,:,iq))
[2102]	404	ENDDO
[2060]	405
[2127]	406	!-------------------------------------------------------------------------------
	407	! Calcul du transport vertical du moment horizontal
	408	!-------------------------------------------------------------------------------
[2060]	409
[2232]	410	! AB: Careful, thermcell_dv2 wasn't checked! It is not sure that it works
	411	! correctly with the plumes stacking (zmin, wmax doesn't make sense).
[2144]	412	IF (dvimpl) THEN
[2232]	413	CALL thermcell_dv2(ngrid,nlay,ptimestep,fm_tot,entr_tot,detr_tot, &
	414	& masse,fraca,zmax,zmin,pu,pv,pduadj,pdvadj,zua,zva)
[2144]	415	ELSE
[2232]	416	CALL thermcell_dq(ngrid,nlay,ptimestep,fm_tot,entr_tot,detr_tot, &
	417	& masse,zu,pduadj,zua)
	418	CALL thermcell_dq(ngrid,nlay,ptimestep,fm_tot,entr_tot,detr_tot, &
	419	& masse,zv,pdvadj,zva)
[2144]	420	ENDIF
[2060]	421
	422
	423	RETURN
	424	END

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

Download in other formats:

Original Format