/[lmdze]/trunk/dyn3d/Dissipation/dissip.f
ViewVC logotype

Diff of /trunk/dyn3d/Dissipation/dissip.f

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

trunk/libf/dyn3d/dissip.f90 revision 26 by guez, Tue Mar 9 15:27:15 2010 UTC trunk/dyn3d/Dissipation/dissip.f revision 254 by guez, Mon Feb 5 10:39:38 2018 UTC
# Line 1  Line 1 
1  SUBROUTINE dissip(vcov,ucov,teta,p,dv,du,dh)  module dissip_m
   
   ! From dyn3d/dissip.F,v 1.1.1.1 2004/05/19 12:53:05  
   ! Avec nouveaux operateurs star :  gradiv2 , divgrad2, nxgraro2  ...  
   ! Auteur:  P. Le Van                                                    
   ! Objet: dissipation horizontale                                              
   
   USE dimens_m, ONLY : llm  
   USE paramet_m, ONLY : iip1, iip2, ip1jm, ip1jmp1, llmp1  
   USE comdissnew, ONLY : lstardis, nitergdiv, nitergrot, niterh  
   USE inidissip_m, ONLY : dtdiss, tetah, tetaudiv, tetaurot  
2    
3    IMPLICIT NONE    IMPLICIT NONE
4    
5    !   Arguments:                                                            contains
6    REAL :: vcov(ip1jm,llm), ucov(ip1jmp1,llm), teta(ip1jmp1,llm)  
7    REAL, INTENT (IN) :: p(ip1jmp1,llmp1)    SUBROUTINE dissip(vcov, ucov, teta, p, dv, du, dh)
8    REAL :: dv(ip1jm,llm), du(ip1jmp1,llm), dh(ip1jmp1,llm)  
9        ! From dyn3d/dissip.F, version 1.1.1.1 2004/05/19 12:53:05
10    !   Local:                                                                    ! Author: P. Le Van
11    REAL :: gdx(ip1jmp1,llm), gdy(ip1jm,llm)      ! Objet : calcul de la dissipation horizontale
12    REAL :: grx(ip1jmp1,llm), gry(ip1jm,llm)      ! Avec opérateurs star : gradiv2, divgrad2, nxgraro2
13    REAL :: te1dt(llm), te2dt(llm), te3dt(llm)  
14    REAL :: deltapres(ip1jmp1,llm)      USE comdissnew, ONLY: nitergdiv, nitergrot, niterh
15        USE dimens_m, ONLY: iim, jjm, llm
16    INTEGER :: l, ij      use divgrad2_m, only: divgrad2
17        use gradiv2_m, only: gradiv2
18    !-----------------------------------------------------------------------      USE inidissip_m, ONLY: dtdiss, tetah, tetaudiv, tetaurot, cdivu, crot, cdivh
19        use nr_util, only: assert
20    !   initialisations:                                                          use nxgraro2_m, only: nxgraro2
21    
22    DO l = 1, llm      REAL, intent(in):: vcov(:, :, :) ! (iim + 1, jjm, llm)
23       te1dt(l) = tetaudiv(l)*dtdiss      REAL, intent(in):: ucov(:, :, :) ! (iim + 1, jjm + 1, llm)
24       te2dt(l) = tetaurot(l)*dtdiss      REAL, intent(in):: teta(:, :, :) ! (iim + 1, jjm + 1, llm)
25       te3dt(l) = tetah(l)*dtdiss      REAL, INTENT(IN):: p(:, :, :) ! (iim + 1, jjm + 1, llm + 1)
26    END DO      REAL, intent(out):: dv(:, :, :) ! (iim + 1, jjm, llm)
27    du = 0.      REAL, intent(out):: du(:, :, :) ! (iim + 1, jjm + 1, llm)
28    dv = 0.      REAL, intent(out):: dh(:, :, :) ! (iim + 1, jjm + 1, llm)
29    dh = 0.  
30        ! Local:
31    !   Calcul de la dissipation:                                                REAL gdx(iim + 1, jjm + 1, llm), gdy(iim + 1, jjm, llm)
32        REAL tedt(llm)
33    !   Calcul de la partie   grad  ( div ) :                                    REAL deltapres(iim + 1, jjm + 1, llm)
34        INTEGER l
35    IF (lstardis) THEN  
36       CALL gradiv2(llm,ucov,vcov,nitergdiv,gdx,gdy)      !-----------------------------------------------------------------------
37    ELSE  
38       CALL gradiv(llm,ucov,vcov,nitergdiv,gdx,gdy)      call assert((/size(vcov, 1), size(ucov, 1), size(teta, 1), size(p, 1), &
39    END IF           size(dv, 1), size(du, 1), size(dh, 1)/) == iim + 1, "dissip iim")
40        call assert((/size(vcov, 2), size(ucov, 2) - 1, size(teta, 2) - 1, &
41    DO l = 1, llm           size(p, 2) - 1, size(dv, 2), size(du, 2) - 1, size(dh, 2) - 1/) &
42             == jjm, "dissip jjm")
43       DO ij = 1, iip1      call assert((/size(vcov, 3), size(ucov, 3), size(teta, 3), size(p, 3) - 1, &
44          gdx(ij,l) = 0.           size(dv, 3), size(du, 3), size(dh, 3)/) == llm, "dissip llm")
45          gdx(ij+ip1jm,l) = 0.  
46       END DO      du(:, 1, :) = 0.
47        du(:, jjm + 1, :) = 0.
48       DO ij = iip2, ip1jm  
49          du(ij,l) = du(ij,l) - te1dt(l)*gdx(ij,l)      ! Calcul de la partie grad(div) :
50       END DO      CALL gradiv2(ucov, vcov, nitergdiv, gdx, gdy, cdivu)
51       DO ij = 1, ip1jm      tedt = tetaudiv * dtdiss
52          dv(ij,l) = dv(ij,l) - te1dt(l)*gdy(ij,l)      forall (l = 1: llm)
53       END DO         du(:, 2: jjm, l) = - tedt(l) * gdx(:, 2: jjm, l)
54    END DO         dv(:, :, l) = - tedt(l) * gdy(:, :, l)
55        END forall
56    !   calcul de la partie   n X grad ( rot ):                              
57        ! Calcul de la partie n X grad(rot) :
58    IF (lstardis) THEN      CALL nxgraro2(ucov, vcov, nitergrot, gdx, gdy, crot)
59       CALL nxgraro2(llm,ucov,vcov,nitergrot,grx,gry)      tedt = tetaurot * dtdiss
60    ELSE      forall (l = 1: llm)
61       CALL nxgrarot(llm,ucov,vcov,nitergrot,grx,gry)         du(:, 2: jjm, l) = du(:, 2: jjm, l) - tedt(l) * gdx(:, 2: jjm, l)
62    END IF         dv(:, :, l) = dv(:, :, l) - tedt(l) * gdy(:, :, l)
63        END forall
64    
65    DO l = 1, llm      ! calcul de la partie div(grad) :
66       DO ij = 1, iip1      forall (l = 1: llm) &
67          grx(ij,l) = 0.           deltapres(:, :, l) = max(0., p(:, :, l) - p(:, :, l + 1))
68       END DO      CALL divgrad2(llm, teta, deltapres, niterh, gdx, cdivh)
69        forall (l = 1: llm) dh(:, :, l) = - tetah(l) * dtdiss * gdx(:, :, l)
70       DO ij = iip2, ip1jm  
71          du(ij,l) = du(ij,l) - te2dt(l)*grx(ij,l)    END SUBROUTINE dissip
      END DO  
      DO ij = 1, ip1jm  
         dv(ij,l) = dv(ij,l) - te2dt(l)*gry(ij,l)  
      END DO  
   END DO  
   
   !   calcul de la partie   div ( grad ):                                  
   
   IF (lstardis) THEN  
   
      DO l = 1, llm  
         DO ij = 1, ip1jmp1  
            deltapres(ij,l) = amax1(0.,p(ij,l)-p(ij,l+1))  
         END DO  
      END DO  
   
      CALL divgrad2(llm,teta,deltapres,niterh,gdx)  
   ELSE  
      CALL divgrad(llm,teta,niterh,gdx)  
   END IF  
   
   DO l = 1, llm  
      DO ij = 1, ip1jmp1  
         dh(ij,l) = dh(ij,l) - te3dt(l)*gdx(ij,l)  
      END DO  
   END DO  
72    
73  END SUBROUTINE dissip  end module dissip_m
Legend:
Removed from v.26  
changed lines
  Added in v.254
  ViewVC Help
Powered by ViewVC 1.1.21