/[lmdze]/trunk/Sources/dyn3d/comgeom.f
ViewVC logotype

Diff of /trunk/Sources/dyn3d/comgeom.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/dyn3d/comgeom.f90 revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC trunk/Sources/dyn3d/comgeom.f revision 161 by guez, Fri Jul 24 14:27:59 2015 UTC
# Line 1  Line 1 
1  module comgeom  module comgeom
2    
3    use dimens_m, only: iim, jjm    use dimens_m, only: iim, jjm
   use paramet_m, only: ip1jmp1, ip1jm  
4    
5    implicit none    implicit none
6    
7    private iim, jjm, ip1jmp1, ip1jm    private iim, jjm
8    
9    real cu_2d(iim + 1, jjm + 1), cv_2d(iim + 1, jjm)    real cu_2d(iim + 1, jjm + 1), cv_2d(iim + 1, jjm) ! in m
10    real cu(ip1jmp1), cv(ip1jm)    real cu((iim + 1) * (jjm + 1)), cv((iim + 1) * jjm) ! in m
11    equivalence (cu, cu_2d), (cv, cv_2d)    equivalence (cu, cu_2d), (cv, cv_2d)
12    
13    real unscu2_2d(iim + 1, jjm + 1)    real unscu2_2d(iim + 1, jjm + 1) ! in m-2
14    real unscu2(ip1jmp1)    real unscu2((iim + 1) * (jjm + 1)) ! in m-2
15    equivalence (unscu2, unscu2_2d)    equivalence (unscu2, unscu2_2d)
16    
17    real unscv2_2d(iim + 1,jjm)    real unscv2_2d(iim + 1, jjm) ! in m-2
18    real unscv2(ip1jm)    real unscv2((iim + 1) * jjm) ! in m-2
19    equivalence (unscv2, unscv2_2d)    equivalence (unscv2, unscv2_2d)
20    
21    real aire_2d(iim + 1,jjm + 1), airesurg_2d(iim + 1,jjm + 1)    real aire((iim + 1) * (jjm + 1)), aire_2d(iim + 1, jjm + 1) ! in m2
22    real aire(ip1jmp1), airesurg(ip1jmp1)    real airesurg_2d(iim + 1, jjm + 1), airesurg((iim + 1) * (jjm + 1))
23    equivalence (aire, aire_2d), (airesurg, airesurg_2d)    equivalence (aire, aire_2d), (airesurg, airesurg_2d)
24    
25    real aireu_2d(iim + 1,jjm + 1)              real aireu_2d(iim + 1, jjm + 1) ! in m2
26    real aireu(ip1jmp1)    real aireu((iim + 1) * (jjm + 1)) ! in m2
27    equivalence (aireu, aireu_2d)    equivalence (aireu, aireu_2d)
28    
29    real airev_2d(iim + 1,jjm),unsaire_2d(iim + 1,jjm + 1)    real airev((iim + 1) * jjm), airev_2d(iim + 1, jjm) ! in m2
30    real airev(ip1jm),unsaire(ip1jmp1)    real unsaire((iim + 1) * (jjm + 1)), unsaire_2d(iim + 1, jjm + 1) ! in m-2
31    equivalence (airev, airev_2d), (unsaire, unsaire_2d)    equivalence (airev, airev_2d), (unsaire, unsaire_2d)
32    
33    real apoln,apols    real apoln, apols ! in m2
34    
35    real unsairez_2d(iim + 1,jjm),airuscv2_2d(iim + 1,jjm)    real unsairez_2d(iim + 1, jjm)
36    real unsairez(ip1jm),airuscv2(ip1jm)    real unsairez((iim + 1) * jjm)
37    equivalence (unsairez, unsairez_2d), (airuscv2, airuscv2_2d)    equivalence (unsairez, unsairez_2d)
38    
39    real airvscu2_2d(iim + 1,jjm)          real alpha1_2d(iim + 1, jjm + 1)
40    real airvscu2(ip1jm)    real alpha1((iim + 1) * (jjm + 1))
41    equivalence (airvscu2, airvscu2_2d)    equivalence (alpha1, alpha1_2d)
   
   real aireij1_2d(iim + 1,jjm + 1),aireij2_2d(iim + 1,jjm + 1)  
   real aireij1(ip1jmp1),aireij2(ip1jmp1)  
   equivalence (aireij1, aireij1_2d), (aireij2, aireij2_2d)  
   
   real aireij3(ip1jmp1)  
   real aireij3_2d(iim + 1,jjm + 1)        
   equivalence (aireij3, aireij3_2d)  
   
   real aireij4_2d(iim + 1,jjm + 1), alpha1_2d(iim + 1,jjm + 1)  
   real aireij4(ip1jmp1), alpha1(ip1jmp1)  
   equivalence (aireij4, aireij4_2d), (alpha1, alpha1_2d)  
