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

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

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

trunk/libf/dyn3d/tau2alpha.f90 revision 37 by guez, Tue Dec 21 15:45:48 2010 UTC trunk/dyn3d/Guide/tau2alpha.f revision 114 by guez, Fri Sep 19 11:41:35 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 comconst, 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(:, :)
20      REAL dxdy_min, dxdy_max  
21        ! Local:
22      !  local :      REAL dxdy
23      REAL alphamin, alphamax, gamma, xi      REAL, save:: dxdy_min, dxdy_max
24      SAVE gamma      REAL alphamin, alphamax, xi
25        REAL, SAVE:: gamma
26      INTEGER i, j, ilon, ilat      INTEGER i, j, ilon, ilat
27        LOGICAL:: first = .TRUE.
28      LOGICAL first      REAL dx(iip1, jjp1), dy(iip1, jjp1)
     SAVE first  
     DATA first/ .TRUE./  
   
     REAL zdx(iip1, jjp1), zdy(iip1, jjp1)  
   
29      REAL zlat      REAL zlat
30      REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)      REAL, save:: 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))               dx(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)            dx(1, j) = dx(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))               dy(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
47            zdx(i, 1) = zdx(i, 2)            dx(i, 1) = dx(i, 2)
48            zdx(i, jjp1) = zdx(i, jjm)            dx(i, jjp1) = dx(i, jjm)
49            zdy(i, 1) = zdy(i, 2)            dy(i, 1) = dy(i, 2)
50            zdy(i, jjp1) = zdy(i, jjm)            dy(i, jjp1) = dy(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(dx(i, j)**2 + dy(i, j)**2)
56            END DO            END DO
57         END DO         END DO
58           CALL writefield("dxdys", dxdys)
59    
60         DO j = 1, jjp1         DO j = 1, jjp1
61            DO i = 1, iim            DO i = 1, iim
62               dxdyu(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))               dxdyu(i, j) = 0.5 * (dxdys(i, j) + dxdys(i + 1, j))
63            END DO            END DO
64            dxdyu(iip1, j) = dxdyu(1, j)            dxdyu(iip1, j) = dxdyu(1, j)
65         END DO         END DO
66    
67         DO j = 1, jjm         DO j = 1, jjm
68            DO i = 1, iip1            DO i = 1, iip1
69               dxdyv(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))               dxdyv(i, j) = 0.5 * (dxdys(i, j) + dxdys(i, j + 1))
70            END DO            END DO
71         END DO         END DO
72    
73         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  
74         CALL coordij(clon, clat, ilon, ilat)         CALL coordij(clon, clat, ilon, ilat)
75         !   aire de la maille au centre du zoom         ! aire de la maille au centre du zoom
76         dxdy_min = dxdys(ilon, ilat)         dxdy_min = dxdys(ilon, ilat)
77         !   dxdy maximale de la maille  
78           ! dxdy maximal de la maille :
79         dxdy_max = 0.         dxdy_max = 0.
80         DO j = 1, jjp1         DO j = 1, jjp1
81            DO i = 1, iip1            DO i = 1, iip1
# Line 88  contains Line 83  contains
83            END DO            END DO
84         END DO         END DO
85    
86         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
87            PRINT *, 'ATTENTION modele peu zoome'            PRINT *, 'Attention : modèle peu zoomé.'
88            PRINT *, 'ATTENTION on prend une constante de guidage cste'            PRINT *, 'On prend une constante de guidage constante.'
           gamma = 0.  
89         ELSE         ELSE
90            gamma = (dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)            gamma = (dxdy_max - 2. * dxdy_min) / (dxdy_max - dxdy_min)
91            PRINT *, 'gamma=', gamma            IF (gamma < 1E-5) THEN
92            IF (gamma<1.E-5) THEN               PRINT *, '(dxdy_max - 2. * dxdy_min) / (dxdy_max - dxdy_min) ' &
93               PRINT *, 'gamma =', gamma, '<1e-5'                    // '< 1e-5'
94               STOP               STOP 1
95            END IF            END IF
96              gamma = log(0.5) / log(gamma)
97            PRINT *, 'gamma=', gamma            PRINT *, 'gamma=', gamma
           gamma = log(0.5)/log(gamma)  
98         END IF         END IF
99           first = .false.
100      END IF      END IF
101    
102      alphamin = factt/taumax      alphamin = factt / taumax
103      alphamax = factt/taumin      alphamax = factt / taumin
104    
105      DO j = 1, pjm      DO j = 1, size(alpha, 2)
106         DO i = 1, pim         DO i = 1, size(alpha, 1)
107            IF (type==1) THEN            IF (type==1) THEN
108               dxdy_ = dxdys(i, j)               dxdy = dxdys(i, j)
109               zlat = rlatu(j)*180./pi               zlat = rlatu(j) * 180. / pi
110            ELSE IF (type==2) THEN            ELSE IF (type==2) THEN
111               dxdy_ = dxdyu(i, j)               dxdy = dxdyu(i, j)
112               zlat = rlatu(j)*180./pi               zlat = rlatu(j) * 180. / pi
113            ELSE IF (type==3) THEN            ELSE IF (type==3) THEN
114               dxdy_ = dxdyv(i, j)               dxdy = dxdyv(i, j)
115               zlat = rlatv(j)*180./pi               zlat = rlatv(j) * 180. / pi
116            END IF            END IF
117            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
118               !  pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin               ! grille regulière
119               alpha(i, j) = alphamin               alpha(i, j) = alphamin
120            ELSE            ELSE
121               xi = ((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma               xi = ((dxdy_max - dxdy) / (dxdy_max - dxdy_min))**gamma
122               xi = min(xi, 1.)               xi = min(xi, 1.)
123               IF (lat_min_guide<=zlat .AND. zlat<=lat_max_guide) THEN               IF (lat_min_guide <= zlat .AND. zlat <= lat_max_guide) THEN
124                  alpha(i, j) = xi*alphamin + (1.-xi)*alphamax                  alpha(i, j) = xi * alphamin + (1. - xi) * alphamax
125               ELSE               ELSE
126                  alpha(i, j) = 0.                  alpha(i, j) = 0.
127               END IF               END IF
Legend:
Removed from v.37  
changed lines
  Added in v.114
  ViewVC Help
Powered by ViewVC 1.1.21