/[lmdze]/trunk/phylmd/Interface_surf/calcul_fluxs.f
ViewVC logotype

Annotation of /trunk/phylmd/Interface_surf/calcul_fluxs.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 104 - (hide annotations)
Thu Sep 4 10:05:52 2014 UTC (9 years, 8 months ago) by guez
File size: 6953 byte(s) 
Removed procedure sortvarc0. Called sortvarc with an additional
argument resetvarc instead. (Following LMDZ.) Moved current time
computations and some printing statements from sortvarc to
caldyn. Could then remove arguments itau and time_0 of sortvarc, and
could remove "use dynetat0". Better to keep "dynetat0.f" as a gcm-only
file.

Moved some variables from module ener to module sortvarc.

Split file "mathelp.f" into single-procedure files.

Removed unused argument nadv of adaptdt. Removed dimension arguments
of bernoui.

Removed unused argument nisurf of interfoce_lim. Changed the size of
argument lmt_sst of interfoce_lim from klon to knon. Removed case when
newlmt is false.

dynredem1 is called only once in each run, either ce0l or gcm. So
variable nb in call to nf95_put_var was always 1. Removed variable nb.

Removed dimension arguments of calcul_fluxs. Removed unused arguments
precip_rain, precip_snow, snow of calcul_fluxs. Changed the size of
all the arrays in calcul_fluxs from klon to knon.

Removed dimension arguments of fonte_neige. Changed the size of all
the arrays in fonte_neige from klon to knon.

