source: trunk/libf/dyn3dpar/cpdet.F @ 37

Last change on this file since 37 was 37, checked in by emillour, 14 years ago

Remise en route chantier compilation -- Ehouarn

  • Modifs et corrections pour pouvoir compiler le gcm (en séentiel, avec

makelmdz_fcm pour l'instant):

  • ajout de fichiers 'arch' pour linux-64 (pour Bellonzi, avec ioipsl et en r8)
  • modification de makelmdz_fcm, ajout de la cléPP_PHYS si on compile avec une physique
  • correction de quelques typos/bugs réléàa compilation:
  • infotrac.F90 : supression des appels àlnblnk' (remplacépar len_trim)
  • bilan_dyn.F : déaration des variables znom3,znom3l,zunites3, planet_type
  • cpdet.F : "use control_mod, ONLY: planet_type" mis aux bons endroits

(idem sur cpdet.F dans dyn3dpar)

  • leapfrog.F : declaration de ztetaec(), dtec, cpdet , itau_w, duspg()
  • diagedyn.F : correction typo; attention dans diagedyn.F il y a du

include "../phylmd/YOMCST.h"
Ca va poser problè qd on change de physique....

  • Avec ces modifs, la compilation marche sans physque, avec ou sans ioipsl et avec la physique terrestre phylmd.
File size: 6.1 KB
Line 
1! ADAPTATION GCM POUR CP(T)
2c======================================================================
3c S. Lebonnois, 10/2010
4c
5c Cp doit être calculé par cpdet(t) pour être valable partout
6c
7c La fonction d'Exner reste pk = RCPD*(play/pref)**RKAPPA
8c (RCPD=cpp, RKAPPA=kappa)
9c
10c On passe de T a teta (temperature potentielle) par t2tpot(t,teta,pk)
11c On passe de teta a T par tpot2t(teta,t,pk)
12c
13c======================================================================
14
15      SUBROUTINE ini_cpdet
16     
17      USE control_mod, ONLY: planet_type
18      IMPLICIT none
19c======================================================================
20c Initialisation de nu_venus et t0_venus
21c======================================================================
22
23! for cpp, nu_venus and t0_venus:
24#include "comconst.h"
25
26      if (planet_type.eq."venus") then
27          nu_venus=0.35
28          t0_venus=460.
29      else
30          nu_venus=0.
31          t0_venus=0.
32      endif
33
34      return
35      end
36
37c======================================================================
38c======================================================================
39
40      FUNCTION cpdet(t)
41
42      USE control_mod, ONLY: planet_type
43      IMPLICIT none
44
45! for cpp, nu_venus and t0_venus:
46#include "comconst.h"
47
48      real cpdet,t
49
50      if (planet_type.eq."venus") then
51          cpdet = cpp*(t/t0_venus)**nu_venus
52      else
53          cpdet = cpp
54      endif
55
56      return
57      end
58     
59c======================================================================
60c======================================================================
61
62      SUBROUTINE t2tpot(npoints, yt, yteta, ypk)
63c======================================================================
64c Arguments:
65c
66c yt   --------input-R- Temperature
67c yteta-------output-R- Temperature potentielle
68c ypk  --------input-R- Fonction d'Exner: RCPD*(pplay/pref)**RKAPPA
69c
70c======================================================================
71
72      USE control_mod, ONLY: planet_type
73      IMPLICIT NONE
74     
75! for cpp, nu_venus and t0_venus:
76#include "comconst.h"
77
78      integer npoints
79      REAL    yt(npoints), yteta(npoints), ypk(npoints)
80     
81      if (planet_type.eq."venus") then
82          yteta = yt**nu_venus                                          &
83     &            - nu_venus * t0_venus**nu_venus * log(ypk/cpp)
84          yteta = yteta**(1./nu_venus)
85      else
86          yteta = yt * cpp/ypk
87      endif
88
89      return
90      end
91
92c======================================================================
93c======================================================================
94
95      SUBROUTINE t2tpot_p(ip1jmp1,llm, yt, yteta, ypk)
96! Parallel version of t2tpot
97      USE parallel
98      USE control_mod, only : planet_type
99      IMPLICIT none
100
101! for cpp, nu_venus and t0_venus:
102#include "comconst.h"
103
104      integer,intent(in) :: ip1jmp1,llm
105      real,intent(in) :: yt(ip1jmp1,llm)
106      real,intent(out) :: yteta(ip1jmp1,llm)
107      real,intent(in) :: ypk(ip1jmp1,llm)
108! local variable:
109      integer :: ij,l,ijb,ije
110     
111      ijb=ij_begin
112      ije=ij_end 
113     
114      if (planet_type.eq."venus") then
115!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
116        do l=1,llm
117          yteta(ijb:ije,l)=yt(ijb:ije,l)**nu_venus                      &
118     &                     -nu_venus*t0_venus**nu_venus*                &
119     &                          log(ypk(ijb:ije,l)/cpp)
120          yteta(ijb:ije,l)=yteta(ijb:ije,l)**(1./nu_venus)
121        enddo
122!$OMP END DO
123      else
124!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
125        do l=1,llm
126          yteta(ijb:ije,l)=yt(ijb:ije,l)*cpp/ypk(ijb:ije,l)
127        enddo
128!$OMP END DO
129      endif ! of if (planet_type.eq."venus")
130
131      END
132
133c======================================================================
134
135      SUBROUTINE tpot2t(npoints,yteta, yt, ypk)
136c======================================================================
137c Arguments:
138c
139c yteta--------input-R- Temperature potentielle
140c yt   -------output-R- Temperature
141c ypk  --------input-R- Fonction d'Exner: RCPD*(pplay/pref)**RKAPPA
142c
143c======================================================================
144
145      USE control_mod, ONLY: planet_type
146      IMPLICIT NONE
147
148! for cpp, nu_venus and t0_venus:
149#include "comconst.h"
150
151      integer npoints
152      REAL    yt(npoints), yteta(npoints), ypk(npoints)
153     
154      if (planet_type.eq."venus") then
155         yt = yteta**nu_venus                                           &
156     &       + nu_venus * t0_venus**nu_venus * log(ypk/cpp)
157         yt = yt**(1./nu_venus)
158      else
159          yt = yteta * ypk/cpp
160      endif
161 
162      return
163      end
164
165c======================================================================
166c======================================================================
167      SUBROUTINE tpot2t_p(ip1jmp1,llm,yteta, yt, ypk)
168! Parallel version of tpot2t
169      USE parallel
170      USE control_mod, only : planet_type
171      IMPLICIT none
172! for cpp, nu_venus and t0_venus:
173#include "comconst.h"
174
175      integer,intent(in) :: ip1jmp1,llm
176      real,intent(out) :: yt(ip1jmp1,llm)
177      real,intent(in) :: yteta(ip1jmp1,llm)
178      real,intent(in) :: ypk(ip1jmp1,llm)
179! local variable:
180      integer :: ij,l,ijb,ije
181     
182      ijb=ij_begin
183      ije=ij_end 
184
185      if (planet_type.eq."venus") then
186!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
187        do l=1,llm
188          yt(ijb:ije,l)=yteta(ijb:ije,l)**nu_venus                      &
189     &                  +nu_venus*t0_venus**nu_venus*                   &
190     &                       log(ypk(ijb:ije,l)/cpp)
191          yt(ijb:ije,l)=yt(ijb:ije,l)**(1./nu_venus)
192        enddo
193!$OMP END DO
194      else
195!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
196        do l=1,llm
197          yt(ijb:ije,l)=yteta(ijb:ije,l)*ypk(ijb:ije,l)/cpp
198        enddo
199!$OMP END DO
200      endif ! of if (planet_type.eq."venus")
201      END
202
203c======================================================================
204c======================================================================
205c
206c ATTENTION
207c
208c Si un jour on a besoin, il faudra coder les routines
209c    dt2dtpot / dtpto2dt
210c
211c======================================================================
212c======================================================================
Note: See TracBrowser for help on using the repository browser.