source: LMDZ4/trunk/libf/phylmd/ajsec.F @ 2221

Last change on this file since 2221 was 878, checked in by Laurent Fairhead, 17 years ago

Bascule de la physique du LMD vers la physique avec thermiques
LF
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 10.5 KB
RevLine 
[524]1!
2! $Header$
3!
[878]4      SUBROUTINE ajsec(paprs, pplay, t,q,limbas,d_t,d_q)
[766]5      USE dimphy
[524]6      IMPLICIT none
7c======================================================================
8c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
9c Objet: ajustement sec (adaptation du GCM du LMD)
10c======================================================================
11c Arguments:
12c t-------input-R- Temperature
13c
14c d_t-----output-R-Incrementation de la temperature
15c======================================================================
[766]16cym#include "dimensions.h"
17cym#include "dimphy.h"
[524]18#include "YOMCST.h"
19      REAL paprs(klon,klev+1), pplay(klon,klev)
20      REAL t(klon,klev), q(klon,klev)
21      REAL d_t(klon,klev), d_q(klon,klev)
22c
[878]23      INTEGER limbas(klon), limhau ! les couches a ajuster
24c
25      LOGICAL mixq
26ccc      PARAMETER (mixq=.TRUE.)
27      PARAMETER (mixq=.FALSE.)
28c
29      REAL zh(klon,klev)
30      REAL zho(klon,klev)
31      REAL zq(klon,klev)
32      REAL zpk(klon,klev)
33      REAL zpkdp(klon,klev)
34      REAL hm, sm, qm
35      LOGICAL modif(klon), down
36      INTEGER i, k, k1, k2
37c
38c Initialisation:
39c
40cym
41      limhau=klev
42 
43      DO k = 1, klev
44      DO i = 1, klon
45         d_t(i,k) = 0.0
46         d_q(i,k) = 0.0
47      ENDDO
48      ENDDO
49c------------------------------------- detection des profils a modifier
50      DO k = 1, limhau
51      DO i = 1, klon
52         zpk(i,k) = pplay(i,k)**RKAPPA
53         zh(i,k) = RCPD * t(i,k)/ zpk(i,k)
54         zho(i,k) = zh(i,k)
55         zq(i,k) = q(i,k)
56      ENDDO
57      ENDDO
58c
59      DO k = 1, limhau
60      DO i = 1, klon
61         zpkdp(i,k) = zpk(i,k) * (paprs(i,k)-paprs(i,k+1))
62      ENDDO
63      ENDDO
64c
65      DO i = 1, klon
66         modif(i) = .FALSE.
67      ENDDO
68      DO k = 2, limhau
69      DO i = 1, klon
70      IF (.NOT.modif(i).and.k-1>limbas(i)) THEN
71         IF ( zh(i,k).LT.zh(i,k-1) ) modif(i) = .TRUE.
72      ENDIF
73      ENDDO
74      ENDDO
75c------------------------------------- correction des profils instables
76      DO 1080 i = 1, klon
77      IF (modif(i)) THEN
78          k2 = limbas(i)
79 8000     CONTINUE
80            k2 = k2 + 1
81            IF (k2 .GT. limhau) goto 8001
82            IF (zh(i,k2) .LT. zh(i,k2-1)) THEN
83              k1 = k2 - 1
84              k = k1
85              sm = zpkdp(i,k2)
86              hm = zh(i,k2)
87              qm = zq(i,k2)
88 8020         CONTINUE
89                sm = sm +zpkdp(i,k)
90                hm = hm +zpkdp(i,k) * (zh(i,k)-hm) / sm
91                qm = qm +zpkdp(i,k) * (zq(i,k)-qm) / sm
92                down = .FALSE.
93                IF (k1 .ne. limbas(i)) THEN
94                  IF (hm .LT. zh(i,k1-1)) down = .TRUE.
95                ENDIF
96                IF (down) THEN
97                  k1 = k1 - 1
98                  k = k1
99                ELSE
100                  IF ((k2 .EQ. limhau)) GOTO 8021
101                  IF ((zh(i,k2+1).GE.hm)) GOTO 8021
102                  k2 = k2 + 1
103                  k = k2
104                ENDIF
105              GOTO 8020
106 8021         CONTINUE
107c------------ nouveau profil : constant (valeur moyenne)
108              DO k = k1, k2
109                zh(i,k) = hm
110                zq(i,k) = qm
111              ENDDO
112              k2 = k2 + 1
113            ENDIF
114          GOTO 8000
115 8001     CONTINUE
116      ENDIF
117 1080 CONTINUE
118c
119      DO k = 1, limhau
120      DO i = 1, klon
121         d_t(i,k) = (zh(i,k)-zho(i,k))*zpk(i,k)/RCPD
122         d_q(i,k) = zq(i,k) - q(i,k)
123      ENDDO
124      ENDDO
125c
126! FH : les d_q et d_t sont maintenant calcules de facon a valoir
127! effectivement 0. si on ne fait rien.
128!
129!     IF (limbas.GT.1) THEN
130!     DO k = 1, limbas-1
131!     DO i = 1, klon
132!        d_t(i,k) = 0.0
133!        d_q(i,k) = 0.0
134!     ENDDO
135!     ENDDO
136!     ENDIF
137c
138!     IF (limhau.LT.klev) THEN
139!     DO k = limhau+1, klev
140!     DO i = 1, klon
141!        d_t(i,k) = 0.0
142!        d_q(i,k) = 0.0
143!     ENDDO
144!     ENDDO
145!     ENDIF
146c
147      IF (.NOT.mixq) THEN
148      DO k = 1, klev
149      DO i = 1, klon
150         d_q(i,k) = 0.0
151      ENDDO
152      ENDDO
153      ENDIF
154c
155      RETURN
156      END
157
158      SUBROUTINE ajsec_convV2(paprs, pplay, t,q, d_t,d_q)
159      USE dimphy
160      IMPLICIT none
161c======================================================================
162c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
163c Objet: ajustement sec (adaptation du GCM du LMD)
164c======================================================================
165c Arguments:
166c t-------input-R- Temperature
167c
168c d_t-----output-R-Incrementation de la temperature
169c======================================================================
170cym#include "dimensions.h"
171cym#include "dimphy.h"
172#include "YOMCST.h"
173      REAL paprs(klon,klev+1), pplay(klon,klev)
174      REAL t(klon,klev), q(klon,klev)
175      REAL d_t(klon,klev), d_q(klon,klev)
176c
[524]177      INTEGER limbas, limhau ! les couches a ajuster
178ccc      PARAMETER (limbas=klev-3, limhau=klev)
[766]179cym      PARAMETER (limbas=1, limhau=klev)
[524]180c
181      LOGICAL mixq
182ccc      PARAMETER (mixq=.TRUE.)
183      PARAMETER (mixq=.FALSE.)
184c
185      REAL zh(klon,klev)
186      REAL zq(klon,klev)
187      REAL zpk(klon,klev)
188      REAL zpkdp(klon,klev)
189      REAL hm, sm, qm
190      LOGICAL modif(klon), down
191      INTEGER i, k, k1, k2
192c
193c Initialisation:
194c
[766]195cym
196      limbas=1
197      limhau=klev
198 
[524]199      DO k = 1, klev
200      DO i = 1, klon
201         d_t(i,k) = 0.0
202         d_q(i,k) = 0.0
203      ENDDO
204      ENDDO
205c------------------------------------- detection des profils a modifier
206      DO k = limbas, limhau
207      DO i = 1, klon
208         zpk(i,k) = pplay(i,k)**RKAPPA
209         zh(i,k) = RCPD * t(i,k)/ zpk(i,k)
210         zq(i,k) = q(i,k)
211      ENDDO
212      ENDDO
213c
214      DO k = limbas, limhau
215      DO i = 1, klon
216         zpkdp(i,k) = zpk(i,k) * (paprs(i,k)-paprs(i,k+1))
217      ENDDO
218      ENDDO
219c
220      DO i = 1, klon
221         modif(i) = .FALSE.
222      ENDDO
223      DO k = limbas+1, limhau
224      DO i = 1, klon
225      IF (.NOT.modif(i)) THEN
226         IF ( zh(i,k).LT.zh(i,k-1) ) modif(i) = .TRUE.
227      ENDIF
228      ENDDO
229      ENDDO
230c------------------------------------- correction des profils instables
231      DO 1080 i = 1, klon
232      IF (modif(i)) THEN
233          k2 = limbas
234 8000     CONTINUE
235            k2 = k2 + 1
236            IF (k2 .GT. limhau) goto 8001
237            IF (zh(i,k2) .LT. zh(i,k2-1)) THEN
238              k1 = k2 - 1
239              k = k1
240              sm = zpkdp(i,k2)
241              hm = zh(i,k2)
242              qm = zq(i,k2)
243 8020         CONTINUE
244                sm = sm +zpkdp(i,k)
245                hm = hm +zpkdp(i,k) * (zh(i,k)-hm) / sm
246                qm = qm +zpkdp(i,k) * (zq(i,k)-qm) / sm
247                down = .FALSE.
248                IF (k1 .ne. limbas) THEN
249                  IF (hm .LT. zh(i,k1-1)) down = .TRUE.
250                ENDIF
251                IF (down) THEN
252                  k1 = k1 - 1
253                  k = k1
254                ELSE
255                  IF ((k2 .EQ. limhau)) GOTO 8021
256                  IF ((zh(i,k2+1).GE.hm)) GOTO 8021
257                  k2 = k2 + 1
258                  k = k2
259                ENDIF
260              GOTO 8020
261 8021         CONTINUE
262c------------ nouveau profil : constant (valeur moyenne)
263              DO k = k1, k2
264                zh(i,k) = hm
265                zq(i,k) = qm
266              ENDDO
267              k2 = k2 + 1
268            ENDIF
269          GOTO 8000
270 8001     CONTINUE
271      ENDIF
272 1080 CONTINUE
273c
274      DO k = limbas, limhau
275      DO i = 1, klon
276         d_t(i,k) = zh(i,k)*zpk(i,k)/RCPD - t(i,k)
277         d_q(i,k) = zq(i,k) - q(i,k)
278      ENDDO
279      ENDDO
280c
281      IF (limbas.GT.1) THEN
282      DO k = 1, limbas-1
283      DO i = 1, klon
284         d_t(i,k) = 0.0
285         d_q(i,k) = 0.0
286      ENDDO
287      ENDDO
288      ENDIF
289c
290      IF (limhau.LT.klev) THEN
291      DO k = limhau+1, klev
292      DO i = 1, klon
293         d_t(i,k) = 0.0
294         d_q(i,k) = 0.0
295      ENDDO
296      ENDDO
297      ENDIF
298c
299      IF (.NOT.mixq) THEN
300      DO k = 1, klev
301      DO i = 1, klon
302         d_q(i,k) = 0.0
303      ENDDO
304      ENDDO
305      ENDIF
306c
307      RETURN
308      END
309      SUBROUTINE ajsec_old(paprs, pplay, t, d_t)
[766]310      USE dimphy
[524]311      IMPLICIT none
312c======================================================================
313c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
314c Objet: ajustement sec (adaptation du GCM du LMD)
315c======================================================================
316c Arguments:
317c t-------input-R- Temperature
318c
319c d_t-----output-R-Incrementation de la temperature
320c======================================================================
[766]321cym#include "dimensions.h"
322cym#include "dimphy.h"
[524]323#include "YOMCST.h"
324      REAL paprs(klon,klev+1), pplay(klon,klev)
325      REAL t(klon,klev)
326      REAL d_t(klon,klev)
327c
328      REAL local_h(klon,klev)
329      REAL hm, sm
330      LOGICAL modif(klon), down
331      INTEGER i, l, l1, l2
332c------------------------------------- detection des profils a modifier
333      DO i = 1, klon
334         modif(i)   = .false.
335      ENDDO
336c
337      DO l = 1, klev
338      DO i = 1, klon
339         local_h(i,l) = RCPD * t(i,l)/ (pplay(i,l)**RKAPPA)
340      ENDDO
341      ENDDO
342c
343      DO l = 2, klev
344      DO i = 1, klon
345         IF ( local_h(i,l).lt.local_h(i,l-1) ) THEN
346            modif(i) = .true.
347         ELSE
348            modif(i) = modif(i)
349         ENDIF
350      ENDDO
351      ENDDO
352c------------------------------------- correction des profils instables
353      do 1080 i = 1, klon
354        if (modif(i)) then
355          l2 = 1
356 8000     continue
357            l2 = l2 + 1
358            if (l2 .gt. klev) goto 8001
359            if (local_h(i, l2) .lt. local_h(i, l2-1)) then
360              l1 = l2 - 1
361              l  = l1
362              sm = pplay(i,l2)**rkappa * (paprs(i,l2)-paprs(i,l2+1))
363              hm = local_h(i, l2)
364 8020         continue
365                sm = sm +pplay(i,l)**rkappa*(paprs(i,l)-paprs(i,l+1))
366                hm = hm +pplay(i,l)**rkappa*(paprs(i,l)-paprs(i,l+1))
367     .                         * (local_h(i, l) - hm) / sm
368                down = .false.
369                if (l1 .ne. 1) then
370                  if (hm .lt. local_h(i, l1-1)) then
371                    down = .true.
372                  end if
373                end if
374                if (down) then
375                  l1 = l1 - 1
376                  l  = l1
377                else
378                  if ((l2 .eq. klev)) GOTO 8021
379                  IF ((local_h(i, l2+1).ge.hm)) goto 8021
380                  l2 = l2 + 1
381                  l  = l2
382                end if
383              go to 8020
384 8021         continue
385c------------ nouveau profil : constant (valeur moyenne)
386              do 1100 l = l1, l2
387                local_h(i, l) = hm
388 1100         continue
389              l2 = l2 + 1
390            end if
391          go to 8000
392 8001     continue
393        end if
394 1080 continue
395c
396      DO l = 1, klev
397      DO i = 1, klon
398         d_t(i,l) = local_h(i,l)*(pplay(i,l)**rkappa)/RCPD - t(i,l)
399      ENDDO
400      ENDDO
401c
402      RETURN
403      END
Note: See TracBrowser for help on using the repository browser.