source: codes/icosagcm/devel/src/kernels_hex/energy_fluxes.k90 @ 724

Last change on this file since 724 was 724, checked in by dubos, 6 years ago

devel : updated generated kernels ; added missing SIMD directives
File size: 7.5 KB
Line 
1   !--------------------------------------------------------------------------
2   !---------------------------- energy_fluxes ----------------------------------
3   ! First diagnose geopotential and temperature, column-wise
4   !$OMP BARRIER
5   !DIR$ SIMD
6   DO ij=ij_omp_begin_ext,ij_omp_end_ext
7      pk(ij,llm) = ptop + .5*g*rhodz(ij,llm)
8   END DO
9   DO l = llm-1,1,-1
10      !DIR$ SIMD
11      DO ij=ij_omp_begin_ext,ij_omp_end_ext
12         pk(ij,l) = pk(ij,l+1) + (.5*g)*( rhodz(ij,l)+rhodz(ij,l+1) )
13      END DO
14   END DO
15   ! NB : at this point pressure is stored in array pk
16   ! pk then serves as buffer to store temperature
17   SELECT CASE(caldyn_thermo)
18   CASE(thermo_theta)
19      DO l = 1,llm
20         !DIR$ SIMD
21         DO ij=ij_omp_begin_ext,ij_omp_end_ext
22            p_ik = pk(ij,l)
23            theta_ik = theta_rhodz(ij,l,1)/rhodz(ij,l)
24            theta(ij,l) = theta_ik
25            temp_ik = theta_ik*(p_ik/preff)**kappa
26            gv = (g*Rd)*temp_ik/p_ik
27            pk(ij,l) = temp_ik
28            geopot(ij,l+1) = geopot(ij,l) + gv*rhodz(ij,l)
29         END DO
30      END DO
31   CASE(thermo_entropy)
32      DO l = 1,llm
33         !DIR$ SIMD
34         DO ij=ij_omp_begin_ext,ij_omp_end_ext
35            p_ik = pk(ij,l)
36            theta_ik = theta_rhodz(ij,l,1)/rhodz(ij,l)
37            temp_ik = Treff*exp((theta_ik + Rd*log(p_ik/preff))/cpp)
38            theta(ij,l) = Treff*exp(theta_ik/cpp)
39            gv = (g*Rd)*temp_ik/p_ik ! specific volume v = Rd*T/p
40            pk(ij,l) = temp_ik
41            geopot(ij,l+1) = geopot(ij,l) + gv*rhodz(ij,l)
42         END DO
43      END DO
44   END SELECT
45   !$OMP BARRIER
46   ! Now accumulate energies and energy fluxes
47   ! NB : at this point temperature is stored in array pk
48   ! pk then serves as buffer to store energy
49   ! enthalpy
50   DO l = ll_begin, ll_end
51      !DIR$ SIMD
52      DO ij=ij_begin_ext, ij_end_ext
53         energy = cpp*pk(ij,l)
54         enthalpy(ij,l) = enthalpy(ij,l) + frac*rhodz(ij,l)*energy
55         pk(ij,l) = energy
56      END DO
57   END DO
58   DO l = ll_begin, ll_end
59      !DIR$ SIMD
60      DO ij=ij_begin_ext, ij_end_ext
61         enthalpy_flux(ij+u_right,l) = enthalpy_flux(ij+u_right,l) + .5*massflux(ij+u_right,l)*(pk(ij,l)+pk(ij+t_right,l))
62         enthalpy_flux(ij+u_lup,l) = enthalpy_flux(ij+u_lup,l) + .5*massflux(ij+u_lup,l)*(pk(ij,l)+pk(ij+t_lup,l))
63         enthalpy_flux(ij+u_ldown,l) = enthalpy_flux(ij+u_ldown,l) + .5*massflux(ij+u_ldown,l)*(pk(ij,l)+pk(ij+t_ldown,l))
64      END DO
65   END DO
66   ! potential energy
67   DO l = ll_begin, ll_end
68      !DIR$ SIMD
69      DO ij=ij_begin_ext, ij_end_ext
70         energy = .5*(geopot(ij,l+1)+geopot(ij,l))
71         epot(ij,l) = epot(ij,l) + frac*rhodz(ij,l)*energy
72         pk(ij,l) = energy
73      END DO
74   END DO
75   DO l = ll_begin, ll_end
76      !DIR$ SIMD
77      DO ij=ij_begin_ext, ij_end_ext
78         epot_flux(ij+u_right,l) = epot_flux(ij+u_right,l) + .5*massflux(ij+u_right,l)*(pk(ij,l)+pk(ij+t_right,l))
79         epot_flux(ij+u_lup,l) = epot_flux(ij+u_lup,l) + .5*massflux(ij+u_lup,l)*(pk(ij,l)+pk(ij+t_lup,l))
80         epot_flux(ij+u_ldown,l) = epot_flux(ij+u_ldown,l) + .5*massflux(ij+u_ldown,l)*(pk(ij,l)+pk(ij+t_ldown,l))
81      END DO
82   END DO
83   ! theta
84   DO l = ll_begin, ll_end
85      !DIR$ SIMD
86      DO ij=ij_begin_ext, ij_end_ext
87         energy = theta(ij,l)
88         thetat(ij,l) = thetat(ij,l) + frac*rhodz(ij,l)*energy
89         pk(ij,l) = energy
90      END DO
91   END DO
92   DO l = ll_begin, ll_end
93      !DIR$ SIMD
94      DO ij=ij_begin_ext, ij_end_ext
95         thetat_flux(ij+u_right,l) = thetat_flux(ij+u_right,l) + .5*massflux(ij+u_right,l)*(pk(ij,l)+pk(ij+t_right,l))
96         thetat_flux(ij+u_lup,l) = thetat_flux(ij+u_lup,l) + .5*massflux(ij+u_lup,l)*(pk(ij,l)+pk(ij+t_lup,l))
97         thetat_flux(ij+u_ldown,l) = thetat_flux(ij+u_ldown,l) + .5*massflux(ij+u_ldown,l)*(pk(ij,l)+pk(ij+t_ldown,l))
98      END DO
99   END DO
100   ! kinetic energy
101   DO l = ll_begin, ll_end
102      !DIR$ SIMD
103      DO ij=ij_begin_ext, ij_end_ext
104         energy=0.d0
105         energy = energy + le(ij+u_rup)*de(ij+u_rup)*ue(ij+u_rup,l)**2
106         energy = energy + le(ij+u_lup)*de(ij+u_lup)*ue(ij+u_lup,l)**2
107         energy = energy + le(ij+u_left)*de(ij+u_left)*ue(ij+u_left,l)**2
108         energy = energy + le(ij+u_ldown)*de(ij+u_ldown)*ue(ij+u_ldown,l)**2
109         energy = energy + le(ij+u_rdown)*de(ij+u_rdown)*ue(ij+u_rdown,l)**2
110         energy = energy + le(ij+u_right)*de(ij+u_right)*ue(ij+u_right,l)**2
111         energy = energy * (.25/Ai(ij))
112         ekin(ij,l) = ekin(ij,l) + frac*rhodz(ij,l)*energy
113         pk(ij,l) = energy
114      END DO
115   END DO
116   DO l = ll_begin, ll_end
117      !DIR$ SIMD
118      DO ij=ij_begin_ext, ij_end_ext
119         ekin_flux(ij+u_right,l) = ekin_flux(ij+u_right,l) + .5*massflux(ij+u_right,l)*(pk(ij,l)+pk(ij+t_right,l))
120         ekin_flux(ij+u_lup,l) = ekin_flux(ij+u_lup,l) + .5*massflux(ij+u_lup,l)*(pk(ij,l)+pk(ij+t_lup,l))
121         ekin_flux(ij+u_ldown,l) = ekin_flux(ij+u_ldown,l) + .5*massflux(ij+u_ldown,l)*(pk(ij,l)+pk(ij+t_ldown,l))
122      END DO
123   END DO
124   ! ulon
125   DO l = ll_begin, ll_end
126      !DIR$ SIMD
127      DO ij=ij_begin_ext, ij_end_ext
128         cx=centroid(ij,1)
129         cy=centroid(ij,2)
130         cz=centroid(ij,3)
131         ux=0. ; uy=0. ; uz=0.
132         ue_le = ne_rup*ue(ij+u_rup,l)*le(ij+u_rup)
133         ux = ux + ue_le*(.5*(xyz_v(ij+z_rup,1)+xyz_v(ij+z_up,1))-cx)
134         uy = uy + ue_le*(.5*(xyz_v(ij+z_rup,2)+xyz_v(ij+z_up,2))-cy)
135         uz = uz + ue_le*(.5*(xyz_v(ij+z_rup,3)+xyz_v(ij+z_up,3))-cz)
136         ue_le = ne_lup*ue(ij+u_lup,l)*le(ij+u_lup)
137         ux = ux + ue_le*(.5*(xyz_v(ij+z_lup,1)+xyz_v(ij+z_up,1))-cx)
138         uy = uy + ue_le*(.5*(xyz_v(ij+z_lup,2)+xyz_v(ij+z_up,2))-cy)
139         uz = uz + ue_le*(.5*(xyz_v(ij+z_lup,3)+xyz_v(ij+z_up,3))-cz)
140         ue_le = ne_left*ue(ij+u_left,l)*le(ij+u_left)
141         ux = ux + ue_le*(.5*(xyz_v(ij+z_lup,1)+xyz_v(ij+z_ldown,1))-cx)
142         uy = uy + ue_le*(.5*(xyz_v(ij+z_lup,2)+xyz_v(ij+z_ldown,2))-cy)
143         uz = uz + ue_le*(.5*(xyz_v(ij+z_lup,3)+xyz_v(ij+z_ldown,3))-cz)
144         ue_le = ne_ldown*ue(ij+u_ldown,l)*le(ij+u_ldown)
145         ux = ux + ue_le*(.5*(xyz_v(ij+z_ldown,1)+xyz_v(ij+z_down,1))-cx)
146         uy = uy + ue_le*(.5*(xyz_v(ij+z_ldown,2)+xyz_v(ij+z_down,2))-cy)
147         uz = uz + ue_le*(.5*(xyz_v(ij+z_ldown,3)+xyz_v(ij+z_down,3))-cz)
148         ue_le = ne_rdown*ue(ij+u_rdown,l)*le(ij+u_rdown)
149         ux = ux + ue_le*(.5*(xyz_v(ij+z_rdown,1)+xyz_v(ij+z_down,1))-cx)
150         uy = uy + ue_le*(.5*(xyz_v(ij+z_rdown,2)+xyz_v(ij+z_down,2))-cy)
151         uz = uz + ue_le*(.5*(xyz_v(ij+z_rdown,3)+xyz_v(ij+z_down,3))-cz)
152         ue_le = ne_right*ue(ij+u_right,l)*le(ij+u_right)
153         ux = ux + ue_le*(.5*(xyz_v(ij+z_rup,1)+xyz_v(ij+z_rdown,1))-cx)
154         uy = uy + ue_le*(.5*(xyz_v(ij+z_rup,2)+xyz_v(ij+z_rdown,2))-cy)
155         uz = uz + ue_le*(.5*(xyz_v(ij+z_rup,3)+xyz_v(ij+z_rdown,3))-cz)
156         ulon_i = ux*elon_i(ij,1) + uy*elon_i(ij,2) + uz*elon_i(ij,3)
157         energy = ulon_i*(1./Ai(ij))
158         ulon(ij,l) = ulon(ij,l) + frac*rhodz(ij,l)*energy
159         pk(ij,l) = energy
160      END DO
161   END DO
162   DO l = ll_begin, ll_end
163      !DIR$ SIMD
164      DO ij=ij_begin_ext, ij_end_ext
165         ulon_flux(ij+u_right,l) = ulon_flux(ij+u_right,l) + .5*massflux(ij+u_right,l)*(pk(ij,l)+pk(ij+t_right,l))
166         ulon_flux(ij+u_lup,l) = ulon_flux(ij+u_lup,l) + .5*massflux(ij+u_lup,l)*(pk(ij,l)+pk(ij+t_lup,l))
167         ulon_flux(ij+u_ldown,l) = ulon_flux(ij+u_ldown,l) + .5*massflux(ij+u_ldown,l)*(pk(ij,l)+pk(ij+t_ldown,l))
168      END DO
169   END DO
170   !---------------------------- energy_fluxes ----------------------------------
171   !--------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.