Context Navigation

surface.F90 @ 4221

Last change on this file since 4221 was 4210, checked in by dubos, 6 years ago
simple_physics : cleanup PRINTs
File size: 8.5 KB

Line
1	MODULE surface
2
3	#include "use_logging.h"
4
5	IMPLICIT NONE
6	PRIVATE
7	SAVE
8
9	REAL, PARAMETER :: pi=2.*ASIN(1.)
10
11	! common variables
12	REAL, PUBLIC :: I_mer,I_ter,Cd_mer,Cd_ter, &
13	& alb_mer,alb_ter,emi_mer,emi_ter
14
15	! local saved variables:
16	! ----------------------
17	REAL :: lambda
18	REAL,ALLOCATABLE :: dz1(:),dz2(:),zc(:,:),zd(:,:)
19	!$OMP THREADPRIVATE(dz1,dz2,zc,zd,lambda)
20
21	PUBLIC :: soil
22
23	CONTAINS
24
25	SUBROUTINE init_soil(ngrid,nsoil)
26	INTEGER, INTENT(IN) :: ngrid, nsoil
27	REAL min_period,dalph_soil
28	REAL fz,rk,fz1,rk1,rk2
29	INTEGER :: jk
30
31	! this is a function definition
32	fz(rk)=fz1(dalph_soil*rk-1.)/(dalph_soil-1.)
33
34	!-----------------------------------------------------------------------
35	! ground levels
36	! grnd=z/l where l is the skin depth of the diurnal cycle:
37	! --------------------------------------------------------
38
39	WRITELOG(,) 'nsoil,ngrid,firstcall=',nsoil,ngrid, .TRUE.
40
41	ALLOCATE(dz1(nsoil),dz2(nsoil))
42	ALLOCATE(zc(ngrid,nsoil),zd(ngrid,nsoil))
43
44	min_period=20000.
45	dalph_soil=2.
46
47	! la premiere couche represente un dixieme de cycle diurne
48	fz1=sqrt(min_period/pi)
49
50	DO jk=1,nsoil
51	rk1=jk
52	rk2=jk-1
53	dz2(jk)=fz(rk1)-fz(rk2)
54	ENDDO
55	DO jk=1,nsoil-1
56	rk1=jk+.5
57	rk2=jk-.5
58	dz1(jk)=1./(fz(rk1)-fz(rk2))
59	ENDDO
60	lambda=fz(.5)*dz1(1)
61	WRITELOG(,) 'full layers, intermediate layers (secoonds)'
62	DO jk=1,nsoil
63	rk=jk
64	rk1=jk+.5
65	rk2=jk-.5
66	WRITELOG(,) fz(rk1)fz(rk2)pi, &
67	& fz(rk)fz(rk)pi
68	ENDDO
69
70	LOG_INFO('init_soil')
71	END SUBROUTINE init_soil
72
73	SUBROUTINE soil(ngrid,nsoil,firstcall,ptherm_i, &
74	& ptimestep,ptsrf,ptsoil, &
75	& pcapcal,pfluxgrd)
76
77	!=======================================================================
78	!
79	! Auteur: Frederic Hourdin 30/01/92
80	! -------
81	!
82	! objet: computation of : the soil temperature evolution
83	! ------ the surfacic heat capacity "Capcal"
84	! the surface conduction flux pcapcal
85	!
86	!
87	! Method: implicit time integration
88	! -------
89	! Consecutive ground temperatures are related by:
90	! T(k+1) = C(k) + D(k)*T(k) (1)
91	! the coefficients C and D are computed at the t-dt time-step.
92	! Routine structure:
93	! 1)new temperatures are computed using (1)
94	! 2)C and D coefficients are computed from the new temperature
95	! profile for the t+dt time-step
96	! 3)the coefficients A and B are computed where the diffusive
97	! fluxes at the t+dt time-step is given by
98	! Fdiff = A + B Ts(t+dt)
99	! or Fdiff = F0 + Capcal (Ts(t+dt)-Ts(t))/dt
100	! with F0 = A + B (Ts(t))
101	! Capcal = B*dt
102	!
103	! Interface:
104	! ----------
105	!
106	! Arguments:
107	! ----------
108	! ngrid number of grid-points
109	! ptimestep physical timestep (s)
110	! pto(ngrid,nsoil) temperature at time-step t (K)
111	! ptn(ngrid,nsoil) temperature at time step t+dt (K)
112	! pcapcal(ngrid) specific heat (Wm-2s*K-1)
113	! pfluxgrd(ngrid) surface diffusive flux from ground (Wm-2)
114	!
115	!=======================================================================
116	! declarations:
117	! -------------
118
119
120	!-----------------------------------------------------------------------
121	! arguments
122	! ---------
123
124	INTEGER ngrid,nsoil
125	REAL ptimestep
126	REAL ptsrf(ngrid),ptsoil(ngrid,nsoil),ptherm_i(ngrid)
127	REAL pcapcal(ngrid),pfluxgrd(ngrid)
128	LOGICAL firstcall
129
130
131	!-----------------------------------------------------------------------
132	! local arrays
133	! ------------
134
135	INTEGER ig,jk
136	REAL za(ngrid),zb(ngrid)
137	REAL zdz2(nsoil),z1(ngrid)
138
139	IF (firstcall) THEN
140	CALL init_soil(ngrid, nsoil)
141	ELSE
142	!-----------------------------------------------------------------------
143	! Computation of the soil temperatures using the Cgrd and Dgrd
144	! coefficient computed at the previous time-step:
145	! -----------------------------------------------
146
147	! surface temperature
148	DO ig=1,ngrid
149	ptsoil(ig,1)=(lambda*zc(ig,1)+ptsrf(ig))/ &
150	& (lambda*(1.-zd(ig,1))+1.)
151	ENDDO
152
153	! other temperatures
154	DO jk=1,nsoil-1
155	DO ig=1,ngrid
156	ptsoil(ig,jk+1)=zc(ig,jk)+zd(ig,jk)*ptsoil(ig,jk)
157	ENDDO
158	ENDDO
159
160	ENDIF
161
162	!-----------------------------------------------------------------------
163	! Computation of the Cgrd and Dgrd coefficient for the next step:
164	! ---------------------------------------------------------------
165
166	DO jk=1,nsoil
167	zdz2(jk)=dz2(jk)/ptimestep
168	ENDDO
169
170	DO ig=1,ngrid
171	z1(ig)=zdz2(nsoil)+dz1(nsoil-1)
172	zc(ig,nsoil-1)=zdz2(nsoil)*ptsoil(ig,nsoil)/z1(ig)
173	zd(ig,nsoil-1)=dz1(nsoil-1)/z1(ig)
174	ENDDO
175
176	DO jk=nsoil-1,2,-1
177	DO ig=1,ngrid
178	z1(ig)=1./(zdz2(jk)+dz1(jk-1)+dz1(jk)*(1.-zd(ig,jk)))
179	zc(ig,jk-1)= &
180	& (ptsoil(ig,jk)zdz2(jk)+dz1(jk)zc(ig,jk))*z1(ig)
181	zd(ig,jk-1)=dz1(jk-1)*z1(ig)
182	ENDDO
183	ENDDO
184
185	!-----------------------------------------------------------------------
186	! computation of the surface diffusive flux from ground and
187	! calorific capacity of the ground:
188	! ---------------------------------
189
190	DO ig=1,ngrid
191	pfluxgrd(ig)=ptherm_i(ig)dz1(1) &
192	& (zc(ig,1)+(zd(ig,1)-1.)*ptsoil(ig,1))
193	z1(ig)=lambda*(1.-zd(ig,1))+1.
194	pcapcal(ig)=ptherm_i(ig)* &
195	& ptimestep(zdz2(1)+(1.-zd(ig,1))dz1(1))/z1(ig)
196	pfluxgrd(ig)=pfluxgrd(ig) &
197	& +pcapcal(ig)(ptsoil(ig,1)z1(ig)-lambda*zc(ig,1)-ptsrf(ig)) &
198	& /ptimestep
199	ENDDO
200
201	END SUBROUTINE soil
202
203	END MODULE surface

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

Download in other formats:

Original Format