Changed the size of arguments tsurf and tsurf_new of interfsurf_hq
from klon to knon. Changed the size of argument ptsrf of soil from
klon to knon.
1 guez 54 module calcul_fluxs_m
2    
3     implicit none
4    
5     contains
6    
7 guez 104 SUBROUTINE calcul_fluxs(nisurf, dtime, tsurf, p1lay, cal, beta, coef1lay, &
8     ps, qsurf, radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, petAcoef, &
9     peqAcoef, petBcoef, peqBcoef, tsurf_new, evap, fluxlat, fluxsens, &
10     dflux_s, dflux_l)
11 guez 54
12 guez 99 ! Cette routine calcule les fluxs en h et q à l'interface et une
13     ! température de surface.
14 guez 54
15 guez 104 ! L. Fairhead April 2000
16 guez 54
17 guez 104 USE abort_gcm_m, ONLY: abort_gcm
18     USE indicesol, ONLY: is_ter
19     USE fcttre, ONLY: dqsatl, dqsats, foede, foeew, qsatl, qsats, thermcep
20     USE interface_surf, ONLY: run_off
21     use nr_util, only: assert_eq
22     USE suphec_m, ONLY: rcpd, rd, retv, rkappa, rlstt, rlvtt, rtt
23     USE yoethf_m, ONLY: r2es, r5ies, r5les, rvtmp2
24    
25     integer, intent(IN):: nisurf ! surface a traiter
26     real, intent(IN):: dtime
27     real, intent(IN):: tsurf(:) ! (knon) temperature de surface
28     real, intent(IN):: p1lay(:) ! (knon) pression 1er niveau (milieu de couche)
29     real, intent(IN):: cal(:) ! (knon) capacité calorifique du sol
30     real, intent(IN):: beta(:) ! (knon) evap reelle
31     real, intent(IN):: coef1lay(:) ! (knon) coefficient d'échange
32     real, intent(IN):: ps(:) ! (knon) pression au sol
33     real, intent(OUT):: qsurf(:) ! (knon) humidite de l'air au dessus du sol
34     real, intent(IN):: radsol(:) ! (knon) rayonnement net au sol (LW + SW)
35    
36     real, intent(IN):: dif_grnd(:) ! (knon)
37     ! coefficient diffusion vers le sol profond
38    
39     real, intent(IN):: t1lay(:), q1lay(:), u1lay(:), v1lay(:) ! (knon)
40    
41     real, intent(IN):: petAcoef(:), peqAcoef(:) ! (knon)
42     ! coefficients A de la résolution de la couche limite pour t et q
43    
44     real, intent(IN):: petBcoef(:), peqBcoef(:) ! (knon)
45 guez 54 ! petBcoef coeff. B de la resolution de la CL pour t
46     ! peqBcoef coeff. B de la resolution de la CL pour q
47    
48 guez 104 real, intent(OUT):: tsurf_new(:) ! (knon) température au sol
49     real, intent(OUT):: evap(:), fluxlat(:), fluxsens(:) ! (knon)
50     ! fluxlat flux de chaleur latente
51 guez 54 ! fluxsens flux de chaleur sensible
52 guez 104 real, intent(OUT):: dflux_s(:), dflux_l(:) ! (knon)
53     ! Dérivées des flux dF/dTs (W m-2 K-1)
54 guez 54 ! dflux_s derivee du flux de chaleur sensible / Ts
55     ! dflux_l derivee du flux de chaleur latente / Ts
56    
57 guez 104 ! Local:
58     integer i
59     real, dimension(size(ps)) :: zx_mh, zx_nh, zx_oh
60     real, dimension(size(ps)) :: zx_mq, zx_nq, zx_oq
61     real, dimension(size(ps)) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
62     real, dimension(size(ps)) :: zx_sl, zx_k1
63     real, dimension(size(ps)) :: zx_q_0 , d_ts
64 guez 103 logical zdelta
65     real zcvm5, zx_qs, zcor, zx_dq_s_dh
66 guez 54 real :: bilan_f, fq_fonte
67     REAL :: subli, fsno
68     REAL :: qsat_new, q1_new
69 guez 104 integer knon ! nombre de points a traiter
70     real, parameter:: t_grnd = 271.35, t_coup = 273.15
71 guez 54
72 guez 104 !---------------------------------------------------------------------
73 guez 54
74 guez 104 knon = assert_eq((/size(tsurf), size(p1lay), size(cal), size(beta), &
75     size(coef1lay), size(ps), size(qsurf), size(radsol), size(dif_grnd), &
76     size(t1lay), size(q1lay), size(u1lay), size(v1lay), size(petAcoef), &
77     size(peqAcoef), size(petBcoef), size(peqBcoef), size(tsurf_new), &
78     size(evap), size(fluxlat), size(fluxsens), size(dflux_s), &
79     size(dflux_l)/), "calcul_fluxs knon")
80 guez 54
81 guez 101 if (size(run_off) /= knon .AND. nisurf == is_ter) then
82 guez 104 print *, 'Bizarre, le nombre de points continentaux'
83     print *, 'a change entre deux appels. J''arrete.'
84     call abort_gcm('calcul_fluxs', 'Pb run_off', 1)
85 guez 54 endif
86    
87 guez 104 ! Traitement humidite du sol
88 guez 54
89     evap = 0.
90     fluxsens=0.
91     fluxlat=0.
92     dflux_s = 0.
93     dflux_l = 0.
94    
95     ! zx_qs = qsat en kg/kg
96    
97     DO i = 1, knon
98     zx_pkh(i) = (ps(i)/ps(i))**RKAPPA
99     IF (thermcep) THEN
100 guez 103 zdelta= rtt >= tsurf(i)
101     zcvm5 = merge(R5IES*RLSTT, R5LES*RLVTT, zdelta)
102 guez 54 zcvm5 = zcvm5 / RCPD / (1.0+RVTMP2*q1lay(i))
103     zx_qs= r2es * FOEEW(tsurf(i), zdelta)/ps(i)
104     zx_qs=MIN(0.5, zx_qs)
105     zcor=1./(1.-retv*zx_qs)
106     zx_qs=zx_qs*zcor
107     zx_dq_s_dh = FOEDE(tsurf(i), zdelta, zcvm5, zx_qs, zcor) &
108     /RLVTT / zx_pkh(i)
109     ELSE
110     IF (tsurf(i).LT.t_coup) THEN
111     zx_qs = qsats(tsurf(i)) / ps(i)
112     zx_dq_s_dh = dqsats(tsurf(i), zx_qs)/RLVTT &
113     / zx_pkh(i)
114     ELSE
115     zx_qs = qsatl(tsurf(i)) / ps(i)
116     zx_dq_s_dh = dqsatl(tsurf(i), zx_qs)/RLVTT &
117     / zx_pkh(i)
118     ENDIF
119     ENDIF
120     zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh
121     zx_qsat(i) = zx_qs
122     zx_coef(i) = coef1lay(i) &
123     * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) &
124     * p1lay(i)/(RD*t1lay(i))
125    
126     ENDDO
127    
128     ! === Calcul de la temperature de surface ===
129    
130     ! zx_sl = chaleur latente d'evaporation ou de sublimation
131    
132     do i = 1, knon
133     zx_sl(i) = RLVTT
134     if (tsurf(i) .LT. RTT) zx_sl(i) = RLSTT
135     zx_k1(i) = zx_coef(i)
136     enddo
137    
138     do i = 1, knon
139     ! Q
140     zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)
141     zx_mq(i) = beta(i) * zx_k1(i) * &
142     (peqAcoef(i) - zx_qsat(i) &
143     + zx_dq_s_dt(i) * tsurf(i)) &
144     / zx_oq(i)
145     zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &
146     / zx_oq(i)
147    
148     ! H
149     zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)
150     zx_mh(i) = zx_k1(i) * petAcoef(i) / zx_oh(i)
151     zx_nh(i) = - (zx_k1(i) * RCPD * zx_pkh(i))/ zx_oh(i)
152    
153     ! Tsurface
154     tsurf_new(i) = (tsurf(i) + cal(i)/(RCPD * zx_pkh(i)) * dtime * &
155     (radsol(i) + zx_mh(i) + zx_sl(i) * zx_mq(i)) &
156     + dif_grnd(i) * t_grnd * dtime)/ &
157     ( 1. - dtime * cal(i)/(RCPD * zx_pkh(i)) * ( &
158     zx_nh(i) + zx_sl(i) * zx_nq(i)) &
159     + dtime * dif_grnd(i))
160    
161    
162     ! Y'a-t-il fonte de neige?
163    
164     ! fonte_neige = (nisurf /= is_oce) .AND. &
165     ! & (snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &
166     ! & .AND. (tsurf_new(i) >= RTT)
167     ! if (fonte_neige) tsurf_new(i) = RTT
168     d_ts(i) = tsurf_new(i) - tsurf(i)
169     ! zx_h_ts(i) = tsurf_new(i) * RCPD * zx_pkh(i)
170     ! zx_q_0(i) = zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
171     !== flux_q est le flux de vapeur d'eau: kg/(m**2 s) positive vers bas
172     !== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
173     evap(i) = - zx_mq(i) - zx_nq(i) * tsurf_new(i)
174     fluxlat(i) = - evap(i) * zx_sl(i)
175     fluxsens(i) = zx_mh(i) + zx_nh(i) * tsurf_new(i)
176     ! Derives des flux dF/dTs (W m-2 K-1):
177     dflux_s(i) = zx_nh(i)
178     dflux_l(i) = (zx_sl(i) * zx_nq(i))
179     ! Nouvelle valeure de l'humidite au dessus du sol
180     qsat_new=zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
181     q1_new = peqAcoef(i) - peqBcoef(i)*evap(i)*dtime
182     qsurf(i)=q1_new*(1.-beta(i)) + beta(i)*qsat_new
183     ENDDO
184    
185     END SUBROUTINE calcul_fluxs
186    
187     end module calcul_fluxs_m
  ViewVC Help
Powered by ViewVC 1.1.21