/[lmdze]/trunk/dyn3d/tau2alpha.f
ViewVC logotype

Diff of /trunk/dyn3d/tau2alpha.f

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

trunk/libf/dyn3d/tau2alpha.f90 revision 39 by guez, Tue Jan 25 15:11:05 2011 UTC trunk/dyn3d/tau2alpha.f revision 112 by guez, Thu Sep 18 13:36:51 2014 UTC
# Line 1  Line 1 
1  module tau2alpha_m  module tau2alpha_m
2    
3      IMPLICIT NONE    IMPLICIT NONE
   
     REAL lat_min_guide, lat_max_guide  
4    
5  contains  contains
6    
7    SUBROUTINE tau2alpha(type, pim, pjm, factt, taumin, taumax, alpha)    SUBROUTINE tau2alpha(type, factt, taumin, taumax, alpha)
8    
9      USE dimens_m, ONLY : iim, jjm      USE comgeom, ONLY: cu_2d, cv_2d, rlatu, rlatv
10      USE paramet_m, ONLY : iip1, jjp1      use conf_guide_m, only: lat_min_guide, lat_max_guide
11      USE nr_util, ONLY : pi      USE dimens_m, ONLY: iim, jjm
12      USE comgeom, ONLY : cu_2d, cv_2d, rlatu, rlatv      USE nr_util, ONLY: pi
13      USE serre, ONLY : clat, clon, grossismx, grossismy      USE paramet_m, ONLY: iip1, jjp1
14        USE serre, ONLY: clat, clon, grossismx, grossismy
15      !   arguments :      use writefield_m, only: writefield
16      INTEGER type  
17      INTEGER pim, pjm      INTEGER, intent(in):: type
18      REAL factt, taumin, taumax      REAL, intent(in):: factt, taumin, taumax
19      REAL dxdy_, alpha(pim, pjm)      real, intent(out):: alpha(:, :)
     REAL dxdy_min, dxdy_max  
20    
21      !  local :      ! Local:
22      REAL alphamin, alphamax, gamma, xi      REAL dxdy
23      SAVE gamma      REAL dxdy_min, dxdy_max
24        REAL alphamin, alphamax, xi
25        REAL, SAVE:: gamma
26      INTEGER i, j, ilon, ilat      INTEGER i, j, ilon, ilat
27        LOGICAL:: first = .TRUE.
     LOGICAL first  
     SAVE first  
     DATA first/ .TRUE./  
   
28      REAL zdx(iip1, jjp1), zdy(iip1, jjp1)      REAL zdx(iip1, jjp1), zdy(iip1, jjp1)
   
29      REAL zlat      REAL zlat
30      REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)      REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)
31      COMMON /comdxdy/dxdys, dxdyu, dxdyv  
32        !------------------------------------------------------------
33    
34      IF (first) THEN      IF (first) THEN
35         DO j = 2, jjm         DO j = 2, jjm
36            DO i = 2, iip1            DO i = 2, iip1
37               zdx(i, j) = 0.5*(cu_2d(i-1, j)+cu_2d(i, j))/cos(rlatu(j))               zdx(i, j) = 0.5 * (cu_2d(i - 1, j) + cu_2d(i, j)) / cos(rlatu(j))
38            END DO            END DO
39            zdx(1, j) = zdx(iip1, j)            zdx(1, j) = zdx(iip1, j)
40         END DO         END DO
41         DO j = 2, jjm         DO j = 2, jjm
42            DO i = 1, iip1            DO i = 1, iip1
43               zdy(i, j) = 0.5*(cv_2d(i, j-1)+cv_2d(i, j))               zdy(i, j) = 0.5 * (cv_2d(i, j - 1) + cv_2d(i, j))
44            END DO            END DO
45         END DO         END DO
46         DO i = 1, iip1         DO i = 1, iip1
# Line 54  contains Line 49  contains
49            zdy(i, 1) = zdy(i, 2)            zdy(i, 1) = zdy(i, 2)
50            zdy(i, jjp1) = zdy(i, jjm)            zdy(i, jjp1) = zdy(i, jjm)
51         END DO         END DO
52    
53         DO j = 1, jjp1         DO j = 1, jjp1
54            DO i = 1, iip1            DO i = 1, iip1
55               dxdys(i, j) = sqrt(zdx(i, j)*zdx(i, j)+zdy(i, j)*zdy(i, j))               dxdys(i, j) = sqrt(zdx(i, j)**2 + zdy(i, j)**2)
56            END DO            END DO
57         END DO         END DO
58         DO j = 1, jjp1         CALL writefield("dxdys", dxdys)
59            DO i = 1, iim  
60               dxdyu(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))         if (type == 2) then
61              DO j = 1, jjp1
62                 DO i = 1, iim
63                    dxdyu(i, j) = 0.5 * (dxdys(i, j) + dxdys(i + 1, j))
64                 END DO
65                 dxdyu(iip1, j) = dxdyu(1, j)
66            END DO            END DO
67            dxdyu(iip1, j) = dxdyu(1, j)         elseif (type == 3) then
68         END DO            DO j = 1, jjm
69         DO j = 1, jjm               DO i = 1, iip1
70            DO i = 1, iip1                  dxdyv(i, j) = 0.5 * (dxdys(i, j) + dxdys(i, j + 1))
71               dxdyv(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))               END DO
72            END DO            END DO
73         END DO         end if
   
        CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL  ')  
        CALL dump2d(iip1, jjp1, dxdyu, 'DX2DY2 U     ')  
        CALL dump2d(iip1, jjp1, dxdyv, 'DX2DY2 v     ')  