42    
43    real alpha2_2d(iim + 1,jjm + 1)            real alpha2_2d(iim + 1, jjm + 1)        
44    real alpha2(ip1jmp1)    real alpha2((iim + 1) * (jjm + 1))
45    equivalence (alpha2, alpha2_2d)    equivalence (alpha2, alpha2_2d)
46    
47    real alpha3_2d(iim + 1,jjm + 1), alpha4_2d(iim + 1,jjm + 1)    real alpha3_2d(iim + 1, jjm + 1), alpha4_2d(iim + 1, jjm + 1)
48    real alpha3(ip1jmp1), alpha4(ip1jmp1)    real alpha3((iim + 1) * (jjm + 1)), alpha4((iim + 1) * (jjm + 1))
49    equivalence (alpha3, alpha3_2d), (alpha4, alpha4_2d)    equivalence (alpha3, alpha3_2d), (alpha4, alpha4_2d)
50    
51    real alpha1p2_2d(iim + 1,jjm + 1)            real alpha1p2_2d(iim + 1, jjm + 1)        
52    real alpha1p2(ip1jmp1)    real alpha1p2((iim + 1) * (jjm + 1))
53    equivalence (alpha1p2, alpha1p2_2d)    equivalence (alpha1p2, alpha1p2_2d)
54    
55    real alpha1p4_2d(iim + 1,jjm + 1),alpha2p3_2d(iim + 1,jjm + 1)    real alpha1p4_2d(iim + 1, jjm + 1), alpha2p3_2d(iim + 1, jjm + 1)
56    real alpha1p4(ip1jmp1),alpha2p3(ip1jmp1)    real alpha1p4((iim + 1) * (jjm + 1)), alpha2p3((iim + 1) * (jjm + 1))
57    equivalence (alpha1p4, alpha1p4_2d), (alpha2p3, alpha2p3_2d)    equivalence (alpha1p4, alpha1p4_2d), (alpha2p3, alpha2p3_2d)
58    
59    real alpha3p4(ip1jmp1)    real alpha3p4((iim + 1) * (jjm + 1))
60    real alpha3p4_2d(iim + 1,jjm + 1)        real alpha3p4_2d(iim + 1, jjm + 1)    
61    equivalence (alpha3p4, alpha3p4_2d)    equivalence (alpha3p4, alpha3p4_2d)
62    
63    real fext_2d(iim + 1,jjm),constang_2d(iim + 1,jjm + 1)    real fext_2d(iim + 1, jjm), constang_2d(iim + 1, jjm + 1)
64    real fext(ip1jm),constang(ip1jmp1)    real fext((iim + 1) * jjm), constang((iim + 1) * (jjm + 1))
65    equivalence (fext, fext_2d), (constang, constang_2d)    equivalence (fext, fext_2d), (constang, constang_2d)
66    
67    real rlatu(jjm + 1)    real cuvsurcv_2d(iim + 1, jjm), cvsurcuv_2d(iim + 1, jjm) ! no dimension
68    ! (latitudes of points of the "scalar" and "u" grid, in rad)    real cuvsurcv((iim + 1) * jjm), cvsurcuv((iim + 1) * jjm) ! no dimension
   
   real rlatv(jjm)  
   ! (latitudes of points of the "v" grid, in rad, in decreasing order)  
   
   real rlonu(iim + 1) ! longitudes of points of the "u" grid, in rad  
   
   real rlonv(iim + 1)  
   ! (longitudes of points of the "scalar" and "v" grid, in rad)  
   
   real cuvsurcv_2d(iim + 1,jjm),cvsurcuv_2d(iim + 1,jjm)    
   real cuvsurcv(ip1jm),cvsurcuv(ip1jm)  
