trunk/dyn3d/tau2alpha.f

module tau2alpha_m

    USE paramet_m, ONLY : iip1, jjp1
    USE dimens_m, ONLY : jjm

    IMPLICIT NONE

    private iip1, jjp1, jjm

    REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)

contains

  SUBROUTINE tau2alpha(type, pim, pjm, factt, taumin, taumax, alpha)

    use conf_guide_m, only: lat_min_guide, lat_max_guide
    USE dimens_m, ONLY : iim
    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
    REAL, intent(in):: 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)
    REAL zlat

    !------------------------------------------------------------

    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

       CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL  ')
       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

       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.
       ELSE
          gamma = (dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)
          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)
       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
    END DO

  END SUBROUTINE tau2alpha

end module tau2alpha_m
1	guez	37	module tau2alpha_m
2
3	guez	44	USE paramet_m, ONLY : iip1, jjp1
4			USE dimens_m, ONLY : jjm
5
6	guez	37	IMPLICIT NONE
7
8	guez	44	private iip1, jjp1, jjm
9
10			REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)
11	guez	37
12			contains
13
14			SUBROUTINE tau2alpha(type, pim, pjm, factt, taumin, taumax, alpha)
15
16	guez	83	use conf_guide_m, only: lat_min_guide, lat_max_guide
17	guez	44	USE dimens_m, ONLY : iim
18	guez	39	USE nr_util, ONLY : pi
19	guez	37	USE comgeom, ONLY : cu_2d, cv_2d, rlatu, rlatv
20			USE serre, ONLY : clat, clon, grossismx, grossismy
21
22			! arguments :
23			INTEGER type
24			INTEGER pim, pjm
25	guez	44	REAL, intent(in):: factt, taumin, taumax
26	guez	37	REAL dxdy_, alpha(pim, pjm)
27			REAL dxdy_min, dxdy_max
28
29			! local :
30			REAL alphamin, alphamax, gamma, xi
31			SAVE gamma
32			INTEGER i, j, ilon, ilat
33
34			LOGICAL first
35			SAVE first
36			DATA first/ .TRUE./
37
38			REAL zdx(iip1, jjp1), zdy(iip1, jjp1)
39			REAL zlat
40
41	guez	44	!------------------------------------------------------------
42
43	guez	37	IF (first) THEN
44			DO j = 2, jjm
45			DO i = 2, iip1
46			zdx(i, j) = 0.5*(cu_2d(i-1, j)+cu_2d(i, j))/cos(rlatu(j))
47			END DO
48			zdx(1, j) = zdx(iip1, j)
49			END DO
50			DO j = 2, jjm
51			DO i = 1, iip1
52			zdy(i, j) = 0.5*(cv_2d(i, j-1)+cv_2d(i, j))
53			END DO
54			END DO
55			DO i = 1, iip1
56			zdx(i, 1) = zdx(i, 2)
57			zdx(i, jjp1) = zdx(i, jjm)
58			zdy(i, 1) = zdy(i, 2)
59			zdy(i, jjp1) = zdy(i, jjm)
60			END DO
61			DO j = 1, jjp1
62			DO i = 1, iip1
63			dxdys(i, j) = sqrt(zdx(i, j)zdx(i, j)+zdy(i, j)zdy(i, j))
64			END DO
65			END DO
66			DO j = 1, jjp1
67			DO i = 1, iim
68			dxdyu(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))
69			END DO
70			dxdyu(iip1, j) = dxdyu(1, j)
71			END DO
72			DO j = 1, jjm
73			DO i = 1, iip1
74			dxdyv(i, j) = 0.5*(dxdys(i, j)+dxdys(i+1, j))
75			END DO
76			END DO
77
78			CALL dump2d(iip1, jjp1, dxdys, 'DX2DY2 SCAL ')
79			CALL dump2d(iip1, jjp1, dxdyu, 'DX2DY2 U ')
80			CALL dump2d(iip1, jjp1, dxdyv, 'DX2DY2 v ')
81
82			! coordonnees du centre du zoom
83			CALL coordij(clon, clat, ilon, ilat)
84			! aire de la maille au centre du zoom
85			dxdy_min = dxdys(ilon, ilat)
86			! dxdy maximale de la maille
87			dxdy_max = 0.
88			DO j = 1, jjp1
89			DO i = 1, iip1
90			dxdy_max = max(dxdy_max, dxdys(i, j))
91			END DO
92			END DO
93
94			IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN
95			PRINT *, 'ATTENTION modele peu zoome'
96			PRINT *, 'ATTENTION on prend une constante de guidage cste'
97			gamma = 0.
98			ELSE
99			gamma = (dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)
100			PRINT *, 'gamma=', gamma
101			IF (gamma<1.E-5) THEN
102			PRINT *, 'gamma =', gamma, '<1e-5'
103			STOP
104			END IF
105			PRINT *, 'gamma=', gamma
106			gamma = log(0.5)/log(gamma)
107			END IF
108			END IF
109
110			alphamin = factt/taumax
111			alphamax = factt/taumin
112
113			DO j = 1, pjm
114			DO i = 1, pim
115			IF (type==1) THEN
116			dxdy_ = dxdys(i, j)
117			zlat = rlatu(j)*180./pi
118			ELSE IF (type==2) THEN
119			dxdy_ = dxdyu(i, j)
120			zlat = rlatu(j)*180./pi
121			ELSE IF (type==3) THEN
122			dxdy_ = dxdyv(i, j)
123			zlat = rlatv(j)*180./pi
124			END IF
125			IF (abs(grossismx-1.)<0.1 .OR. abs(grossismy-1.)<0.1) THEN
126			! pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin
127			alpha(i, j) = alphamin
128			ELSE
129			xi = ((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma
130			xi = min(xi, 1.)
131			IF (lat_min_guide<=zlat .AND. zlat<=lat_max_guide) THEN
132			alpha(i, j) = xialphamin + (1.-xi)alphamax
133			ELSE
134			alpha(i, j) = 0.
135			END IF
136			END IF
137			END DO
138			END DO
139
140			END SUBROUTINE tau2alpha
141
142			end module tau2alpha_m