source: trunk/LMDZ.COMMON/libf/dyn3d/cpdet_mod.F90 @ 1086

Last change on this file since 1086 was 1086, checked in by slebonnois, 11 years ago

SL: new way to deal with cpdet routines in dyn3d*... This one should work.
File size: 9.5 KB
Line 
1module cpdet_mod
2
3implicit none
4
5! ADAPTATION OF GCM TO CP(T)
6!======================================================================
7! S. Lebonnois, 10/2010
8!
9! Cp must be computed using cpdet(t) to be valid
10!
11! The Exner function is still pk = RCPD*(play/pref)**RKAPPA
12! (RCPD=cpp, RKAPPA=kappa)
13!
14! One goes from T to teta (potential temperature) using t2tpot(t,teta,pk)
15! One goes from teta to T using tpot2t(teta,t,pk)
16!
17!======================================================================
18
19contains
20
21      SUBROUTINE ini_cpdet
22     
23      USE control_mod, ONLY: planet_type
24      IMPLICIT none
25!======================================================================
26! Initialization of nu_venus and t0_venus
27!======================================================================
28
29! for cpp, nu_venus and t0_venus:
30#include "comconst.h"
31
32      if (planet_type.eq."venus") then
33          nu_venus=0.35
34          t0_venus=460.
35      else
36          nu_venus=0.
37          t0_venus=0.
38      endif
39
40      return
41      end subroutine ini_cpdet
42
43!======================================================================
44!======================================================================
45
46      FUNCTION cpdet(t)
47
48      USE control_mod, ONLY: planet_type
49      IMPLICIT none
50
51! for cpp, nu_venus and t0_venus:
52#include "comconst.h"
53
54      real,intent(in) :: t
55      real cpdet
56
57      if (planet_type.eq."venus") then
58          cpdet = cpp*(t/t0_venus)**nu_venus
59      else
60          cpdet = cpp
61      endif
62
63      return
64      end function cpdet
65     
66!======================================================================
67!======================================================================
68
69      SUBROUTINE t2tpot(npoints, yt, yteta, ypk)
70!======================================================================
71! Arguments:
72!
73! yt   --------input-R- Temperature
74! yteta-------output-R- Temperature potentielle
75! ypk  --------input-R- Fonction d'Exner: RCPD*(pplay/pref)**RKAPPA
76!
77!======================================================================
78
79      USE control_mod, ONLY: planet_type
80      IMPLICIT NONE
81     
82! for cpp, nu_venus and t0_venus:
83#include "comconst.h"
84
85      integer,intent(in) :: npoints
86      REAL,intent(in) :: yt(npoints), ypk(npoints)
87      REAL,intent(out) :: yteta(npoints)
88     
89      if (planet_type.eq."venus") then
90          yteta = yt**nu_venus                                          &
91     &            - nu_venus * t0_venus**nu_venus * log(ypk/cpp)
92          yteta = yteta**(1./nu_venus)
93      else
94          yteta = yt * cpp/ypk
95      endif
96
97      return
98      end subroutine t2tpot
99
100!======================================================================
101!======================================================================
102
103      SUBROUTINE tpot2t(npoints,yteta, yt, ypk)
104!======================================================================
105! Arguments:
106!
107! yteta--------input-R- Temperature potentielle
108! yt   -------output-R- Temperature
109! ypk  --------input-R- Fonction d'Exner: RCPD*(pplay/pref)**RKAPPA
110!
111!======================================================================
112
113      USE control_mod, ONLY: planet_type
114      IMPLICIT NONE
115
116! for cpp, nu_venus and t0_venus:
117#include "comconst.h"
118
119      integer,intent(in) :: npoints
120      REAL,intent(in) :: yteta(npoints), ypk(npoints)
121      REAL,intent(out) :: yt(npoints)
122     
123      if (planet_type.eq."venus") then
124         yt = yteta**nu_venus                                           &
125     &       + nu_venus * t0_venus**nu_venus * log(ypk/cpp)
126         yt = yt**(1./nu_venus)
127      else
128          yt = yteta * ypk/cpp
129      endif
130 
131      return
132      end subroutine tpot2t
133
134!======================================================================
135!======================================================================
136! Routines pour les calculs paralleles
137!======================================================================
138#ifdef CPP_PARA
139!======================================================================
140
141      SUBROUTINE t2tpot_p(nlon,nlev, yt, yteta, ypk)
142! Parallel version of t2tpot, for an arbitrary number of columns
143      USE control_mod, only : planet_type
144      IMPLICIT none
145
146! for cpp, nu_venus and t0_venus:
147#include "comconst.h"
148
149      integer,intent(in) :: nlon,nlev
150      real,intent(in) :: yt(nlon,nlev)
151      real,intent(out) :: yteta(nlon,nlev)
152      real,intent(in) :: ypk(nlon,nlev)
153! local variable:
154      integer :: l
155     
156      if (planet_type.eq."venus") then
157!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
158        do l=1,nlev
159          yteta(:,l)=yt(:,l)**nu_venus                                  &
160     &                     -nu_venus*t0_venus**nu_venus*                &
161     &                          log(ypk(:,l)/cpp)
162          yteta(:,l)=yteta(:,l)**(1./nu_venus)
163        enddo
164!$OMP END DO
165      else
166!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
167        do l=1,nlev
168          yteta(:,l)=yt(:,l)*cpp/ypk(:,l)
169        enddo
170!$OMP END DO
171      endif ! of if (planet_type.eq."venus")
172
173      end subroutine t2tpot_p
174
175!======================================================================
176!======================================================================
177
178      SUBROUTINE t2tpot_glo_p(yt, yteta, ypk)
179! Parallel version of t2tpot, over the full dynamics (scalar) grid
180! (more efficient than multiple calls to t2tpot_p() with slices of data)
181      USE parallel_lmdz, only : jj_begin,jj_end
182      USE control_mod, only : planet_type
183      IMPLICIT none
184! for iip1, jjp1 and llm
185#include "dimensions.h"
186#include "paramet.h"
187! for cpp, nu_venus and t0_venus:
188#include "comconst.h"
189
190      real,intent(in) :: yt(iip1,jjp1,llm)
191      real,intent(out) :: yteta(iip1,jjp1,llm)
192      real,intent(in) :: ypk(iip1,jjp1,llm)
193! local variable:
194      integer :: j,l
195      integer :: jjb,jje
196     
197      jjb=jj_begin
198      jje=jj_end
199
200      if (planet_type.eq."venus") then
201!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
202        do l=1,llm
203          yteta(:,jjb:jje,l)=yt(:,jjb:jje,l)**nu_venus                  &
204     &                     -nu_venus*t0_venus**nu_venus*                &
205     &                          log(ypk(:,jjb:jje,l)/cpp)
206          yteta(:,jjb:jje,l)=yteta(:,jjb:jje,l)**(1./nu_venus)
207        enddo
208!$OMP END DO
209      else
210!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
211        do l=1,llm
212          yteta(:,jjb:jje,l)=yt(:,jjb:jje,l)*cpp/ypk(:,jjb:jje,l)
213        enddo
214!$OMP END DO
215      endif ! of if (planet_type.eq."venus")
216
217      end subroutine t2tpot_glo_p
218
219!======================================================================
220!======================================================================
221
222      SUBROUTINE tpot2t_p(nlon,nlev,yteta,yt,ypk)
223! Parallel version of tpot2t, for an arbitrary number of columns
224      USE control_mod, only : planet_type
225      IMPLICIT none
226! for cpp, nu_venus and t0_venus:
227#include "comconst.h"
228
229      integer,intent(in) :: nlon,nlev
230      real,intent(out) :: yt(nlon,nlev)
231      real,intent(in) :: yteta(nlon,nlev)
232      real,intent(in) :: ypk(nlon,nlev)
233
234! local variable:
235      integer :: l
236
237      if (planet_type.eq."venus") then
238!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
239        do l=1,nlev
240          yt(:,l)=yteta(:,l)**nu_venus                                  &
241     &                  +nu_venus*t0_venus**nu_venus*                   &
242     &                       log(ypk(:,l)/cpp)
243          yt(:,l)=yt(:,l)**(1./nu_venus)
244        enddo
245!$OMP END DO
246      else
247!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
248        do l=1,nlev
249          yt(:,l)=yteta(:,l)*ypk(:,l)/cpp
250        enddo
251!$OMP END DO
252      endif ! of if (planet_type.eq."venus")
253      end subroutine tpot2t_p
254
255!======================================================================
256!======================================================================
257
258      SUBROUTINE tpot2t_glo_p(yteta,yt,ypk)
259! Parallel version of tpot2t, over the full dynamics (scalar) grid
260! (more efficient than multiple calls to tpot2t_p() with slices of data)
261      USE parallel_lmdz, only : jj_begin,jj_end
262      USE control_mod, only : planet_type
263      IMPLICIT none
264! for iip1, jjp1 and llm
265#include "dimensions.h"
266#include "paramet.h"
267! for cpp, nu_venus and t0_venus:
268#include "comconst.h"
269
270      real,intent(out) :: yt(iip1,jjp1,llm)
271      real,intent(in) :: yteta(iip1,jjp1,llm)
272      real,intent(in) :: ypk(iip1,jjp1,llm)
273! local variable:
274      integer :: j,l
275      integer :: jjb,jje
276     
277      jjb=jj_begin
278      jje=jj_end
279
280      if (planet_type.eq."venus") then
281!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
282        do l=1,llm
283          yt(:,jjb:jje,l)=yteta(:,jjb:jje,l)**nu_venus                  &
284     &                  +nu_venus*t0_venus**nu_venus*                   &
285     &                       log(ypk(:,jjb:jje,l)/cpp)
286          yt(:,jjb:jje,l)=yt(:,jjb:jje,l)**(1./nu_venus)
287        enddo
288!$OMP END DO
289      else
290!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
291        do l=1,llm
292          yt(:,jjb:jje,l)=yteta(:,jjb:jje,l)*ypk(:,jjb:jje,l)/cpp
293        enddo
294!$OMP END DO
295      endif ! of if (planet_type.eq."venus")
296      end subroutine tpot2t_glo_p
297
298!======================================================================
299#endif
300!======================================================================
301! Fin routines specifiques parallele
302!======================================================================
303!======================================================================
304!
305! ATTENTION
306!
307! Si un jour on a besoin, il faudra coder les routines
308!    dt2dtpot / dtpto2dt
309!
310!======================================================================
311!======================================================================
312end module cpdet_mod
Note: See TracBrowser for help on using the repository browser.