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

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

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

revision 112 by guez, Thu Sep 18 13:36:51 2014 UTC revision 120 by guez, Tue Jan 13 14:56:15 2015 UTC
# Line 4  module fyhyp_m Line 4  module fyhyp_m
4    
5  contains  contains
6    
7    SUBROUTINE fyhyp(yzoomdeg, grossism, dzooma, tau, rrlatu, yyprimu, rrlatv, &    SUBROUTINE fyhyp(rlatu, yyprimu, rlatv, rlatu2, yprimu2, rlatu1, yprimu1)
        yyprimv, rlatu2, yprimu2, rlatu1, yprimu1, champmin, champmax)  
8    
9      ! From LMDZ4/libf/dyn3d/fyhyp.F, version 1.2, 2005/06/03 09:11:32      ! From LMDZ4/libf/dyn3d/fyhyp.F, version 1.2, 2005/06/03 09:11:32
10    
11      ! Author: P. Le Van, from analysis by R. Sadourny      ! Author: P. Le Van, from analysis by R. Sadourny
12    
13      ! Calcule les latitudes et dérivées dans la grille du GCM pour une      ! Calcule les latitudes et dérivées dans la grille du GCM pour une
14      ! fonction f(y) à tangente hyperbolique.      ! fonction f(y) à dérivée tangente hyperbolique.
15    
16      ! Nota bene : il vaut mieux avoir grossism * dzoom < pi / 2 (rad),      ! Nota bene : il vaut mieux avoir grossismy * dzoomy < pi / 2 (radians).
     ! en latitude.  
17    
18        use coefpoly_m, only: coefpoly
19      USE dimens_m, only: jjm      USE dimens_m, only: jjm
20      USE paramet_m, only: JJP1      use serre, only: clat, grossismy, dzoomy, tauy
21    
22      REAL, intent(in):: yzoomdeg      REAL, intent(out):: rlatu(jjm + 1), yyprimu(jjm + 1)
23        REAL, intent(out):: rlatv(jjm)
24      REAL, intent(in):: grossism      real, intent(out):: rlatu2(jjm), yprimu2(jjm), rlatu1(jjm), yprimu1(jjm)
     ! grossissement (= 2 si 2 fois, = 3 si 3 fois, etc.)  
   
     REAL, intent(in):: dzooma  
   
     REAL, intent(in):: tau  
     ! raideur de la transition de l'intérieur à l'extérieur du zoom  
   
     ! arguments de sortie  
   
     REAL rrlatu(jjp1), yyprimu(jjp1), rrlatv(jjm), yyprimv(jjm)  
     real rlatu2(jjm), yprimu2(jjm), rlatu1(jjm), yprimu1(jjm)  
     DOUBLE PRECISION champmin, champmax  