74    
75         !   coordonnees du centre du zoom         ! coordonnees du centre du zoom
76         CALL coordij(clon, clat, ilon, ilat)         CALL coordij(clon, clat, ilon, ilat)
77         !   aire de la maille au centre du zoom         ! aire de la maille au centre du zoom
78         dxdy_min = dxdys(ilon, ilat)         dxdy_min = dxdys(ilon, ilat)
79         !   dxdy maximale de la maille  
80           ! dxdy maximal de la maille :
81         dxdy_max = 0.         dxdy_max = 0.
82         DO j = 1, jjp1         DO j = 1, jjp1
83            DO i = 1, iip1            DO i = 1, iip1
# Line 87  contains Line 85  contains
85            END DO            END DO
86         END DO         END DO
87    
88         IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN         IF (abs(grossismx - 1.) < 0.1 .OR. abs(grossismy - 1.) < 0.1) THEN
89            PRINT *, 'ATTENTION modele peu zoome'            PRINT *, 'Attention : modèle peu zoomé.'
90            PRINT *, 'ATTENTION on prend une constante de guidage cste'            PRINT *, 'On prend une constante de guidage constante.'
           gamma = 0.  
91         ELSE         ELSE
92            gamma = (dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)            gamma = (dxdy_max - 2. * dxdy_min) / (dxdy_max - dxdy_min)
93            PRINT *, 'gamma=', gamma            PRINT *, 'gamma=', gamma
94            IF (gamma<1.E-5) THEN            IF (gamma<1.E-5) THEN
95               PRINT *, 'gamma =', gamma, '<1e-5'               PRINT *, 'gamma =', gamma, '<1e-5'
96               STOP               STOP
97            END IF            END IF
98            PRINT *, 'gamma=', gamma            PRINT *, 'gamma=', gamma
99            gamma = log(0.5)/log(gamma)            gamma = log(0.5) / log(gamma)
100         END IF         END IF
101           first = .false.
102      END IF      END IF
103    
104      alphamin = factt/taumax      alphamin = factt / taumax
105      alphamax = factt/taumin      alphamax = factt / taumin
106    
107      DO j = 1, pjm      DO j = 1, size(alpha, 2)
108         DO i = 1, pim         DO i = 1, size(alpha, 1)
109            IF (type==1) THEN            IF (type==1) THEN
110               dxdy_ = dxdys(i, j)               dxdy = dxdys(i, j)
111               zlat = rlatu(j)*180./pi               zlat = rlatu(j) * 180. / pi
112            ELSE IF (type==2) THEN            ELSE IF (type==2) THEN
113               dxdy_ = dxdyu(i, j)               dxdy = dxdyu(i, j)
114               zlat = rlatu(j)*180./pi               zlat = rlatu(j) * 180. / pi
115            ELSE IF (type==3) THEN            ELSE IF (type==3) THEN
116               dxdy_ = dxdyv(i, j)               dxdy = dxdyv(i, j)
117               zlat = rlatv(j)*180./pi               zlat = rlatv(j) * 180. / pi
118            END IF            END IF
119            IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN            IF (abs(grossismx - 1.) < 0.1 .OR. abs(grossismy - 1.) < 0.1) THEN
120               !  pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin               ! grille regulière
121               alpha(i, j) = alphamin               alpha(i, j) = alphamin
122            ELSE            ELSE
123               xi = ((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma               xi = ((dxdy_max - dxdy) / (dxdy_max - dxdy_min))**gamma
124               xi = min(xi, 1.)               xi = min(xi, 1.)
125               IF (lat_min_guide<=zlat .AND. zlat<=lat_max_guide) THEN               IF (lat_min_guide <= zlat .AND. zlat <= lat_max_guide) THEN
126                  alpha(i, j) = xi*alphamin + (1.-xi)*alphamax                  alpha(i, j) = xi * alphamin + (1. - xi) * alphamax
127               ELSE               ELSE
128                  alpha(i, j) = 0.                  alpha(i, j) = 0.
129               END IF               END IF
Legend:
Removed from v.39  
changed lines
  Added in v.112
  ViewVC Help
Powered by ViewVC 1.1.21