69    equivalence (cuvsurcv, cuvsurcv_2d), (cvsurcuv, cvsurcuv_2d)    equivalence (cuvsurcv, cuvsurcv_2d), (cvsurcuv, cvsurcuv_2d)
70    
71    real cvusurcu_2d(iim + 1,jjm + 1),cusurcvu_2d(iim + 1,jjm + 1)    real cvusurcu_2d(iim + 1, jjm + 1), cusurcvu_2d(iim + 1, jjm + 1)
72    real cvusurcu(ip1jmp1),cusurcvu(ip1jmp1)    ! no dimension
73      real cvusurcu((iim + 1) * (jjm + 1)), cusurcvu((iim + 1) * (jjm + 1))
74      ! no dimension
75    equivalence (cvusurcu, cvusurcu_2d), (cusurcvu, cusurcvu_2d)    equivalence (cvusurcu, cvusurcu_2d), (cusurcvu, cusurcvu_2d)
76    
77    real cuvscvgam1_2d(iim + 1,jjm)    real cuvscvgam1_2d(iim + 1, jjm)
78    real cuvscvgam1(ip1jm)    real cuvscvgam1((iim + 1) * jjm)
79    equivalence (cuvscvgam1, cuvscvgam1_2d)    equivalence (cuvscvgam1, cuvscvgam1_2d)
80    
81    real cuvscvgam2_2d(iim + 1,jjm),cvuscugam1_2d(iim + 1,jjm + 1)    real cuvscvgam2_2d(iim + 1, jjm), cvuscugam1_2d(iim + 1, jjm + 1)
82    real cuvscvgam2(ip1jm),cvuscugam1(ip1jmp1)    real cuvscvgam2((iim + 1) * jjm), cvuscugam1((iim + 1) * (jjm + 1))
83    equivalence (cuvscvgam2, cuvscvgam2_2d), (cvuscugam1, cvuscugam1_2d)    equivalence (cuvscvgam2, cuvscvgam2_2d), (cvuscugam1, cvuscugam1_2d)
84    
85    real cvuscugam2_2d(iim + 1,jjm + 1),cvscuvgam_2d(iim + 1,jjm)    real cvuscugam2_2d(iim + 1, jjm + 1), cvscuvgam_2d(iim + 1, jjm)
86    real cvuscugam2(ip1jmp1),cvscuvgam(ip1jm)    real cvuscugam2((iim + 1) * (jjm + 1)), cvscuvgam((iim + 1) * jjm)
87    equivalence (cvuscugam2, cvuscugam2_2d), (cvscuvgam, cvscuvgam_2d)    equivalence (cvuscugam2, cvuscugam2_2d), (cvscuvgam, cvscuvgam_2d)
88    
89    real cuscvugam(ip1jmp1)    real cuscvugam((iim + 1) * (jjm + 1))
90    real cuscvugam_2d(iim + 1,jjm + 1)    real cuscvugam_2d(iim + 1, jjm + 1)
91    equivalence (cuscvugam, cuscvugam_2d)    equivalence (cuscvugam, cuscvugam_2d)
92    
93    real unsapolnga1,unsapolnga2,unsapolsga1,unsapolsga2                    real unsapolnga1, unsapolnga2, unsapolsga1, unsapolsga2                
94    
95    real unsair_gam1_2d(iim + 1,jjm + 1),unsair_gam2_2d(iim + 1,jjm + 1)    real unsair_gam1_2d(iim + 1, jjm + 1), unsair_gam2_2d(iim + 1, jjm + 1)
96    real unsair_gam1(ip1jmp1),unsair_gam2(ip1jmp1)    real unsair_gam1((iim + 1) * (jjm + 1)), unsair_gam2((iim + 1) * (jjm + 1))
97    equivalence (unsair_gam1, unsair_gam1_2d), (unsair_gam2, unsair_gam2_2d)    equivalence (unsair_gam1, unsair_gam1_2d), (unsair_gam2, unsair_gam2_2d)
98    
99    real unsairz_gam_2d(iim + 1,jjm)    real unsairz_gam_2d(iim + 1, jjm)
100    real unsairz_gam(ip1jm)    real unsairz_gam((iim + 1) * jjm)
101    equivalence (unsairz_gam, unsairz_gam_2d)    equivalence (unsairz_gam, unsairz_gam_2d)
102    
103    real aivscu2gam_2d(iim + 1,jjm),aiuscv2gam_2d(iim + 1,jjm)    save
104    real aivscu2gam(ip1jm),aiuscv2gam(ip1jm)  
105    equivalence (aivscu2gam, aivscu2gam_2d), (aiuscv2gam, aiuscv2gam_2d)  contains
106    
107    real xprimu(iim + 1),xprimv(iim + 1)    SUBROUTINE inigeom
108    
109    save      ! Auteur : P. Le Van
110    
111        ! Calcul des élongations cuij1, ..., cuij4, cvij1, ..., cvij4 aux mêmes
112        ! endroits que les aires aireij1_2d, ..., aireij4_2d.
113    
114        ! Calcul des coefficients cu_2d, cv_2d, 1. / cu_2d**2, 1. /
115        ! cv_2d**2. Les coefficients cu_2d et cv_2d permettent de passer
116        ! des vitesses naturelles aux vitesses covariantes et
117        ! contravariantes, ou vice-versa.
118    
119        ! On a :
120        ! u(covariant) = cu_2d * u(naturel), u(contravariant) = u(naturel) / cu_2d
121        ! v(covariant) = cv_2d * v(naturel), v(contravariant) = v(naturel) / cv_2d
122    
123        ! On en tire :
124        ! u(covariant) = cu_2d * cu_2d * u(contravariant)
125        ! v(covariant) = cv_2d * cv_2d * v(contravariant)
126    
127        ! x est la longitude du point en radians.
128        ! y est la latitude du point en radians.
129        !
130        ! On a : cu_2d(i, j) = rad * cos(y) * dx / dX
131        ! cv(j) = rad * dy / dY
132        ! aire_2d(i, j) = cu_2d(i, j) * cv(j)
133        !
134        ! y, dx / dX, dy / dY calculés aux points concernés. cv, bien que
135        ! dépendant de j uniquement, sera ici indicé aussi en i pour un
136        ! adressage plus facile en ij.
137    
138        ! cv_2d est aux points v. cu_2d est aux points u. Cf. "inigeom.txt".
139    
140        USE comconst, ONLY : g, omeg, rad
141        USE comdissnew, ONLY : coefdis, nitergdiv, nitergrot, niterh
142        use dynetat0_m, only: xprimp025, xprimm025, rlatu1, rlatu2, rlatu, rlatv, &
143             yprimu1, yprimu2
144        use nr_util, only: pi
145        USE paramet_m, ONLY : iip1, jjp1
146    
147        ! Local:
148        INTEGER i, j
149        REAL ai14, ai23, airez, un4rad2
150        REAL coslatm, coslatp, radclatm, radclatp
151        REAL, dimension(iip1, jjp1):: cuij1, cuij2, cuij3, cuij4 ! in m
152        REAL, dimension(iip1, jjp1):: cvij1, cvij2, cvij3, cvij4 ! in m
153        REAL gamdi_gdiv, gamdi_grot, gamdi_h
154        real, dimension(iim + 1, jjm + 1):: aireij1_2d, aireij2_2d, aireij3_2d, &
155             aireij4_2d ! in m2
156    
157        !------------------------------------------------------------------
158    
159        PRINT *, 'Call sequence information: inigeom'
160    
161        IF (nitergdiv /= 2) THEN
162           gamdi_gdiv = coefdis / (nitergdiv - 2)
163        ELSE
164           gamdi_gdiv = 0.
165        END IF
166    
167        IF (nitergrot /= 2) THEN
168           gamdi_grot = coefdis / (nitergrot - 2)
169        ELSE
170           gamdi_grot = 0.
171        END IF
172    
173        IF (niterh /= 2) THEN
174           gamdi_h = coefdis / (niterh - 2)
175        ELSE
176           gamdi_h = 0.
177        END IF
178    
179        print *, 'gamdi_gdiv = ', gamdi_gdiv
180        print *, "gamdi_grot = ", gamdi_grot
181        print *, "gamdi_h = ", gamdi_h
182    
183        un4rad2 = 0.25 * rad * rad
184    
185        ! Cf. "inigeom.txt". Calcul des quatre aires élémentaires
186        ! aireij1_2d, aireij2_2d, aireij3_2d, aireij4_2d qui entourent
187        ! chaque aire_2d(i, j), ainsi que les quatre élongations
188        ! élémentaires cuij et les quatre élongations cvij qui sont
189        ! calculées aux mêmes endroits que les aireij.
190    
191        coslatm = cos(rlatu1(1))
192        radclatm = 0.5 * rad * coslatm
193    
194        aireij1_2d(:iim, 1) = 0.
195        aireij2_2d(:iim, 1) = un4rad2 * coslatm * xprimp025(:iim) * yprimu1(1)
196        aireij3_2d(:iim, 1) = un4rad2 * coslatm * xprimm025(:iim) * yprimu1(1)
197        aireij4_2d(:iim, 1) = 0.
198    
199        cuij1(:iim, 1) = 0.
200        cuij2(:iim, 1) = radclatm * xprimp025(:iim)
201        cuij3(:iim, 1) = radclatm * xprimm025(:iim)
202        cuij4(:iim, 1) = 0.
203    
204        cvij1(:iim, 1) = 0.
205        cvij2(:iim, 1) = 0.5 * rad * yprimu1(1)
206        cvij3(:iim, 1) = cvij2(:iim, 1)
207        cvij4(:iim, 1) = 0.
208    
209        do j = 2, jjm
210           coslatm = cos(rlatu1(j))
211           coslatp = cos(rlatu2(j-1))
212           radclatp = 0.5 * rad * coslatp
213           radclatm = 0.5 * rad * coslatm
214           ai14 = un4rad2 * coslatp * yprimu2(j-1)
215           ai23 = un4rad2 * coslatm * yprimu1(j)
216    
217           aireij1_2d(:iim, j) = ai14 * xprimp025(:iim)
218           aireij2_2d(:iim, j) = ai23 * xprimp025(:iim)
219           aireij3_2d(:iim, j) = ai23 * xprimm025(:iim)
220           aireij4_2d(:iim, j) = ai14 * xprimm025(:iim)
221           cuij1(:iim, j) = radclatp * xprimp025(:iim)
222           cuij2(:iim, j) = radclatm * xprimp025(:iim)
223           cuij3(:iim, j) = radclatm * xprimm025(:iim)
224           cuij4(:iim, j) = radclatp * xprimm025(:iim)
225           cvij1(:iim, j) = 0.5 * rad * yprimu2(j-1)
226           cvij2(:iim, j) = 0.5 * rad * yprimu1(j)
227           cvij3(:iim, j) = cvij2(:iim, j)
228           cvij4(:iim, j) = cvij1(:iim, j)
229        end do
230    
231        coslatp = cos(rlatu2(jjm))
232        radclatp = 0.5 * rad * coslatp
233    
234        aireij1_2d(:iim, jjp1) = un4rad2 * coslatp * xprimp025(:iim) * yprimu2(jjm)
235        aireij2_2d(:iim, jjp1) = 0.
236        aireij3_2d(:iim, jjp1) = 0.
237        aireij4_2d(:iim, jjp1) = un4rad2 * coslatp * xprimm025(:iim) * yprimu2(jjm)
238    
239        cuij1(:iim, jjp1) = radclatp * xprimp025(:iim)
240        cuij2(:iim, jjp1) = 0.
241        cuij3(:iim, jjp1) = 0.
242        cuij4(:iim, jjp1) = radclatp * xprimm025(:iim)
243    
244        cvij1(:iim, jjp1) = 0.5 * rad * yprimu2(jjm)
245        cvij2(:iim, jjp1) = 0.
246        cvij3(:iim, jjp1) = 0.
247        cvij4(:iim, jjp1) = cvij1(:iim, jjp1)
248    
249        ! Périodicité :
250    
251        cvij1(iip1, :) = cvij1(1, :)
252        cvij2(iip1, :) = cvij2(1, :)
253        cvij3(iip1, :) = cvij3(1, :)
254        cvij4(iip1, :) = cvij4(1, :)
255    
256        cuij1(iip1, :) = cuij1(1, :)
257        cuij2(iip1, :) = cuij2(1, :)
258        cuij3(iip1, :) = cuij3(1, :)
259        cuij4(iip1, :) = cuij4(1, :)
260    
261        aireij1_2d(iip1, :) = aireij1_2d(1, :)
262        aireij2_2d(iip1, :) = aireij2_2d(1, :)
263        aireij3_2d(iip1, :) = aireij3_2d(1, :)
264        aireij4_2d(iip1, :) = aireij4_2d(1, :)
265    
266        DO j = 1, jjp1
267           DO i = 1, iim
268              aire_2d(i, j) = aireij1_2d(i, j) + aireij2_2d(i, j) &
269                   + aireij3_2d(i, j) + aireij4_2d(i, j)
270              alpha1_2d(i, j) = aireij1_2d(i, j) / aire_2d(i, j)
271              alpha2_2d(i, j) = aireij2_2d(i, j) / aire_2d(i, j)
272              alpha3_2d(i, j) = aireij3_2d(i, j) / aire_2d(i, j)
273              alpha4_2d(i, j) = aireij4_2d(i, j) / aire_2d(i, j)
274              alpha1p2_2d(i, j) = alpha1_2d(i, j) + alpha2_2d(i, j)
275              alpha1p4_2d(i, j) = alpha1_2d(i, j) + alpha4_2d(i, j)
276              alpha2p3_2d(i, j) = alpha2_2d(i, j) + alpha3_2d(i, j)
277              alpha3p4_2d(i, j) = alpha3_2d(i, j) + alpha4_2d(i, j)
278           END DO
279    
280           aire_2d(iip1, j) = aire_2d(1, j)
281           alpha1_2d(iip1, j) = alpha1_2d(1, j)
282           alpha2_2d(iip1, j) = alpha2_2d(1, j)
283           alpha3_2d(iip1, j) = alpha3_2d(1, j)
284           alpha4_2d(iip1, j) = alpha4_2d(1, j)
285           alpha1p2_2d(iip1, j) = alpha1p2_2d(1, j)
286           alpha1p4_2d(iip1, j) = alpha1p4_2d(1, j)
287           alpha2p3_2d(iip1, j) = alpha2p3_2d(1, j)
288           alpha3p4_2d(iip1, j) = alpha3p4_2d(1, j)
289        END DO
290    
291        DO j = 1, jjp1
292           DO i = 1, iim
293              aireu_2d(i, j) = aireij1_2d(i, j) + aireij2_2d(i, j) + &
294                   aireij4_2d(i + 1, j) + aireij3_2d(i + 1, j)
295              unsaire_2d(i, j) = 1. / aire_2d(i, j)
296              unsair_gam1_2d(i, j) = unsaire_2d(i, j)**(-gamdi_gdiv)
297              unsair_gam2_2d(i, j) = unsaire_2d(i, j)**(-gamdi_h)
298              airesurg_2d(i, j) = aire_2d(i, j) / g
299           END DO
300           aireu_2d(iip1, j) = aireu_2d(1, j)
301           unsaire_2d(iip1, j) = unsaire_2d(1, j)
302           unsair_gam1_2d(iip1, j) = unsair_gam1_2d(1, j)
303           unsair_gam2_2d(iip1, j) = unsair_gam2_2d(1, j)
304           airesurg_2d(iip1, j) = airesurg_2d(1, j)
305        END DO
306    
307        DO j = 1, jjm
308           DO i = 1, iim
309              airev_2d(i, j) = aireij2_2d(i, j) + aireij3_2d(i, j) + &
310                   aireij1_2d(i, j + 1) + aireij4_2d(i, j + 1)
311           END DO
312           DO i = 1, iim
313              airez = aireij2_2d(i, j) + aireij1_2d(i, j + 1) &
314                   + aireij3_2d(i + 1, j) + aireij4_2d(i + 1, j + 1)
315              unsairez_2d(i, j) = 1. / airez
316              unsairz_gam_2d(i, j) = unsairez_2d(i, j)**(-gamdi_grot)
317              fext_2d(i, j) = airez * sin(rlatv(j)) * 2. * omeg
318           END DO
319           airev_2d(iip1, j) = airev_2d(1, j)
320           unsairez_2d(iip1, j) = unsairez_2d(1, j)
321           fext_2d(iip1, j) = fext_2d(1, j)
322           unsairz_gam_2d(iip1, j) = unsairz_gam_2d(1, j)
323        END DO
324    
325        ! Calcul des élongations cu_2d, cv_2d
326    
327        DO j = 1, jjm
328           DO i = 1, iim
329              cv_2d(i, j) = 0.5 * &
330                   (cvij2(i, j) + cvij3(i, j) + cvij1(i, j + 1) + cvij4(i, j + 1))
331              unscv2_2d(i, j) = 1. / cv_2d(i, j)**2
332           END DO
333           DO i = 1, iim
334              cuvsurcv_2d(i, j) = airev_2d(i, j) * unscv2_2d(i, j)
335              cvsurcuv_2d(i, j) = 1. / cuvsurcv_2d(i, j)
336              cuvscvgam1_2d(i, j) = cuvsurcv_2d(i, j)**(-gamdi_gdiv)
337              cuvscvgam2_2d(i, j) = cuvsurcv_2d(i, j)**(-gamdi_h)
338              cvscuvgam_2d(i, j) = cvsurcuv_2d(i, j)**(-gamdi_grot)
339           END DO
340           cv_2d(iip1, j) = cv_2d(1, j)
341           unscv2_2d(iip1, j) = unscv2_2d(1, j)
342           cuvsurcv_2d(iip1, j) = cuvsurcv_2d(1, j)
343           cvsurcuv_2d(iip1, j) = cvsurcuv_2d(1, j)
344           cuvscvgam1_2d(iip1, j) = cuvscvgam1_2d(1, j)
345           cuvscvgam2_2d(iip1, j) = cuvscvgam2_2d(1, j)
346           cvscuvgam_2d(iip1, j) = cvscuvgam_2d(1, j)
347        END DO
348    
349        DO j = 2, jjm
350           DO i = 1, iim
351              cu_2d(i, j) = 0.5 * (cuij1(i, j) + cuij4(i + 1, j) + cuij2(i, j) &
352                   + cuij3(i + 1, j))
353              unscu2_2d(i, j) = 1. / cu_2d(i, j)**2
354              cvusurcu_2d(i, j) = aireu_2d(i, j) * unscu2_2d(i, j)
355              cusurcvu_2d(i, j) = 1. / cvusurcu_2d(i, j)
356              cvuscugam1_2d(i, j) = cvusurcu_2d(i, j)**(-gamdi_gdiv)
357              cvuscugam2_2d(i, j) = cvusurcu_2d(i, j)**(-gamdi_h)
358              cuscvugam_2d(i, j) = cusurcvu_2d(i, j)**(-gamdi_grot)
359           END DO
360           cu_2d(iip1, j) = cu_2d(1, j)
361           unscu2_2d(iip1, j) = unscu2_2d(1, j)
362           cvusurcu_2d(iip1, j) = cvusurcu_2d(1, j)
363           cusurcvu_2d(iip1, j) = cusurcvu_2d(1, j)
364           cvuscugam1_2d(iip1, j) = cvuscugam1_2d(1, j)
365           cvuscugam2_2d(iip1, j) = cvuscugam2_2d(1, j)
366           cuscvugam_2d(iip1, j) = cuscvugam_2d(1, j)
367        END DO
368    
369        ! Calcul aux pôles
370    
371        cu_2d(:, 1) = 0.
372        unscu2_2d(:, 1) = 0.
373    
374        cu_2d(:, jjp1) = 0.
375        unscu2_2d(:, jjp1) = 0.
376    
377        ! Calcul des aires aux pôles :
378    
379        apoln = sum(aire_2d(:iim, 1))
380        apols = sum(aire_2d(:iim, jjp1))
381        unsapolnga1 = 1. / (apoln**(-gamdi_gdiv))
382        unsapolsga1 = 1. / (apols**(-gamdi_gdiv))
383        unsapolnga2 = 1. / (apoln**(-gamdi_h))
384        unsapolsga2 = 1. / (apols**(-gamdi_h))
385    
386        ! Changement F. Hourdin calcul conservatif pour fext_2d
387        ! constang_2d contient le produit a * cos (latitude) * omega
388    
389        DO i = 1, iim
390           constang_2d(i, 1) = 0.
391        END DO
392        DO j = 1, jjm - 1
393           DO i = 1, iim
394              constang_2d(i, j + 1) = rad * omeg * cu_2d(i, j + 1) &
395                   * cos(rlatu(j + 1))
396           END DO
397        END DO
398        DO i = 1, iim
399           constang_2d(i, jjp1) = 0.
400        END DO
401    
402        ! Périodicité en longitude
403        DO j = 1, jjp1
404           constang_2d(iip1, j) = constang_2d(1, j)
405        END DO
406    
407      END SUBROUTINE inigeom
408    
409  end module comgeom  end module comgeom
Legend:
Removed from v.3  
changed lines
  Added in v.161
  ViewVC Help
Powered by ViewVC 1.1.21