/[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/libf/dyn3d/tau2alpha.f90 revision 37 by guez, Tue Dec 21 15:45:48 2010 UTC trunk/dyn3d/Guide/tau2alpha.f revision 115 by guez, Fri Sep 19 17:36:20 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(dxdy, rlat, taumin, taumax, alpha)
     !=======================================================================  
   
     USE dimens_m, ONLY : iim, jjm  
     USE paramet_m, ONLY : iip1, jjp1  
     USE comconst, ONLY : pi  
     USE comgeom, ONLY : cu_2d, cv_2d, rlatu, rlatv  
     USE serre, ONLY : clat, clon, grossismx, grossismy  
   
     !   arguments :  
     INTEGER type  
     INTEGER pim, pjm  
     REAL factt, taumin, taumax  
     REAL dxdy_, alpha(pim, pjm)  
     REAL dxdy_min, dxdy_max  
   
     !  local :  
     REAL alphamin, alphamax, gamma, xi  
     SAVE gamma  
     INTEGER i, j, ilon, ilat  
   
     LOGICAL first  
     SAVE first  
     DATA first/ .TRUE./  
   
     REAL zdx(iip1, jjp1), zdy(iip1, jjp1)  
8    
9        use conf_guide_m, only: lat_min_guide, lat_max_guide, factt
10        use init_tau2alpha_m, only: dxdy_min, dxdy_max, gamma
11        USE nr_util, ONLY: pi, assert_eq
12    
13        REAL, intent(in):: dxdy(:, :) ! (n_lon, n_lat)
14        REAL, intent(in):: rlat(:) ! (n_lat)
15        REAL, intent(in):: taumin, taumax
16        real, intent(out):: alpha(:, :) ! (n_lon, n_lat)
17    
18        ! Local:
19        REAL alphamin, alphamax, xi
20        INTEGER i, j, n_lon, n_lat
21      REAL zlat      REAL zlat
     REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)  
     COMMON /comdxdy/dxdys, dxdyu, dxdyv  
22    
23      IF (first) THEN      !------------------------------------------------------------
        DO j = 2, jjm  
           DO i = 2, iip1  
              zdx(i, j) = 0.5*(cu_2d(i-1, j)+cu_2d(i, j))/cos(rlatu(j))  
           END DO  
           zdx(1, j) = zdx(iip1, j)  
        END DO  
        DO j = 2, jjm  
           DO i = 1, iip1  
              zdy(i, j) = 0.5*(cv_2d(i, j-1)+cv_2d(i, j))  
           END DO  
        END DO  
        DO i = 1, iip1  
           zdx(i, 1) = zdx(i, 2)  
           zdx(i, jjp1) = zdx(i, jjm)  
           zdy(i, 1) = zdy(i, 2)  
           zdy(i, jjp1) = zdy(i, jjm)  
        END DO  
        DO j = 1, jjp1  
           DO i = 1, iip1  
              dxdys(i, j) = sqrt(zdx(i, j)*zdx(i, j)+zdy(i, j)*zdy(i, j))  
           END DO  
        END DO  
        DO j = 1, jjp1  
           DO i = 1, iim  
              dxdyu(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))  
           END DO  
           dxdyu(iip1, j) = dxdyu(1, j)  
        END DO  
        DO j = 1, jjm  
           DO i = 1, iip1  
              dxdyv(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))  
           END DO  
        END DO  
24    
25         CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL  ')      PRINT *, 'Call sequence information: tau2alpha'
        CALL dump2d(iip1, jjp1, dxdyu, 'DX2DY2 U     ')  
        CALL dump2d(iip1, jjp1, dxdyv, 'DX2DY2 v     ')  
   
        !   coordonnees du centre du zoom  
        CALL coordij(clon, clat, ilon, ilat)  
        !   aire de la maille au centre du zoom  
        dxdy_min = dxdys(ilon, ilat)  
        !   dxdy maximale de la maille  
        dxdy_max = 0.  
        DO j = 1, jjp1  
           DO i = 1, iip1  
              dxdy_max = max(dxdy_max, dxdys(i, j))  
           END DO  
        END DO  
26    
27         IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN      n_lon = assert_eq(size(alpha, 1), size(dxdy, 1), "tau2alpha n_lon")
28            PRINT *, 'ATTENTION modele peu zoome'      n_lat = assert_eq(size(alpha, 2), size(dxdy, 2), size(rlat), &
29            PRINT *, 'ATTENTION on prend une constante de guidage cste'           "tau2alpha n_lat")
30            gamma = 0.  
31        alphamin = factt / taumax
32        alphamax = factt / taumin
33    
34        DO j = 1, n_lat
35           zlat = rlat(j) * 180. / pi
36           IF (lat_min_guide <= zlat .AND. zlat <= lat_max_guide) THEN
37              DO i = 1, n_lon
38                 xi = min(((dxdy_max - dxdy(i, j)) &
39                      / (dxdy_max - dxdy_min))**gamma, 1.)
40                 alpha(i, j) = xi * alphamin + (1. - xi) * alphamax
41              END DO
42         ELSE         ELSE
43            gamma = (dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)            alpha(:, j) = 0.
           PRINT *, 'gamma=', gamma  
           IF (gamma<1.E-5) THEN  
              PRINT *, 'gamma =', gamma, '<1e-5'  
              STOP  
           END IF  
           PRINT *, 'gamma=', gamma  
           gamma = log(0.5)/log(gamma)  
44         END IF         END IF
     END IF  
   
     alphamin = factt/taumax  
     alphamax = factt/taumin  
   
     DO j = 1, pjm  
        DO i = 1, pim  
           IF (type==1) THEN  
              dxdy_ = dxdys(i, j)  
              zlat = rlatu(j)*180./pi  
           ELSE IF (type==2) THEN  
              dxdy_ = dxdyu(i, j)  
              zlat = rlatu(j)*180./pi  
           ELSE IF (type==3) THEN  
              dxdy_ = dxdyv(i, j)  
              zlat = rlatv(j)*180./pi  
           END IF  
           IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN  
              !  pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin  
              alpha(i, j) = alphamin  
           ELSE  
              xi = ((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma  
              xi = min(xi, 1.)  
              IF (lat_min_guide<=zlat .AND. zlat<=lat_max_guide) THEN  
                 alpha(i, j) = xi*alphamin + (1.-xi)*alphamax  
              ELSE  
                 alpha(i, j) = 0.  
              END IF  
           END IF  
        END DO  
45      END DO      END DO
46    
47    END SUBROUTINE tau2alpha    END SUBROUTINE tau2alpha
Legend:
Removed from v.37  
changed lines
  Added in v.115
  ViewVC Help
Powered by ViewVC 1.1.21