Context Navigation

iniacademic.F @ 1455

Last change on this file since 1455 was 1454, checked in by Laurent Fairhead, 14 years ago

Merge of LMDZ5V1.0-dev branch r1453 into LMDZ5 trunk r1434

Fusion entre la version r1453 de la branche de développement LMDZ5V1.0-dev
et le tronc LMDZ5 (r1434)

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

File size: 8.4 KB

Line
1	!
2	! $Id: iniacademic.F 1454 2010-11-18 12:01:24Z fairhead $
3	!
4	c
5	c
6	SUBROUTINE iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
7
8	USE filtreg_mod
9	USE infotrac, ONLY : nqtot
10	USE control_mod, ONLY: day_step,planet_type
11	#ifdef CPP_IOIPSL
12	USE IOIPSL
13	#else
14	! if not using IOIPSL, we still need to use (a local version of) getin
15	USE ioipsl_getincom
16	#endif
17	USE Write_Field
18
19
20	c%W% %G%
21	c=======================================================================
22	c
23	c Author: Frederic Hourdin original: 15/01/93
24	c -------
25	c
26	c Subject:
27	c ------
28	c
29	c Method:
30	c --------
31	c
32	c Interface:
33	c ----------
34	c
35	c Input:
36	c ------
37	c
38	c Output:
39	c -------
40	c
41	c=======================================================================
42	IMPLICIT NONE
43	c-----------------------------------------------------------------------
44	c Declararations:
45	c ---------------
46
47	#include "dimensions.h"
48	#include "paramet.h"
49	#include "comvert.h"
50	#include "comconst.h"
51	#include "comgeom.h"
52	#include "academic.h"
53	#include "ener.h"
54	#include "temps.h"
55	#include "iniprint.h"
56	#include "logic.h"
57
58	c Arguments:
59	c ----------
60
61	real time_0
62
63	c variables dynamiques
64	REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
65	REAL teta(ip1jmp1,llm) ! temperature potentielle
66	REAL q(ip1jmp1,llm,nqtot) ! champs advectes
67	REAL ps(ip1jmp1) ! pression au sol
68	REAL masse(ip1jmp1,llm) ! masse d'air
69	REAL phis(ip1jmp1) ! geopotentiel au sol
70
71	c Local:
72	c ------
73
74	REAL p (ip1jmp1,llmp1 ) ! pression aux interfac.des couches
75	REAL pks(ip1jmp1) ! exner au sol
76	REAL pk(ip1jmp1,llm) ! exner au milieu des couches
77	REAL pkf(ip1jmp1,llm) ! exner filt.au milieu des couches
78	REAL phi(ip1jmp1,llm) ! geopotentiel
79	REAL ddsin,tetastrat,zsig,tetapv,w_pv ! variables auxiliaires
80	real tetajl(jjp1,llm)
81	INTEGER i,j,l,lsup,ij
82
83	REAL teta0,ttp,delt_y,delt_z,eps ! Constantes pour profil de T
84	REAL k_f,k_c_a,k_c_s ! Constantes de rappel
85	LOGICAL ok_geost ! Initialisation vent geost. ou nul
86	LOGICAL ok_pv ! Polar Vortex
87	REAL phi_pv,dphi_pv,gam_pv ! Constantes pour polar vortex
88
89	real zz,ran1
90	integer idum
91
92	REAL alpha(ip1jmp1,llm),beta(ip1jmp1,llm),zdtvr
93
94	c-----------------------------------------------------------------------
95	! 1. Initializations for Earth-like case
96	! --------------------------------------
97	c
98	! initialize planet radius, rotation rate,...
99	call conf_planete
100
101	time_0=0.
102	day_ref=1
103	annee_ref=0
104
105	im = iim
106	jm = jjm
107	day_ini = 1
108	dtvr = daysec/REAL(day_step)
109	zdtvr=dtvr
110	etot0 = 0.
111	ptot0 = 0.
112	ztot0 = 0.
113	stot0 = 0.
114	ang0 = 0.
115
116	if (llm.eq.1) then
117	! specific initializations for the shallow water case
118	kappa=1
119	endif
120
121	CALL iniconst
122	CALL inigeom
123	CALL inifilr
124
125	if (llm.eq.1) then
126	! initialize fields for the shallow water case, if required
127	if (.not.read_start) then
128	phis(:)=0.
129	q(:,:,:)=0
130	CALL sw_case_williamson91_6(vcov,ucov,teta,masse,ps)
131	endif
132	endif
133
134	if (iflag_phys.eq.2) then
135	! initializations for the academic case
136
137	! if (planet_type=="earth") then
138
139	! 1. local parameters
140	! by convention, winter is in the southern hemisphere
141	! Geostrophic wind or no wind?
142	ok_geost=.TRUE.
143	CALL getin('ok_geost',ok_geost)
144	! Constants for Newtonian relaxation and friction
145	k_f=1. !friction
146	CALL getin('k_j',k_f)
147	k_f=1./(daysec*k_f)
148	k_c_s=4. !cooling surface
149	CALL getin('k_c_s',k_c_s)
150	k_c_s=1./(daysec*k_c_s)
151	k_c_a=40. !cooling free atm
152	CALL getin('k_c_a',k_c_a)
153	k_c_a=1./(daysec*k_c_a)
154	! Constants for Teta equilibrium profile
155	teta0=315. ! mean Teta (S.H. 315K)
156	CALL getin('teta0',teta0)
157	ttp=200. ! Tropopause temperature (S.H. 200K)
158	CALL getin('ttp',ttp)
159	eps=0. ! Deviation to N-S symmetry(~0-20K)
160	CALL getin('eps',eps)
161	delt_y=60. ! Merid Temp. Gradient (S.H. 60K)
162	CALL getin('delt_y',delt_y)
163	delt_z=10. ! Vertical Gradient (S.H. 10K)
164	CALL getin('delt_z',delt_z)
165	! Polar vortex
166	ok_pv=.false.
167	CALL getin('ok_pv',ok_pv)
168	phi_pv=-50. ! Latitude of edge of vortex
169	CALL getin('phi_pv',phi_pv)
170	phi_pv=phi_pv*pi/180.
171	dphi_pv=5. ! Width of the edge
172	CALL getin('dphi_pv',dphi_pv)
173	dphi_pv=dphi_pv*pi/180.
174	gam_pv=4. ! -dT/dz vortex (in K/km)
175	CALL getin('gam_pv',gam_pv)
176
177	! 2. Initialize fields towards which to relax
178	! Friction
179	knewt_g=k_c_a
180	DO l=1,llm
181	zsig=presnivs(l)/preff
182	knewt_t(l)=(k_c_s-k_c_a)*MAX(0.,(zsig-0.7)/0.3)
183	kfrict(l)=k_f*MAX(0.,(zsig-0.7)/0.3)
184	ENDDO
185	DO j=1,jjp1
186	clat4((j-1)iip1+1:jiip1)=cos(rlatu(j))**4
187	ENDDO
188
189	! Potential temperature
190	DO l=1,llm
191	zsig=presnivs(l)/preff
192	tetastrat=ttpzsig*(-kappa)
193	tetapv=tetastrat
194	IF ((ok_pv).AND.(zsig.LT.0.1)) THEN
195	tetapv=tetastrat(zsig10.)*(kappacpp*gam_pv/1000./g)
196	ENDIF
197	DO j=1,jjp1
198	! Troposphere
199	ddsin=sin(rlatu(j))
200	tetajl(j,l)=teta0-delt_yddsinddsin+eps*ddsin
201	& -delt_z(1.-ddsinddsin)*log(zsig)
202	! Profil stratospherique isotherme (+vortex)
203	w_pv=(1.-tanh((rlatu(j)-phi_pv)/dphi_pv))/2.
204	tetastrat=tetastrat(1.-w_pv)+tetapvw_pv
205	tetajl(j,l)=MAX(tetajl(j,l),tetastrat)
206	ENDDO
207	ENDDO ! of DO l=1,llm
208
209	! CALL writefield('theta_eq',tetajl)
210
211	do l=1,llm
212	do j=1,jjp1
213	do i=1,iip1
214	ij=(j-1)*iip1+i
215	tetarappel(ij,l)=tetajl(j,l)
216	enddo
217	enddo
218	enddo
219
220
221	! else
222	! write(lunout,*)"iniacademic: planet types other than earth",
223	! & " not implemented (yet)."
224	! stop
225	! endif ! of if (planet_type=="earth")
226
227	! 3. Initialize fields (if necessary)
228	IF (.NOT. read_start) THEN
229	! surface pressure
230	ps(:)=preff
231	! ground geopotential
232	phis(:)=0.
233
234	CALL pression ( ip1jmp1, ap, bp, ps, p )
235	CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
236	CALL massdair(p,masse)
237
238	! bulk initialization of temperature
239	teta(:,:)=tetarappel(:,:)
240
241	! geopotential
242	CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
243
244	! winds
245	if (ok_geost) then
246	call ugeostr(phi,ucov)
247	else
248	ucov(:,:)=0.
249	endif
250	vcov(:,:)=0.
251
252	! bulk initialization of tracers
253	if (planet_type=="earth") then
254	! Earth: first two tracers will be water
255	do i=1,nqtot
256	if (i.eq.1) q(:,:,i)=1.e-10
257	if (i.eq.2) q(:,:,i)=1.e-15
258	if (i.gt.2) q(:,:,i)=0.
259	enddo
260	else
261	q(:,:,:)=0
262	endif ! of if (planet_type=="earth")
263
264	! add random perturbation to temperature
265	idum = -1
266	zz = ran1(idum)
267	idum = 0
268	do l=1,llm
269	do ij=iip2,ip1jm
270	teta(ij,l)=teta(ij,l)(1.+0.005ran1(idum))
271	enddo
272	enddo
273
274	! maintain periodicity in longitude
275	do l=1,llm
276	do ij=1,ip1jmp1,iip1
277	teta(ij+iim,l)=teta(ij,l)
278	enddo
279	enddo
280
281	ENDIF ! of IF (.NOT. read_start)
282	endif ! of if (iflag_phys.eq.2)
283
284	END
285	c-----------------------------------------------------------------------

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

Download in other formats:

Original Format