25    
26      ! Local:      ! Local:
27    
28        DOUBLE PRECISION champmin, champmax
29      INTEGER, PARAMETER:: nmax=30000, nmax2=2*nmax      INTEGER, PARAMETER:: nmax=30000, nmax2=2*nmax
30      REAL dzoom ! distance totale de la zone du zoom (en radians)      REAL dzoom ! distance totale de la zone du zoom (en radians)
31      DOUBLE PRECISION ylat(jjp1), yprim(jjp1)      DOUBLE PRECISION ylat(jjm + 1), yprim(jjm + 1)
32      DOUBLE PRECISION yuv      DOUBLE PRECISION yuv
33      DOUBLE PRECISION, save:: yt(0:nmax2)      DOUBLE PRECISION, save:: yt(0:nmax2)
34      DOUBLE PRECISION fhyp(0:nmax2), beta      DOUBLE PRECISION fhyp(0:nmax2), beta
# Line 48  contains Line 36  contains
36      DOUBLE PRECISION fxm(0:nmax2)      DOUBLE PRECISION fxm(0:nmax2)
37      DOUBLE PRECISION, save:: yf(0:nmax2)      DOUBLE PRECISION, save:: yf(0:nmax2)
38      DOUBLE PRECISION yypr(0:nmax2)      DOUBLE PRECISION yypr(0:nmax2)
39      DOUBLE PRECISION yvrai(jjp1), yprimm(jjp1), ylatt(jjp1)      DOUBLE PRECISION yvrai(jjm + 1), yprimm(jjm + 1), ylatt(jjm + 1)
40      DOUBLE PRECISION pi, pis2, epsilon, y0, pisjm      DOUBLE PRECISION pi, pis2, epsilon, y0, pisjm
41      DOUBLE PRECISION yo1, yi, ylon2, ymoy, yprimin      DOUBLE PRECISION yo1, yi, ylon2, ymoy, yprimin
42      DOUBLE PRECISION yfi, yf1, ffdy      DOUBLE PRECISION yfi, yf1, ffdy
# Line 66  contains Line 54  contains
54    
55      !-------------------------------------------------------------------      !-------------------------------------------------------------------
56    
57        print *, "Call sequence information: fyhyp"
58    
59      pi = 2.*asin(1.)      pi = 2.*asin(1.)
60      pis2 = pi/2.      pis2 = pi/2.
61      pisjm = pi/real(jjm)      pisjm = pi/real(jjm)
62      epsilon = 1e-3      epsilon = 1e-3
63      y0 = yzoomdeg*pi/180.      y0 = clat*pi/180.
64    
65      IF (dzooma<1.) THEN      IF (dzoomy<1.) THEN
66         dzoom = dzooma*pi         dzoom = dzoomy*pi
67      ELSE IF (dzooma<12.) THEN      ELSE IF (dzoomy<12.) THEN
68         print *, "Le paramètre dzoomy pour fyhyp est trop petit. L'augmenter " &         print *, "Le paramètre dzoomy pour fyhyp est trop petit. L'augmenter " &
69              // "et relancer."              // "et relancer."
70         STOP 1         STOP 1
71      ELSE      ELSE
72         dzoom = dzooma * pi/180.         dzoom = dzoomy * pi/180.
73      END IF      END IF
74    
75      print *, 'yzoom(rad), grossism, tau, dzoom (rad):'      print *, 'yzoom(rad), grossismy, tauy, dzoom (rad):'
76      print *, y0, grossism, tau, dzoom      print *, y0, grossismy, tauy, dzoom
77    
78      DO i = 0, nmax2      DO i = 0, nmax2
79         yt(i) = -pis2 + real(i)*pi/nmax2         yt(i) = -pis2 + real(i)*pi/nmax2
# Line 99  contains Line 89  contains
89    
90      DO i = 0, nmax2      DO i = 0, nmax2
91         IF (yt(i)<y0) THEN         IF (yt(i)<y0) THEN
92            fa(i) = tau*(yt(i)-y0 + dzoom/2.)            fa(i) = tauy*(yt(i)-y0 + dzoom/2.)
93            fb(i) = (yt(i)-2.*y0*heavyy0m + pis2)*(y0-yt(i))            fb(i) = (yt(i)-2.*y0*heavyy0m + pis2)*(y0-yt(i))
94         ELSE IF (yt(i)>y0) THEN         ELSE IF (yt(i)>y0) THEN
95            fa(i) = tau*(y0-yt(i) + dzoom/2.)            fa(i) = tauy*(y0-yt(i) + dzoom/2.)
96            fb(i) = (2.*y0*heavyy0-yt(i) + pis2)*(yt(i)-y0)            fb(i) = (2.*y0*heavyy0-yt(i) + pis2)*(yt(i)-y0)
97         END IF         END IF
98    
# Line 125  contains Line 115  contains
115      DO i = 1, nmax2      DO i = 1, nmax2
116         ymoy = 0.5*(yt(i-1) + yt(i))         ymoy = 0.5*(yt(i-1) + yt(i))
117         IF (ymoy<y0) THEN         IF (ymoy<y0) THEN
118            fa(i) = tau*(ymoy-y0 + dzoom/2.)            fa(i) = tauy*(ymoy-y0 + dzoom/2.)
119            fb(i) = (ymoy-2.*y0*heavyy0m + pis2)*(y0-ymoy)            fb(i) = (ymoy-2.*y0*heavyy0m + pis2)*(y0-ymoy)
120         ELSE IF (ymoy>y0) THEN         ELSE IF (ymoy>y0) THEN
121            fa(i) = tau*(y0-ymoy + dzoom/2.)            fa(i) = tauy*(y0-ymoy + dzoom/2.)
122            fb(i) = (2.*y0*heavyy0-ymoy + pis2)*(ymoy-y0)            fb(i) = (2.*y0*heavyy0-ymoy + pis2)*(ymoy-y0)
123         END IF         END IF
124    
# Line 144  contains Line 134  contains
134         ffdy = ffdy + fxm(i)*(yt(i)-yt(i-1))         ffdy = ffdy + fxm(i)*(yt(i)-yt(i-1))
135      END DO      END DO
136    
137      beta = (grossism*ffdy-pi)/(ffdy-pi)      beta = (grossismy*ffdy-pi)/(ffdy-pi)
138    
139      IF (2. * beta - grossism <= 0.) THEN      IF (2. * beta - grossismy <= 0.) THEN
140         print *, 'Attention ! La valeur beta calculee dans la routine fyhyp ' &         print *, 'Attention ! La valeur beta calculee dans la routine fyhyp ' &
141              // 'est mauvaise. Modifier les valeurs de grossismy, tauy ou ' &              // 'est mauvaise. Modifier les valeurs de grossismy, tauy ou ' &
142              // 'dzoomy et relancer.'              // 'dzoomy et relancer.'
# Line 156  contains Line 146  contains
146      ! calcul de Ytprim      ! calcul de Ytprim
147    
148      DO i = 0, nmax2      DO i = 0, nmax2
149         ytprim(i) = beta + (grossism-beta)*fhyp(i)         ytprim(i) = beta + (grossismy-beta)*fhyp(i)
150      END DO      END DO
151    
152      ! Calcul de Yf      ! Calcul de Yf
153    
154      yf(0) = -pis2      yf(0) = -pis2
155      DO i = 1, nmax2      DO i = 1, nmax2
156         yypr(i) = beta + (grossism-beta)*fxm(i)         yypr(i) = beta + (grossismy-beta)*fxm(i)
157      END DO      END DO
158    
159      DO i = 1, nmax2      DO i = 1, nmax2
# Line 245  contains Line 235  contains
235    
236         IF (ik==1) THEN         IF (ik==1) THEN
237            ypn = pis2            ypn = pis2
238            DO j = jlat, 1, -1            DO j = jjm + 1, 1, -1
239               IF (yvrai(j)<=ypn) exit               IF (yvrai(j)<=ypn) exit
240            END DO            END DO
241    
# Line 279  contains Line 269  contains
269         END DO         END DO
270    
271         IF (ik==1) THEN         IF (ik==1) THEN
272            DO j = 1, jlat            DO j = 1, jjm + 1
273               rrlatu(j) = ylat(j)               rlatu(j) = ylat(j)
274               yyprimu(j) = yprim(j)               yyprimu(j) = yprim(j)
275            END DO            END DO
276         ELSE IF (ik==2) THEN         ELSE IF (ik==2) THEN
277            DO j = 1, jlat            DO j = 1, jjm
278               rrlatv(j) = ylat(j)               rlatv(j) = ylat(j)
              yyprimv(j) = yprim(j)  
