source: LMDZ5/branches/testing/libf/dyn3d/calfis.F @ 1665

Last change on this file since 1665 was 1665, checked in by Laurent Fairhead, 12 years ago

Version testing basée sur la r1628

http://lmdz.lmd.jussieu.fr/utilisateurs/distribution-du-modele/versions-intermediaires


Testing release based on r1628

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 18.2 KB
Line 
1!
2! $Id: calfis.F 1665 2012-10-09 13:35:26Z fairhead $
3!
4C
5C
6      SUBROUTINE calfis(lafin,
7     $                  jD_cur, jH_cur,
8     $                  pucov,
9     $                  pvcov,
10     $                  pteta,
11     $                  pq,
12     $                  pmasse,
13     $                  pps,
14     $                  pp,
15     $                  ppk,
16     $                  pphis,
17     $                  pphi,
18     $                  pducov,
19     $                  pdvcov,
20     $                  pdteta,
21     $                  pdq,
22     $                  flxw,
23     $                  clesphy0,
24     $                  pdufi,
25     $                  pdvfi,
26     $                  pdhfi,
27     $                  pdqfi,
28     $                  pdpsfi)
29c
30c    Auteur :  P. Le Van, F. Hourdin
31c   .........
32      USE infotrac
33      USE control_mod
34 
35
36      IMPLICIT NONE
37c=======================================================================
38c
39c   1. rearrangement des tableaux et transformation
40c      variables dynamiques  >  variables physiques
41c   2. calcul des termes physiques
42c   3. retransformation des tendances physiques en tendances dynamiques
43c
44c   remarques:
45c   ----------
46c
47c    - les vents sont donnes dans la physique par leurs composantes
48c      naturelles.
49c    - la variable thermodynamique de la physique est une variable
50c      intensive :   T
51c      pour la dynamique on prend    T * ( preff / p(l) ) **kappa
52c    - les deux seules variables dependant de la geometrie necessaires
53c      pour la physique sont la latitude pour le rayonnement et
54c      l'aire de la maille quand on veut integrer une grandeur
55c      horizontalement.
56c    - les points de la physique sont les points scalaires de la
57c      la dynamique; numerotation:
58c          1 pour le pole nord
59c          (jjm-1)*iim pour l'interieur du domaine
60c          ngridmx pour le pole sud
61c      ---> ngridmx=2+(jjm-1)*iim
62c
63c     Input :
64c     -------
65c       pucov           covariant zonal velocity
66c       pvcov           covariant meridional velocity
67c       pteta           potential temperature
68c       pps             surface pressure
69c       pmasse          masse d'air dans chaque maille
70c       pts             surface temperature  (K)
71c       callrad         clef d'appel au rayonnement
72c
73c    Output :
74c    --------
75c        pdufi          tendency for the natural zonal velocity (ms-1)
76c        pdvfi          tendency for the natural meridional velocity
77c        pdhfi          tendency for the potential temperature
78c        pdtsfi         tendency for the surface temperature
79c
80c        pdtrad         radiative tendencies  \  both input
81c        pfluxrad       radiative fluxes      /  and output
82c
83c=======================================================================
84c
85c-----------------------------------------------------------------------
86c
87c    0.  Declarations :
88c    ------------------
89
90#include "dimensions.h"
91#include "paramet.h"
92#include "temps.h"
93
94      INTEGER ngridmx
95      PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm   )
96
97#include "comconst.h"
98#include "comvert.h"
99#include "comgeom2.h"
100#include "iniprint.h"
101
102c    Arguments :
103c    -----------
104      LOGICAL  lafin
105
106
107      REAL pvcov(iip1,jjm,llm)
108      REAL pucov(iip1,jjp1,llm)
109      REAL pteta(iip1,jjp1,llm)
110      REAL pmasse(iip1,jjp1,llm)
111      REAL pq(iip1,jjp1,llm,nqtot)
112      REAL pphis(iip1,jjp1)
113      REAL pphi(iip1,jjp1,llm)
114c
115      REAL pdvcov(iip1,jjm,llm)
116      REAL pducov(iip1,jjp1,llm)
117      REAL pdteta(iip1,jjp1,llm)
118      REAL pdq(iip1,jjp1,llm,nqtot)
119c
120      REAL pps(iip1,jjp1)
121      REAL pp(iip1,jjp1,llmp1)
122      REAL ppk(iip1,jjp1,llm)
123c
124      REAL pdvfi(iip1,jjm,llm)
125      REAL pdufi(iip1,jjp1,llm)
126      REAL pdhfi(iip1,jjp1,llm)
127      REAL pdqfi(iip1,jjp1,llm,nqtot)
128      REAL pdpsfi(iip1,jjp1)
129
130      INTEGER        longcles
131      PARAMETER    ( longcles = 20 )
132      REAL clesphy0( longcles )
133
134
135c    Local variables :
136c    -----------------
137
138      INTEGER i,j,l,ig0,ig,iq,iiq
139      REAL zpsrf(ngridmx)
140      REAL zplev(ngridmx,llm+1),zplay(ngridmx,llm)
141      REAL zphi(ngridmx,llm),zphis(ngridmx)
142c
143      REAL zufi(ngridmx,llm), zvfi(ngridmx,llm)
144      REAL ztfi(ngridmx,llm),zqfi(ngridmx,llm,nqtot)
145c
146      REAL pcvgu(ngridmx,llm), pcvgv(ngridmx,llm)
147      REAL pcvgt(ngridmx,llm), pcvgq(ngridmx,llm,2)
148c
149      REAL zdufi(ngridmx,llm),zdvfi(ngridmx,llm)
150      REAL zdtfi(ngridmx,llm),zdqfi(ngridmx,llm,nqtot)
151      REAL zdpsrf(ngridmx)
152c
153      REAL zdufic(ngridmx,llm),zdvfic(ngridmx,llm)
154      REAL zdtfic(ngridmx,llm),zdqfic(ngridmx,llm,nqtot)
155      REAL jH_cur_split,zdt_split
156      LOGICAL debut_split,lafin_split
157      INTEGER isplit
158
159      REAL zsin(iim),zcos(iim),z1(iim)
160      REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
161      REAL unskap, pksurcp
162c
163cIM diagnostique PVteta, Amip2
164      INTEGER ntetaSTD
165      PARAMETER(ntetaSTD=3)
166      REAL rtetaSTD(ntetaSTD)
167      DATA rtetaSTD/350., 380., 405./ ! Earth-specific values, beware !!
168      REAL PVteta(ngridmx,ntetaSTD)
169c
170      REAL flxw(iip1,jjp1,llm)  ! Flux de masse verticale sur la grille dynamique
171      REAL flxwfi(ngridmx,llm)  ! Flux de masse verticale sur la grille physiq
172c
173     
174      REAL SSUM
175
176      LOGICAL firstcal, debut
177      DATA firstcal/.true./
178      SAVE firstcal,debut
179!      REAL rdayvrai
180      REAL, intent(in):: jD_cur, jH_cur
181c
182c-----------------------------------------------------------------------
183c
184c    1. Initialisations :
185c    --------------------
186c
187c
188      IF ( firstcal )  THEN
189        debut = .TRUE.
190        IF (ngridmx.NE.2+(jjm-1)*iim) THEN
191         write(lunout,*) 'STOP dans calfis'
192         write(lunout,*)
193     &   'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim'
194         write(lunout,*) '  ngridmx  jjm   iim   '
195         write(lunout,*) ngridmx,jjm,iim
196         STOP
197        ENDIF
198      ELSE
199        debut = .FALSE.
200      ENDIF ! of IF (firstcal)
201
202c
203c
204c-----------------------------------------------------------------------
205c   40. transformation des variables dynamiques en variables physiques:
206c   ---------------------------------------------------------------
207
208c   41. pressions au sol (en Pascals)
209c   ----------------------------------
210
211       
212      zpsrf(1) = pps(1,1)
213
214      ig0  = 2
215      DO j = 2,jjm
216         CALL SCOPY( iim,pps(1,j),1,zpsrf(ig0), 1 )
217         ig0 = ig0+iim
218      ENDDO
219
220      zpsrf(ngridmx) = pps(1,jjp1)
221
222
223c   42. pression intercouches :
224c
225c   -----------------------------------------------------------------
226c     .... zplev  definis aux (llm +1) interfaces des couches  ....
227c     .... zplay  definis aux (  llm )    milieux des couches  ....
228c   -----------------------------------------------------------------
229
230c    ...    Exner = cp * ( p(l) / preff ) ** kappa     ....
231c
232       unskap   = 1./ kappa
233c
234      DO l = 1, llmp1
235        zplev( 1,l ) = pp(1,1,l)
236        ig0 = 2
237          DO j = 2, jjm
238             DO i =1, iim
239              zplev( ig0,l ) = pp(i,j,l)
240              ig0 = ig0 +1
241             ENDDO
242          ENDDO
243        zplev( ngridmx,l ) = pp(1,jjp1,l)
244      ENDDO
245c
246c
247
248c   43. temperature naturelle (en K) et pressions milieux couches .
249c   ---------------------------------------------------------------
250
251      DO l=1,llm
252
253         pksurcp     =  ppk(1,1,l) / cpp
254         zplay(1,l)  =  preff * pksurcp ** unskap
255         ztfi(1,l)   =  pteta(1,1,l) *  pksurcp
256         pcvgt(1,l)  =  pdteta(1,1,l) * pksurcp / pmasse(1,1,l)
257         ig0         = 2
258
259         DO j = 2, jjm
260            DO i = 1, iim
261              pksurcp        = ppk(i,j,l) / cpp
262              zplay(ig0,l)   = preff * pksurcp ** unskap
263              ztfi(ig0,l)    = pteta(i,j,l)  * pksurcp
264              pcvgt(ig0,l)   = pdteta(i,j,l) * pksurcp / pmasse(i,j,l)
265              ig0            = ig0 + 1
266            ENDDO
267         ENDDO
268
269         pksurcp       = ppk(1,jjp1,l) / cpp
270         zplay(ig0,l)  = preff * pksurcp ** unskap
271         ztfi (ig0,l)  = pteta(1,jjp1,l)  * pksurcp
272         pcvgt(ig0,l)  = pdteta(1,jjp1,l) * pksurcp/ pmasse(1,jjp1,l)
273
274      ENDDO
275
276c   43.bis traceurs
277c   ---------------
278c
279      DO iq=1,nqtot
280          iiq=niadv(iq)
281         DO l=1,llm
282            zqfi(1,l,iq) = pq(1,1,l,iiq)
283            ig0          = 2
284            DO j=2,jjm
285               DO i = 1, iim
286                  zqfi(ig0,l,iq)  = pq(i,j,l,iiq)
287                  ig0             = ig0 + 1
288               ENDDO
289            ENDDO
290            zqfi(ig0,l,iq) = pq(1,jjp1,l,iiq)
291         ENDDO
292      ENDDO
293
294c   convergence dynamique pour les traceurs "EAU"
295! Earth-specific treatment of first 2 tracers (water)
296       if (planet_type=="earth") then
297        DO iq=1,2
298         DO l=1,llm
299            pcvgq(1,l,iq)= pdq(1,1,l,iq) / pmasse(1,1,l)
300            ig0          = 2
301            DO j=2,jjm
302               DO i = 1, iim
303                  pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l)
304                  ig0             = ig0 + 1
305               ENDDO
306            ENDDO
307            pcvgq(ig0,l,iq)= pdq(1,jjp1,l,iq) / pmasse(1,jjp1,l)
308         ENDDO
309        ENDDO
310       endif ! of if (planet_type=="earth")
311
312
313c   Geopotentiel calcule par rapport a la surface locale:
314c   -----------------------------------------------------
315
316      CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pphi,zphi)
317      CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,pphis,zphis)
318      DO l=1,llm
319         DO ig=1,ngridmx
320           zphi(ig,l)=zphi(ig,l)-zphis(ig)
321         ENDDO
322      ENDDO
323
324c   ....  Calcul de la vitesse  verticale  ( en Pa*m*s  ou Kg/s )  ....
325c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux
326c    de masse est calclue dans advtrac.F 
327c      DO l=1,llm
328c        pvervel(1,l)=pw(1,1,l) * g /apoln
329c        ig0=2
330c       DO j=2,jjm
331c           DO i = 1, iim
332c              pvervel(ig0,l) = pw(i,j,l) * g * unsaire(i,j)
333c              ig0 = ig0 + 1
334c           ENDDO
335c       ENDDO
336c        pvervel(ig0,l)=pw(1,jjp1,l) * g /apols
337c      ENDDO
338
339c
340c   45. champ u:
341c   ------------
342
343      DO 50 l=1,llm
344
345         DO 25 j=2,jjm
346            ig0 = 1+(j-2)*iim
347            zufi(ig0+1,l)= 0.5 *
348     $      ( pucov(iim,j,l)/cu(iim,j) + pucov(1,j,l)/cu(1,j) )
349            pcvgu(ig0+1,l)= 0.5 *
350     $      ( pducov(iim,j,l)/cu(iim,j) + pducov(1,j,l)/cu(1,j) )
351            DO 10 i=2,iim
352               zufi(ig0+i,l)= 0.5 *
353     $         ( pucov(i-1,j,l)/cu(i-1,j) + pucov(i,j,l)/cu(i,j) )
354               pcvgu(ig0+i,l)= 0.5 *
355     $         ( pducov(i-1,j,l)/cu(i-1,j) + pducov(i,j,l)/cu(i,j) )
35610         CONTINUE
35725      CONTINUE
358
35950    CONTINUE
360
361
362c   46.champ v:
363c   -----------
364
365      DO l=1,llm
366         DO j=2,jjm
367            ig0=1+(j-2)*iim
368            DO i=1,iim
369               zvfi(ig0+i,l)= 0.5 *
370     $         ( pvcov(i,j-1,l)/cv(i,j-1) + pvcov(i,j,l)/cv(i,j) )
371               pcvgv(ig0+i,l)= 0.5 *
372     $         ( pdvcov(i,j-1,l)/cv(i,j-1) + pdvcov(i,j,l)/cv(i,j) )
373            ENDDO
374         ENDDO
375      ENDDO
376
377
378c   47. champs de vents aux pole nord   
379c   ------------------------------
380c        U = 1 / pi  *  integrale [ v * cos(long) * d long ]
381c        V = 1 / pi  *  integrale [ v * sin(long) * d long ]
382
383      DO l=1,llm
384
385         z1(1)   =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1)
386         z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,1,l)/cv(1,1)
387         DO i=2,iim
388            z1(i)   =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
389            z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1)
390         ENDDO
391
392         DO i=1,iim
393            zcos(i)   = COS(rlonv(i))*z1(i)
394            zcosbis(i)= COS(rlonv(i))*z1bis(i)
395            zsin(i)   = SIN(rlonv(i))*z1(i)
396            zsinbis(i)= SIN(rlonv(i))*z1bis(i)
397         ENDDO
398
399         zufi(1,l)  = SSUM(iim,zcos,1)/pi
400         pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi
401         zvfi(1,l)  = SSUM(iim,zsin,1)/pi
402         pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi
403
404      ENDDO
405
406
407c   48. champs de vents aux pole sud:
408c   ---------------------------------
409c        U = 1 / pi  *  integrale [ v * cos(long) * d long ]
410c        V = 1 / pi  *  integrale [ v * sin(long) * d long ]
411
412      DO l=1,llm
413
414         z1(1)   =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm)
415         z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,jjm,l)/cv(1,jjm)
416         DO i=2,iim
417            z1(i)   =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
418            z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm)
419         ENDDO
420
421         DO i=1,iim
422            zcos(i)    = COS(rlonv(i))*z1(i)
423            zcosbis(i) = COS(rlonv(i))*z1bis(i)
424            zsin(i)    = SIN(rlonv(i))*z1(i)
425            zsinbis(i) = SIN(rlonv(i))*z1bis(i)
426         ENDDO
427
428         zufi(ngridmx,l)  = SSUM(iim,zcos,1)/pi
429         pcvgu(ngridmx,l) = SSUM(iim,zcosbis,1)/pi
430         zvfi(ngridmx,l)  = SSUM(iim,zsin,1)/pi
431         pcvgv(ngridmx,l) = SSUM(iim,zsinbis,1)/pi
432
433      ENDDO
434c
435      if (planet_type=="earth") then
436#ifdef CPP_PHYS
437! PVtheta calls tetalevel, which is in the physics
438cIM calcul PV a teta=350, 380, 405K
439      CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta,
440     $           ztfi,zplay,zplev,
441     $           ntetaSTD,rtetaSTD,PVteta)
442#endif
443      endif
444c
445c On change de grille, dynamique vers physiq, pour le flux de masse verticale
446      CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,flxw,flxwfi)
447
448c-----------------------------------------------------------------------
449c   Appel de la physique:
450c   ---------------------
451
452
453
454!      write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys
455      zdt_split=dtphys/nsplit_phys
456      zdufic(:,:)=0.
457      zdvfic(:,:)=0.
458      zdtfic(:,:)=0.
459      zdqfic(:,:,:)=0.
460
461      if (planet_type=="earth") then
462#ifdef CPP_PHYS
463
464       do isplit=1,nsplit_phys
465
466         jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys)
467         debut_split=debut.and.isplit==1
468         lafin_split=lafin.and.isplit==nsplit_phys
469
470         CALL physiq (ngridmx,
471     .             llm,
472     .             debut_split,
473     .             lafin_split,
474     .             jD_cur,
475     .             jH_cur_split,
476     .             zdt_split,
477     .             zplev,
478     .             zplay,
479     .             zphi,
480     .             zphis,
481     .             presnivs,
482     .             clesphy0,
483     .             zufi,
484     .             zvfi,
485     .             ztfi,
486     .             zqfi,
487     .             flxwfi,
488     .             zdufi,
489     .             zdvfi,
490     .             zdtfi,
491     .             zdqfi,
492     .             zdpsrf,
493cIM diagnostique PVteta, Amip2         
494     .             pducov,
495     .             PVteta)
496
497         zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split
498         zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split
499         ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split
500         zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split
501
502         zdufic(:,:)=zdufic(:,:)+zdufi(:,:)
503         zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:)
504         zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:)
505         zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:)
506
507       enddo ! of do isplit=1,nsplit_phys
508
509#endif
510! of #ifdef CPP_PHYS
511      endif ! of if (planet_type=="earth")
512
513      zdufi(:,:)=zdufic(:,:)/nsplit_phys
514      zdvfi(:,:)=zdvfic(:,:)/nsplit_phys
515      zdtfi(:,:)=zdtfic(:,:)/nsplit_phys
516      zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys
517
518
519500   CONTINUE
520
521c-----------------------------------------------------------------------
522c   transformation des tendances physiques en tendances dynamiques:
523c   ---------------------------------------------------------------
524
525c  tendance sur la pression :
526c  -----------------------------------
527
528      CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zdpsrf,pdpsfi)
529c
530c   62. enthalpie potentielle
531c   ---------------------
532
533      DO l=1,llm
534
535         DO i=1,iip1
536          pdhfi(i,1,l)    = cpp *  zdtfi(1,l)      / ppk(i, 1  ,l)
537          pdhfi(i,jjp1,l) = cpp *  zdtfi(ngridmx,l)/ ppk(i,jjp1,l)
538         ENDDO
539
540         DO j=2,jjm
541            ig0=1+(j-2)*iim
542            DO i=1,iim
543               pdhfi(i,j,l) = cpp * zdtfi(ig0+i,l) / ppk(i,j,l)
544            ENDDO
545               pdhfi(iip1,j,l) =  pdhfi(1,j,l)
546         ENDDO
547
548      ENDDO
549
550
551c   62. humidite specifique
552c   ---------------------
553! Ehouarn: removed this useless bit: was overwritten at step 63 anyways
554!      DO iq=1,nqtot
555!         DO l=1,llm
556!            DO i=1,iip1
557!               pdqfi(i,1,l,iq)    = zdqfi(1,l,iq)
558!               pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,iq)
559!            ENDDO
560!            DO j=2,jjm
561!               ig0=1+(j-2)*iim
562!               DO i=1,iim
563!                  pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,iq)
564!               ENDDO
565!               pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,iq)
566!            ENDDO
567!         ENDDO
568!      ENDDO
569
570c   63. traceurs
571c   ------------
572C     initialisation des tendances
573      pdqfi(:,:,:,:)=0.
574C
575      DO iq=1,nqtot
576         iiq=niadv(iq)
577         DO l=1,llm
578            DO i=1,iip1
579               pdqfi(i,1,l,iiq)    = zdqfi(1,l,iq)
580               pdqfi(i,jjp1,l,iiq) = zdqfi(ngridmx,l,iq)
581            ENDDO
582            DO j=2,jjm
583               ig0=1+(j-2)*iim
584               DO i=1,iim
585                  pdqfi(i,j,l,iiq) = zdqfi(ig0+i,l,iq)
586               ENDDO
587               pdqfi(iip1,j,l,iiq) = pdqfi(1,j,l,iq)
588            ENDDO
589         ENDDO
590      ENDDO
591
592c   65. champ u:
593c   ------------
594
595      DO l=1,llm
596
597         DO i=1,iip1
598            pdufi(i,1,l)    = 0.
599            pdufi(i,jjp1,l) = 0.
600         ENDDO
601
602         DO j=2,jjm
603            ig0=1+(j-2)*iim
604            DO i=1,iim-1
605               pdufi(i,j,l)=
606     $         0.5*(zdufi(ig0+i,l)+zdufi(ig0+i+1,l))*cu(i,j)
607            ENDDO
608            pdufi(iim,j,l)=
609     $      0.5*(zdufi(ig0+1,l)+zdufi(ig0+iim,l))*cu(iim,j)
610            pdufi(iip1,j,l)=pdufi(1,j,l)
611         ENDDO
612
613      ENDDO
614
615
616c   67. champ v:
617c   ------------
618
619      DO l=1,llm
620
621         DO j=2,jjm-1
622            ig0=1+(j-2)*iim
623            DO i=1,iim
624               pdvfi(i,j,l)=
625     $         0.5*(zdvfi(ig0+i,l)+zdvfi(ig0+i+iim,l))*cv(i,j)
626            ENDDO
627            pdvfi(iip1,j,l) = pdvfi(1,j,l)
628         ENDDO
629      ENDDO
630
631
632c   68. champ v pres des poles:
633c   ---------------------------
634c      v = U * cos(long) + V * SIN(long)
635
636      DO l=1,llm
637
638         DO i=1,iim
639            pdvfi(i,1,l)=
640     $      zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i))
641            pdvfi(i,jjm,l)=zdufi(ngridmx,l)*COS(rlonv(i))
642     $      +zdvfi(ngridmx,l)*SIN(rlonv(i))
643            pdvfi(i,1,l)=
644     $      0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1)
645            pdvfi(i,jjm,l)=
646     $      0.5*(pdvfi(i,jjm,l)+zdvfi(ngridmx-iip1+i,l))*cv(i,jjm)
647          ENDDO
648
649         pdvfi(iip1,1,l)  = pdvfi(1,1,l)
650         pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
651
652      ENDDO
653
654c-----------------------------------------------------------------------
655
656700   CONTINUE
657 
658      firstcal = .FALSE.
659
660      RETURN
661      END
Note: See TracBrowser for help on using the repository browser.