source: codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_hydro.F90 @ 849

Last change on this file since 849 was 849, checked in by jisesh, 5 years ago

devel: separate module for compute_caldyn_slow_hydro

File size: 4.5 KB
Line 
1MODULE compute_caldyn_slow_hydro_mod
2  USE grid_param, ONLY : llm
3  IMPLICIT NONE
4  PRIVATE
5
6  PUBLIC :: compute_caldyn_slow_hydro
7
8CONTAINS
9
10  SUBROUTINE compute_caldyn_slow_hydro(u,rhodz,hv, hflux,Kv,du, zero)
11    USE icosa
12    USE omp_para, ONLY : ll_begin, ll_end
13    USE caldyn_vars_mod
14    LOGICAL, INTENT(IN) :: zero
15    REAL(rstd),INTENT(IN)  :: u(3*iim*jjm,llm)    ! prognostic "velocity"
16    REAL(rstd),INTENT(IN)  :: Kv(2*iim*jjm,llm)   ! kinetic energy at vertices
17    REAL(rstd),INTENT(IN)  :: hv(2*iim*jjm,llm)   ! height/mass averaged to vertices
18    REAL(rstd),INTENT(IN)  :: rhodz(iim*jjm,llm)
19    REAL(rstd),INTENT(OUT) :: hflux(3*iim*jjm,llm) ! hflux in kg/s
20    REAL(rstd),INTENT(INOUT) :: du(3*iim*jjm,llm)
21   
22    REAL(rstd) :: berni(iim*jjm,llm)  ! Bernoulli function
23    REAL(rstd) :: berni1(iim*jjm)  ! Bernoulli function
24    REAL(rstd) :: uu_right, uu_lup, uu_ldown, ke, uu
25    INTEGER :: ij,l
26
27    CALL trace_start("compute_caldyn_slow_hydro")
28
29    IF(dysl_slow_hydro) THEN
30
31#define BERNI(ij,l) berni(ij,l)
32#include "../kernels_hex/caldyn_slow_hydro.k90"
33#undef BERNI
34
35     ELSE
36
37#define BERNI(ij) berni1(ij)
38
39    DO l = ll_begin, ll_end
40       !  Compute mass fluxes
41       IF (caldyn_conserv==conserv_energy) CALL test_message(req_qu) 
42
43       IF(caldyn_kinetic==kinetic_trisk) THEN
44          !DIR$ SIMD
45          DO ij=ij_begin_ext,ij_end_ext
46             uu_right=0.5*(rhodz(ij,l)+rhodz(ij+t_right,l))*u(ij+u_right,l)
47             uu_lup=0.5*(rhodz(ij,l)+rhodz(ij+t_lup,l))*u(ij+u_lup,l)
48             uu_ldown=0.5*(rhodz(ij,l)+rhodz(ij+t_ldown,l))*u(ij+u_ldown,l)
49             uu_right= uu_right*le_de(ij+u_right)
50             uu_lup  = uu_lup  *le_de(ij+u_lup)
51             uu_ldown= uu_ldown*le_de(ij+u_ldown)
52             hflux(ij+u_right,l)=uu_right
53             hflux(ij+u_lup,l)  =uu_lup
54             hflux(ij+u_ldown,l)=uu_ldown
55          ENDDO
56       ELSE ! mass flux deriving from consistent kinetic energy
57          !DIR$ SIMD
58          DO ij=ij_begin_ext,ij_end_ext
59             uu_right=0.5*(hv(ij+z_rup,l)+hv(ij+z_rdown,l))*u(ij+u_right,l)
60             uu_lup=0.5*(hv(ij+z_up,l)+hv(ij+z_lup,l))*u(ij+u_lup,l)
61             uu_ldown=0.5*(hv(ij+z_ldown,l)+hv(ij+z_down,l))*u(ij+u_ldown,l)
62             uu_right= uu_right*le_de(ij+u_right)
63             uu_lup  = uu_lup  *le_de(ij+u_lup)
64             uu_ldown= uu_ldown*le_de(ij+u_ldown)
65             hflux(ij+u_right,l)=uu_right
66             hflux(ij+u_lup,l)  =uu_lup
67             hflux(ij+u_ldown,l)=uu_ldown
68          ENDDO
69       END IF
70
71       ! Compute Bernoulli=kinetic energy
72       IF(caldyn_kinetic==kinetic_trisk) THEN
73          !DIR$ SIMD
74          DO ij=ij_begin,ij_end         
75             BERNI(ij) = &
76                  1/(4*Ai(ij))*(le_de(ij+u_right)*u(ij+u_right,l)**2 +    &
77                                le_de(ij+u_rup)*u(ij+u_rup,l)**2     +    &
78                                le_de(ij+u_lup)*u(ij+u_lup,l)**2     +    &
79                                le_de(ij+u_left)*u(ij+u_left,l)**2   +    &
80                                le_de(ij+u_ldown)*u(ij+u_ldown,l)**2 +    &
81                                le_de(ij+u_rdown)*u(ij+u_rdown,l)**2 ) 
82          ENDDO
83       ELSE
84          !DIR$ SIMD
85          DO ij=ij_begin,ij_end
86             BERNI(ij) = Riv(ij,vup)   *Kv(ij+z_up,l)    + &
87                         Riv(ij,vlup)  *Kv(ij+z_lup,l)   + &
88                         Riv(ij,vldown)*Kv(ij+z_ldown,l) + &
89                         Riv(ij,vdown) *Kv(ij+z_down,l)  + &
90                         Riv(ij,vrdown)*Kv(ij+z_rdown,l) + &
91                         Riv(ij,vrup)  *Kv(ij+z_rup,l)
92          END DO
93       END IF
94       ! Compute du=-grad(Bernoulli)
95       IF(zero) THEN
96          !DIR$ SIMD
97          DO ij=ij_begin,ij_end
98             du(ij+u_right,l) = ne_right*(BERNI(ij)-BERNI(ij+t_right))
99             du(ij+u_lup,l)   = ne_lup*(BERNI(ij)-BERNI(ij+t_lup))
100             du(ij+u_ldown,l) = ne_ldown*(BERNI(ij)-BERNI(ij+t_ldown))
101          END DO
102       ELSE
103          !DIR$ SIMD
104          DO ij=ij_begin,ij_end
105             du(ij+u_right,l) = du(ij+u_right,l) + &
106                  ne_right*(BERNI(ij)-BERNI(ij+t_right))
107             du(ij+u_lup,l)   = du(ij+u_lup,l) + &
108                  ne_lup*(BERNI(ij)-BERNI(ij+t_lup))
109             du(ij+u_ldown,l) = du(ij+u_ldown,l) + &
110                  ne_ldown*(BERNI(ij)-BERNI(ij+t_ldown))
111          END DO
112       END IF
113    END DO
114
115#undef BERNI
116
117    END IF ! dysl
118    CALL trace_end("compute_caldyn_slow_hydro")   
119  END SUBROUTINE compute_caldyn_slow_hydro
120
121END MODULE compute_caldyn_slow_hydro_mod
Note: See TracBrowser for help on using the repository browser.