Context Navigation

yamada.F @ 965

Last change on this file since 965 was 868, checked in by Laurent Fairhead, 17 years ago
Preparation du remplacement de la physique utilisee pour l'exercice IPCC_AR4 par la version de la physique avec thermique. On garde le repertoire phylmd pour un petit moment pour que les utilisateurs ne soient pas trop perdus ... phy_IPCC_AR4 = phylmd LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 5.5 KB

Rev	Line
[868]	1	!
	2	! $Header$
	3	!
	4	SUBROUTINE yamada(ngrid,dt,g,rconst,plev,temp
	5	s ,zlev,zlay,u,v,teta,cd,q2,km,kn,ustar
	6	s ,l_mix)
	7	use dimphy
	8	IMPLICIT NONE
	9	c.......................................................................
	10	cym#include "dimensions.h"
	11	cym#include "dimphy.h"
	12	c.......................................................................
	13	c
	14	c dt : pas de temps
	15	c g : g
	16	c zlev : altitude a chaque niveau (interface inferieure de la couche
	17	c de meme indice)
	18	c zlay : altitude au centre de chaque couche
	19	c u,v : vitesse au centre de chaque couche
	20	c (en entree : la valeur au debut du pas de temps)
	21	c teta : temperature potentielle au centre de chaque couche
	22	c (en entree : la valeur au debut du pas de temps)
	23	c cd : cdrag
	24	c (en entree : la valeur au debut du pas de temps)
	25	c q2 : $q^2$ au bas de chaque couche
	26	c (en entree : la valeur au debut du pas de temps)
	27	c (en sortie : la valeur a la fin du pas de temps)
	28	c km : diffusivite turbulente de quantite de mouvement (au bas de chaque
	29	c couche)
	30	c (en sortie : la valeur a la fin du pas de temps)
	31	c kn : diffusivite turbulente des scalaires (au bas de chaque couche)
	32	c (en sortie : la valeur a la fin du pas de temps)
	33	c
	34	c.......................................................................
	35	REAL dt,g,rconst
	36	real plev(klon,klev+1),temp(klon,klev)
	37	real ustar(klon),snstable
	38	REAL zlev(klon,klev+1)
	39	REAL zlay(klon,klev)
	40	REAL u(klon,klev)
	41	REAL v(klon,klev)
	42	REAL teta(klon,klev)
	43	REAL cd(klon)
	44	REAL q2(klon,klev+1)
	45	REAL km(klon,klev+1)
	46	REAL kn(klon,klev+1)
	47	integer l_mix,ngrid
	48
	49
	50	integer nlay,nlev
	51	cym PARAMETER (nlay=klev)
	52	cym PARAMETER (nlev=klev+1)
	53
	54	logical first
	55	save first
	56	data first/.true./
	57
	58
	59	integer ig,k
	60
	61	real ri,zrif,zalpha,zsm
	62	real rif(klon,klev+1),sm(klon,klev+1),alpha(klon,klev)
	63
	64	real m2(klon,klev+1),dz(klon,klev+1),zq,n2(klon,klev+1)
	65	real l(klon,klev+1),l0(klon)
	66
	67	real sq(klon),sqz(klon),zz(klon,klev+1)
	68	integer iter
	69
	70	real ric,rifc,b1,kap
	71	save ric,rifc,b1,kap
	72	data ric,rifc,b1,kap/0.195,0.191,16.6,0.3/
	73
	74	real frif,falpha,fsm
	75
	76	frif(ri)=0.6588(ri+0.1776-sqrt(riri-0.3221*ri+0.03156))
	77	falpha(ri)=1.318*(0.2231-ri)/(0.2341-ri)
	78	fsm(ri)=1.96(0.1912-ri)(0.2341-ri)/((1.-ri)*(0.2231-ri))
	79
	80	nlay=klev
	81	nlev=klev+1
	82
	83	if (0.eq.1.and.first) then
	84	do ig=1,1000
	85	ri=(ig-800.)/500.
	86	if (ri.lt.ric) then
	87	zrif=frif(ri)
	88	else
	89	zrif=rifc
	90	endif
	91	if(zrif.lt.0.16) then
	92	zalpha=falpha(zrif)
	93	zsm=fsm(zrif)
	94	else
	95	zalpha=1.12
	96	zsm=0.085
	97	endif
	98	print*,ri,rif,zalpha,zsm
	99	enddo
	100	first=.false.
	101	endif
	102
	103	c Correction d'un bug sauvage a verifier.
	104	c do k=2,nlev
	105	do k=2,nlay
	106	do ig=1,ngrid
	107	dz(ig,k)=zlay(ig,k)-zlay(ig,k-1)
	108	m2(ig,k)=((u(ig,k)-u(ig,k-1))2+(v(ig,k)-v(ig,k-1))2)
	109	s /(dz(ig,k)*dz(ig,k))
	110	n2(ig,k)=g2.(teta(ig,k)-teta(ig,k-1))
	111	s /(teta(ig,k-1)+teta(ig,k)) /dz(ig,k)
	112	ri=n2(ig,k)/max(m2(ig,k),1.e-10)
	113	if (ri.lt.ric) then
	114	rif(ig,k)=frif(ri)
	115	else
	116	rif(ig,k)=rifc
	117	endif
	118	if(rif(ig,k).lt.0.16) then
	119	alpha(ig,k)=falpha(rif(ig,k))
	120	sm(ig,k)=fsm(rif(ig,k))
	121	else
	122	alpha(ig,k)=1.12
	123	sm(ig,k)=0.085
	124	endif
	125	zz(ig,k)=b1m2(ig,k)(1.-rif(ig,k))*sm(ig,k)
	126	enddo
	127	enddo
	128
	129	c iterration pour determiner la longueur de melange
	130
	131	do ig=1,ngrid
	132	l0(ig)=100.
	133	enddo
	134	do k=2,klev-1
	135	do ig=1,ngrid
	136	l(ig,k)=l0(ig)kapzlev(ig,k)/(kap*zlev(ig,k)+l0(ig))
	137	enddo
	138	enddo
	139
	140	do iter=1,10
	141	do ig=1,ngrid
	142	sq(ig)=1.e-10
	143	sqz(ig)=1.e-10
	144	enddo
	145	do k=2,klev-1
	146	do ig=1,ngrid
	147	q2(ig,k)=l(ig,k)*2zz(ig,k)
	148	l(ig,k)=min(l0(ig)kapzlev(ig,k)/(kap*zlev(ig,k)+l0(ig))
	149	s ,0.5*sqrt(q2(ig,k))/sqrt(max(n2(ig,k),1.e-10)))
	150	zq=sqrt(q2(ig,k))
	151	sqz(ig)=sqz(ig)+zqzlev(ig,k)(zlay(ig,k)-zlay(ig,k-1))
	152	sq(ig)=sq(ig)+zq*(zlay(ig,k)-zlay(ig,k-1))
	153	enddo
	154	enddo
	155	do ig=1,ngrid
	156	l0(ig)=0.2*sqz(ig)/sq(ig)
	157	enddo
	158	c(abd 3 5 2) print*,'ITER=',iter,' L0=',l0
	159
	160	enddo
	161
	162	do k=2,klev
	163	do ig=1,ngrid
	164	l(ig,k)=min(l0(ig)kapzlev(ig,k)/(kap*zlev(ig,k)+l0(ig))
	165	s ,0.5*sqrt(q2(ig,k))/sqrt(max(n2(ig,k),1.e-10)))
	166	q2(ig,k)=l(ig,k)*2zz(ig,k)
	167	km(ig,k)=l(ig,k)sqrt(q2(ig,k))sm(ig,k)
	168	kn(ig,k)=km(ig,k)*alpha(ig,k)
	169	enddo
	170	enddo
	171
	172	return
	173	end

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

Download in other formats:

Original Format