source: codes/icosagcm/devel/src/kernels_unst/pvort_only.k90 @ 913

Last change on this file since 913 was 658, checked in by dubos, 7 years ago

devel/unstructured : updated kernels

File size: 2.9 KB
Line 
1   !--------------------------------------------------------------------------
2   !---------------------------- pvort_only ----------------------------------
3   !$OMP DO SCHEDULE(STATIC)
4   DO ij = 1, dual_num
5      ! this VLOOP iterates over dual cell edges
6      SELECT CASE(dual_deg(ij))
7      CASE(3)
8         edge1 = dual_edge(1,ij)
9         edge2 = dual_edge(2,ij)
10         edge3 = dual_edge(3,ij)
11         sign1 = dual_ne(1,ij)
12         sign2 = dual_ne(2,ij)
13         sign3 = dual_ne(3,ij)
14         vertex1 = dual_vertex(1,ij)
15         vertex2 = dual_vertex(2,ij)
16         vertex3 = dual_vertex(3,ij)
17         !DIR$ SIMD
18         DO l = 1, llm
19            etav = 0.d0
20            etav = etav + sign1*u(l,edge1)
21            etav = etav + sign2*u(l,edge2)
22            etav = etav + sign3*u(l,edge3)
23            hv=0.
24            hv = hv + Riv2(1,ij)*rhodz(l,vertex1)
25            hv = hv + Riv2(2,ij)*rhodz(l,vertex2)
26            hv = hv + Riv2(3,ij)*rhodz(l,vertex3)
27            qv(l,ij) = (etav + fv(ij)*Av(ij) )/(hv*Av(ij))
28         END DO
29      CASE(4)
30         edge1 = dual_edge(1,ij)
31         edge2 = dual_edge(2,ij)
32         edge3 = dual_edge(3,ij)
33         edge4 = dual_edge(4,ij)
34         sign1 = dual_ne(1,ij)
35         sign2 = dual_ne(2,ij)
36         sign3 = dual_ne(3,ij)
37         sign4 = dual_ne(4,ij)
38         vertex1 = dual_vertex(1,ij)
39         vertex2 = dual_vertex(2,ij)
40         vertex3 = dual_vertex(3,ij)
41         vertex4 = dual_vertex(4,ij)
42         !DIR$ SIMD
43         DO l = 1, llm
44            etav = 0.d0
45            etav = etav + sign1*u(l,edge1)
46            etav = etav + sign2*u(l,edge2)
47            etav = etav + sign3*u(l,edge3)
48            etav = etav + sign4*u(l,edge4)
49            hv=0.
50            hv = hv + Riv2(1,ij)*rhodz(l,vertex1)
51            hv = hv + Riv2(2,ij)*rhodz(l,vertex2)
52            hv = hv + Riv2(3,ij)*rhodz(l,vertex3)
53            hv = hv + Riv2(4,ij)*rhodz(l,vertex4)
54            qv(l,ij) = (etav + fv(ij)*Av(ij) )/(hv*Av(ij))
55         END DO
56      CASE DEFAULT
57         !DIR$ SIMD
58         DO l = 1, llm
59            etav = 0.d0
60            DO iedge = 1, dual_deg(ij)
61               edge = dual_edge(iedge,ij)
62               etav = etav + dual_ne(iedge,ij)*u(l,edge)
63            END DO
64            hv=0.
65            DO ivertex = 1, dual_deg(ij)
66               vertex = dual_vertex(ivertex,ij)
67               hv = hv + Riv2(ivertex,ij)*rhodz(l,vertex)
68            END DO
69            qv(l,ij) = (etav + fv(ij)*Av(ij) )/(hv*Av(ij))
70         END DO
71      END SELECT
72   END DO
73   !$OMP END DO
74   !$OMP DO SCHEDULE(STATIC)
75   DO edge = 1, edge_num
76      ij_up = up(edge)
77      ij_down = down(edge)
78      !DIR$ SIMD
79      DO l = 1, llm
80         qu(l,edge)=0.5d0*(qv(l,ij_down)+qv(l,ij_up))
81      END DO
82   END DO
83   !$OMP END DO
84   !---------------------------- pvort_only ----------------------------------
85   !--------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.