/[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.f revision 103 by guez, Fri Aug 29 13:00:05 2014 UTC trunk/Sources/dyn3d/Guide/tau2alpha.f revision 211 by guez, Tue Dec 13 17:23:09 2016 UTC
# Line 1  Line 1 
1  module tau2alpha_m  module tau2alpha_m
2    
   USE paramet_m, ONLY : iip1, jjp1  
   USE dimens_m, ONLY : jjm  
   
3    IMPLICIT NONE    IMPLICIT NONE
4    
   private iip1, jjp1, jjm  
   
   REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)  
   
5  contains  contains
6    
7    SUBROUTINE tau2alpha(type, factt, taumin, taumax, alpha)    SUBROUTINE tau2alpha(lat_min_guide, lat_max_guide, factt, dxdy, rlat, &
8           taumin, taumax, alpha)
9    
10      USE comgeom, ONLY : cu_2d, cv_2d, rlatu, rlatv      use init_tau2alpha_m, only: dxdy_min, dxdy_max, gamma
11      use conf_guide_m, only: lat_min_guide, lat_max_guide      USE nr_util, ONLY: assert_eq
12      USE dimens_m, ONLY : iim  
13      use dump2d_m, only: dump2d      ! Dans le cas où on n'a les analyses que sur une bande de latitudes :
14      USE nr_util, ONLY : pi      REAL, intent(in):: lat_min_guide ! minimum latitude for nudging, in rad
15      USE serre, ONLY : clat, clon, grossismx, grossismy      real, intent(in):: lat_max_guide ! maximum latitude for nudging, in rad
16    
17      INTEGER, intent(in):: type      REAL, intent(in):: factt
18      REAL, intent(in):: factt, taumin, taumax      ! pas de temps entre deux appels au guidage, en jours
19      real, intent(out):: alpha(:, :)      
20        REAL, intent(in):: dxdy(:, :) ! (n_lon, n_lat)
21        REAL, intent(in):: rlat(:) ! (n_lat)
22        REAL, intent(in):: taumin, taumax
23        real, intent(out):: alpha(:, :) ! (n_lon, n_lat)
24    
25      ! Local:      ! Local:
26      REAL dxdy      REAL a_min, a_max, xi
27      REAL dxdy_min, dxdy_max      INTEGER i, j, n_lon, n_lat
     REAL alphamin, alphamax, xi  
     REAL, SAVE:: gamma  
     INTEGER i, j, ilon, ilat  
     LOGICAL:: first = .TRUE.  
     REAL zdx(iip1, jjp1), zdy(iip1, jjp1)  
     REAL zlat  
28    
29      !------------------------------------------------------------      !------------------------------------------------------------
30    
31      IF (first) THEN      PRINT *, 'Call sequence information: tau2alpha'
        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)**2 + zdy(i, j)**2)  
           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, j + 1))  
           END DO  
        END DO  
32    
33         CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL ')      n_lon = assert_eq(size(alpha, 1), size(dxdy, 1), "tau2alpha n_lon")
34         CALL dump2d(iip1, jjp1, dxdyu, 'DX2DY2 U ')      n_lat = assert_eq(size(alpha, 2), size(dxdy, 2), size(rlat), &
35         CALL dump2d(iip1, jjp1, dxdyv, 'DX2DY2 v ')           "tau2alpha n_lat")
36    
37         ! coordonnees du centre du zoom      a_min = factt / taumax
38         CALL coordij(clon, clat, ilon, ilat)      a_max = factt / taumin
39         ! aire de la maille au centre du zoom  
40         dxdy_min = dxdys(ilon, ilat)      DO j = 1, n_lat
41           IF (lat_min_guide <= rlat(j) .AND. rlat(j) <= lat_max_guide) THEN
42         ! dxdy maximal de la maille :            DO i = 1, n_lon
43         dxdy_max = 0.               xi = min(((dxdy_max - dxdy(i, j)) &
44         DO j = 1, jjp1                    / (dxdy_max - dxdy_min))**gamma, 1.)
45            DO i = 1, iip1               alpha(i, j) = 1. - exp(- xi * a_min - (1. - xi) * a_max)
              dxdy_max = max(dxdy_max, dxdys(i, j))  
46            END DO            END DO
        END DO  
   
        IF (abs(grossismx - 1.)<0.1 .OR. abs(grossismy - 1.)<0.1) THEN  
           PRINT *, 'ATTENTION modele peu zoome'  
           PRINT *, 'ATTENTION on prend une constante de guidage cste'  
           gamma = 0.  
47         ELSE         ELSE
48            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)  
49         END IF         END IF
        first = .false.  
     END IF  
   
     alphamin = factt / taumax  
     alphamax = factt / taumin  
   
     DO j = 1, size(alpha, 2)  
        DO i = 1, size(alpha, 1)  
           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  
50      END DO      END DO
51    
52    END SUBROUTINE tau2alpha    END SUBROUTINE tau2alpha
Legend:
Removed from v.103  
changed lines
  Added in v.211
  ViewVC Help
Powered by ViewVC 1.1.21