1 | !-------------------------------------------------------------------------- |
---|
2 | !---------------------------- wind_centered ---------------------------------- |
---|
3 | ! Perot reconstruction based on Gauss theorem |
---|
4 | ! u = sum( u.edge_normal * edge_length * (edge_midpoint-cell_centroid) ) /cell_area |
---|
5 | !$OMP DO SCHEDULE(STATIC) |
---|
6 | DO ij = 1, primal_num |
---|
7 | ! this VLOOP iterates over primal cell edges |
---|
8 | SELECT CASE(primal_deg(ij)) |
---|
9 | CASE(4) |
---|
10 | edge1 = primal_edge(1,ij) |
---|
11 | edge2 = primal_edge(2,ij) |
---|
12 | edge3 = primal_edge(3,ij) |
---|
13 | edge4 = primal_edge(4,ij) |
---|
14 | sign1 = primal_ne(1,ij) |
---|
15 | sign2 = primal_ne(2,ij) |
---|
16 | sign3 = primal_ne(3,ij) |
---|
17 | sign4 = primal_ne(4,ij) |
---|
18 | ij_up1 = up(edge1) |
---|
19 | ij_up2 = up(edge2) |
---|
20 | ij_up3 = up(edge3) |
---|
21 | ij_up4 = up(edge4) |
---|
22 | ij_down1 = down(edge1) |
---|
23 | ij_down2 = down(edge2) |
---|
24 | ij_down3 = down(edge3) |
---|
25 | ij_down4 = down(edge4) |
---|
26 | !DIR$ SIMD |
---|
27 | DO l = 1, llm |
---|
28 | cx=centroid(ij,1) |
---|
29 | cy=centroid(ij,2) |
---|
30 | cz=centroid(ij,3) |
---|
31 | ux=0. ; uy=0. ; uz=0. |
---|
32 | ue_le = sign1*ue(l,edge1)*le(edge1) |
---|
33 | ux = ux + ue_le*(.5*(xyz_v(ij_up1,1)+xyz_v(ij_down1,1))-cx) |
---|
34 | uy = uy + ue_le*(.5*(xyz_v(ij_up1,2)+xyz_v(ij_down1,2))-cy) |
---|
35 | uz = uz + ue_le*(.5*(xyz_v(ij_up1,3)+xyz_v(ij_down1,3))-cz) |
---|
36 | ue_le = sign2*ue(l,edge2)*le(edge2) |
---|
37 | ux = ux + ue_le*(.5*(xyz_v(ij_up2,1)+xyz_v(ij_down2,1))-cx) |
---|
38 | uy = uy + ue_le*(.5*(xyz_v(ij_up2,2)+xyz_v(ij_down2,2))-cy) |
---|
39 | uz = uz + ue_le*(.5*(xyz_v(ij_up2,3)+xyz_v(ij_down2,3))-cz) |
---|
40 | ue_le = sign3*ue(l,edge3)*le(edge3) |
---|
41 | ux = ux + ue_le*(.5*(xyz_v(ij_up3,1)+xyz_v(ij_down3,1))-cx) |
---|
42 | uy = uy + ue_le*(.5*(xyz_v(ij_up3,2)+xyz_v(ij_down3,2))-cy) |
---|
43 | uz = uz + ue_le*(.5*(xyz_v(ij_up3,3)+xyz_v(ij_down3,3))-cz) |
---|
44 | ue_le = sign4*ue(l,edge4)*le(edge4) |
---|
45 | ux = ux + ue_le*(.5*(xyz_v(ij_up4,1)+xyz_v(ij_down4,1))-cx) |
---|
46 | uy = uy + ue_le*(.5*(xyz_v(ij_up4,2)+xyz_v(ij_down4,2))-cy) |
---|
47 | uz = uz + ue_le*(.5*(xyz_v(ij_up4,3)+xyz_v(ij_down4,3))-cz) |
---|
48 | fac = scale*(1./Ai(ij)) |
---|
49 | ucenter(l,ij,1)=ux*fac |
---|
50 | ucenter(l,ij,2)=uy*fac |
---|
51 | ucenter(l,ij,3)=uz*fac |
---|
52 | END DO |
---|
53 | CASE(5) |
---|
54 | edge1 = primal_edge(1,ij) |
---|
55 | edge2 = primal_edge(2,ij) |
---|
56 | edge3 = primal_edge(3,ij) |
---|
57 | edge4 = primal_edge(4,ij) |
---|
58 | edge5 = primal_edge(5,ij) |
---|
59 | sign1 = primal_ne(1,ij) |
---|
60 | sign2 = primal_ne(2,ij) |
---|
61 | sign3 = primal_ne(3,ij) |
---|
62 | sign4 = primal_ne(4,ij) |
---|
63 | sign5 = primal_ne(5,ij) |
---|
64 | ij_up1 = up(edge1) |
---|
65 | ij_up2 = up(edge2) |
---|
66 | ij_up3 = up(edge3) |
---|
67 | ij_up4 = up(edge4) |
---|
68 | ij_up5 = up(edge5) |
---|
69 | ij_down1 = down(edge1) |
---|
70 | ij_down2 = down(edge2) |
---|
71 | ij_down3 = down(edge3) |
---|
72 | ij_down4 = down(edge4) |
---|
73 | ij_down5 = down(edge5) |
---|
74 | !DIR$ SIMD |
---|
75 | DO l = 1, llm |
---|
76 | cx=centroid(ij,1) |
---|
77 | cy=centroid(ij,2) |
---|
78 | cz=centroid(ij,3) |
---|
79 | ux=0. ; uy=0. ; uz=0. |
---|
80 | ue_le = sign1*ue(l,edge1)*le(edge1) |
---|
81 | ux = ux + ue_le*(.5*(xyz_v(ij_up1,1)+xyz_v(ij_down1,1))-cx) |
---|
82 | uy = uy + ue_le*(.5*(xyz_v(ij_up1,2)+xyz_v(ij_down1,2))-cy) |
---|
83 | uz = uz + ue_le*(.5*(xyz_v(ij_up1,3)+xyz_v(ij_down1,3))-cz) |
---|
84 | ue_le = sign2*ue(l,edge2)*le(edge2) |
---|
85 | ux = ux + ue_le*(.5*(xyz_v(ij_up2,1)+xyz_v(ij_down2,1))-cx) |
---|
86 | uy = uy + ue_le*(.5*(xyz_v(ij_up2,2)+xyz_v(ij_down2,2))-cy) |
---|
87 | uz = uz + ue_le*(.5*(xyz_v(ij_up2,3)+xyz_v(ij_down2,3))-cz) |
---|
88 | ue_le = sign3*ue(l,edge3)*le(edge3) |
---|
89 | ux = ux + ue_le*(.5*(xyz_v(ij_up3,1)+xyz_v(ij_down3,1))-cx) |
---|
90 | uy = uy + ue_le*(.5*(xyz_v(ij_up3,2)+xyz_v(ij_down3,2))-cy) |
---|
91 | uz = uz + ue_le*(.5*(xyz_v(ij_up3,3)+xyz_v(ij_down3,3))-cz) |
---|
92 | ue_le = sign4*ue(l,edge4)*le(edge4) |
---|
93 | ux = ux + ue_le*(.5*(xyz_v(ij_up4,1)+xyz_v(ij_down4,1))-cx) |
---|
94 | uy = uy + ue_le*(.5*(xyz_v(ij_up4,2)+xyz_v(ij_down4,2))-cy) |
---|
95 | uz = uz + ue_le*(.5*(xyz_v(ij_up4,3)+xyz_v(ij_down4,3))-cz) |
---|
96 | ue_le = sign5*ue(l,edge5)*le(edge5) |
---|
97 | ux = ux + ue_le*(.5*(xyz_v(ij_up5,1)+xyz_v(ij_down5,1))-cx) |
---|
98 | uy = uy + ue_le*(.5*(xyz_v(ij_up5,2)+xyz_v(ij_down5,2))-cy) |
---|
99 | uz = uz + ue_le*(.5*(xyz_v(ij_up5,3)+xyz_v(ij_down5,3))-cz) |
---|
100 | fac = scale*(1./Ai(ij)) |
---|
101 | ucenter(l,ij,1)=ux*fac |
---|
102 | ucenter(l,ij,2)=uy*fac |
---|
103 | ucenter(l,ij,3)=uz*fac |
---|
104 | END DO |
---|
105 | CASE(6) |
---|
106 | edge1 = primal_edge(1,ij) |
---|
107 | edge2 = primal_edge(2,ij) |
---|
108 | edge3 = primal_edge(3,ij) |
---|
109 | edge4 = primal_edge(4,ij) |
---|
110 | edge5 = primal_edge(5,ij) |
---|
111 | edge6 = primal_edge(6,ij) |
---|
112 | sign1 = primal_ne(1,ij) |
---|
113 | sign2 = primal_ne(2,ij) |
---|
114 | sign3 = primal_ne(3,ij) |
---|
115 | sign4 = primal_ne(4,ij) |
---|
116 | sign5 = primal_ne(5,ij) |
---|
117 | sign6 = primal_ne(6,ij) |
---|
118 | ij_up1 = up(edge1) |
---|
119 | ij_up2 = up(edge2) |
---|
120 | ij_up3 = up(edge3) |
---|
121 | ij_up4 = up(edge4) |
---|
122 | ij_up5 = up(edge5) |
---|
123 | ij_up6 = up(edge6) |
---|
124 | ij_down1 = down(edge1) |
---|
125 | ij_down2 = down(edge2) |
---|
126 | ij_down3 = down(edge3) |
---|
127 | ij_down4 = down(edge4) |
---|
128 | ij_down5 = down(edge5) |
---|
129 | ij_down6 = down(edge6) |
---|
130 | !DIR$ SIMD |
---|
131 | DO l = 1, llm |
---|
132 | cx=centroid(ij,1) |
---|
133 | cy=centroid(ij,2) |
---|
134 | cz=centroid(ij,3) |
---|
135 | ux=0. ; uy=0. ; uz=0. |
---|
136 | ue_le = sign1*ue(l,edge1)*le(edge1) |
---|
137 | ux = ux + ue_le*(.5*(xyz_v(ij_up1,1)+xyz_v(ij_down1,1))-cx) |
---|
138 | uy = uy + ue_le*(.5*(xyz_v(ij_up1,2)+xyz_v(ij_down1,2))-cy) |
---|
139 | uz = uz + ue_le*(.5*(xyz_v(ij_up1,3)+xyz_v(ij_down1,3))-cz) |
---|
140 | ue_le = sign2*ue(l,edge2)*le(edge2) |
---|
141 | ux = ux + ue_le*(.5*(xyz_v(ij_up2,1)+xyz_v(ij_down2,1))-cx) |
---|
142 | uy = uy + ue_le*(.5*(xyz_v(ij_up2,2)+xyz_v(ij_down2,2))-cy) |
---|
143 | uz = uz + ue_le*(.5*(xyz_v(ij_up2,3)+xyz_v(ij_down2,3))-cz) |
---|
144 | ue_le = sign3*ue(l,edge3)*le(edge3) |
---|
145 | ux = ux + ue_le*(.5*(xyz_v(ij_up3,1)+xyz_v(ij_down3,1))-cx) |
---|
146 | uy = uy + ue_le*(.5*(xyz_v(ij_up3,2)+xyz_v(ij_down3,2))-cy) |
---|
147 | uz = uz + ue_le*(.5*(xyz_v(ij_up3,3)+xyz_v(ij_down3,3))-cz) |
---|
148 | ue_le = sign4*ue(l,edge4)*le(edge4) |
---|
149 | ux = ux + ue_le*(.5*(xyz_v(ij_up4,1)+xyz_v(ij_down4,1))-cx) |
---|
150 | uy = uy + ue_le*(.5*(xyz_v(ij_up4,2)+xyz_v(ij_down4,2))-cy) |
---|
151 | uz = uz + ue_le*(.5*(xyz_v(ij_up4,3)+xyz_v(ij_down4,3))-cz) |
---|
152 | ue_le = sign5*ue(l,edge5)*le(edge5) |
---|
153 | ux = ux + ue_le*(.5*(xyz_v(ij_up5,1)+xyz_v(ij_down5,1))-cx) |
---|
154 | uy = uy + ue_le*(.5*(xyz_v(ij_up5,2)+xyz_v(ij_down5,2))-cy) |
---|
155 | uz = uz + ue_le*(.5*(xyz_v(ij_up5,3)+xyz_v(ij_down5,3))-cz) |
---|
156 | ue_le = sign6*ue(l,edge6)*le(edge6) |
---|
157 | ux = ux + ue_le*(.5*(xyz_v(ij_up6,1)+xyz_v(ij_down6,1))-cx) |
---|
158 | uy = uy + ue_le*(.5*(xyz_v(ij_up6,2)+xyz_v(ij_down6,2))-cy) |
---|
159 | uz = uz + ue_le*(.5*(xyz_v(ij_up6,3)+xyz_v(ij_down6,3))-cz) |
---|
160 | fac = scale*(1./Ai(ij)) |
---|
161 | ucenter(l,ij,1)=ux*fac |
---|
162 | ucenter(l,ij,2)=uy*fac |
---|
163 | ucenter(l,ij,3)=uz*fac |
---|
164 | END DO |
---|
165 | CASE DEFAULT |
---|
166 | !DIR$ SIMD |
---|
167 | DO l = 1, llm |
---|
168 | cx=centroid(ij,1) |
---|
169 | cy=centroid(ij,2) |
---|
170 | cz=centroid(ij,3) |
---|
171 | ux=0. ; uy=0. ; uz=0. |
---|
172 | DO iedge = 1, primal_deg(ij) |
---|
173 | edge = primal_edge(iedge,ij) |
---|
174 | ij_up = up(edge) |
---|
175 | ij_down = down(edge) |
---|
176 | ue_le = primal_ne(iedge,ij)*ue(l,edge)*le(edge) |
---|
177 | ux = ux + ue_le*(.5*(xyz_v(ij_up,1)+xyz_v(ij_down,1))-cx) |
---|
178 | uy = uy + ue_le*(.5*(xyz_v(ij_up,2)+xyz_v(ij_down,2))-cy) |
---|
179 | uz = uz + ue_le*(.5*(xyz_v(ij_up,3)+xyz_v(ij_down,3))-cz) |
---|
180 | END DO |
---|
181 | fac = scale*(1./Ai(ij)) |
---|
182 | ucenter(l,ij,1)=ux*fac |
---|
183 | ucenter(l,ij,2)=uy*fac |
---|
184 | ucenter(l,ij,3)=uz*fac |
---|
185 | END DO |
---|
186 | END SELECT |
---|
187 | END DO |
---|
188 | !$OMP END DO |
---|
189 | !---------------------------- wind_centered ---------------------------------- |
---|
190 | !-------------------------------------------------------------------------- |
---|