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

Diff of /trunk/Sources/dyn3d/Guide/tau2alpha.f

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

trunk/dyn3d/tau2alpha.f90 revision 76 by guez, Fri Nov 15 18:45:49 2013 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      USE paramet_m, ONLY : iip1, jjp1    IMPLICIT NONE
     USE dimens_m, ONLY : jjm  
   
     IMPLICIT NONE  
   
     private iip1, jjp1, jjm  
   
     REAL lat_min_guide, lat_max_guide  
     REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)  
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 comgeom, ONLY: cu_2d, cv_2d, rlatu, rlatv
10        use conf_guide_m, only: lat_min_guide, lat_max_guide
11        USE dimens_m, ONLY: iim, jjm
12        USE nr_util, ONLY: pi
13        USE paramet_m, ONLY: iip1, jjp1
14        USE serre, ONLY: clat, clon, grossismx, grossismy
15        use writefield_m, only: writefield
16    
17      USE dimens_m, ONLY : iim      INTEGER, intent(in):: type
     USE nr_util, ONLY : pi  
     USE comgeom, ONLY : cu_2d, cv_2d, rlatu, rlatv  
     USE serre, ONLY : clat, clon, grossismx, grossismy  
   
     !   arguments :  
     INTEGER type  
     INTEGER pim, pjm  
18      REAL, intent(in):: 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)
31    
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 58  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
74    
75         CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL  ')         ! coordonnees du centre du zoom
        CALL dump2d(iip1, jjp1, dxdyu, 'DX2DY2 U     ')  
        CALL dump2d(iip1, jjp1, dxdyv, 'DX2DY2 v     ')  
   
        !   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 91  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.76  
changed lines
  Added in v.112
  ViewVC Help
Powered by ViewVC 1.1.21