279            END DO            END DO
280         ELSE IF (ik==3) THEN         ELSE IF (ik==3) THEN
281            DO j = 1, jlat            DO j = 1, jjm
282               rlatu2(j) = ylat(j)               rlatu2(j) = ylat(j)
283               yprimu2(j) = yprim(j)               yprimu2(j) = yprim(j)
284            END DO            END DO
285         ELSE IF (ik==4) THEN         ELSE IF (ik==4) THEN
286            DO j = 1, jlat            DO j = 1, jjm
287               rlatu1(j) = ylat(j)               rlatu1(j) = ylat(j)
288               yprimu1(j) = yprim(j)               yprimu1(j) = yprim(j)
289            END DO            END DO
# Line 302  contains Line 291  contains
291      END DO loop_ik      END DO loop_ik
292    
293      DO j = 1, jjm      DO j = 1, jjm
294         ylat(j) = rrlatu(j) - rrlatu(j + 1)         ylat(j) = rlatu(j) - rlatu(j + 1)
295      END DO      END DO
296      champmin = 1e12      champmin = 1e12
297      champmax = -1e12      champmax = -1e12
# Line 313  contains Line 302  contains
302      champmin = champmin*180./pi      champmin = champmin*180./pi
303      champmax = champmax*180./pi      champmax = champmax*180./pi
304    
305        DO j = 1, jjm
306           IF (rlatu1(j) <= rlatu2(j)) THEN
307              print *, 'Attention ! rlatu1 < rlatu2 ', rlatu1(j), rlatu2(j), j
308              STOP 13
309           ENDIF
310    
311           IF (rlatu2(j) <= rlatu(j+1)) THEN
312              print *, 'Attention ! rlatu2 < rlatup1 ', rlatu2(j), rlatu(j+1), j
313              STOP 14
314           ENDIF
315    
316           IF (rlatu(j) <= rlatu1(j)) THEN
317              print *, ' Attention ! rlatu < rlatu1 ', rlatu(j), rlatu1(j), j
318              STOP 15
319           ENDIF
320    
321           IF (rlatv(j) <= rlatu2(j)) THEN
322              print *, ' Attention ! rlatv < rlatu2 ', rlatv(j), rlatu2(j), j
323              STOP 16
324           ENDIF
325    
326           IF (rlatv(j) >= rlatu1(j)) THEN
327              print *, ' Attention ! rlatv > rlatu1 ', rlatv(j), rlatu1(j), j
328              STOP 17
329           ENDIF
330    
331           IF (rlatv(j) >= rlatu(j)) THEN
332              print *, ' Attention ! rlatv > rlatu ', rlatv(j), rlatu(j), j
333              STOP 18
334           ENDIF
335        ENDDO
336    
337        print *, 'Latitudes'
338        print 3, champmin, champmax
339    
340    3   Format(1x, ' Au centre du zoom, la longueur de la maille est', &
341             ' d environ ', f0.2, ' degres ', /, &
342             ' alors que la maille en dehors de la zone du zoom est ', &
343             "d'environ ", f0.2, ' degres ')
344    
345    END SUBROUTINE fyhyp    END SUBROUTINE fyhyp
346    
347  end module fyhyp_m  end module fyhyp_m
Legend:
Removed from v.112  
changed lines
  Added in v.120
  ViewVC Help
Powered by ViewVC 1.1.21