Nektar++
Loading...
Searching...
No Matches
TestHexCollection.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: TestHexCollection.cpp
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20/// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description:
32//
33///////////////////////////////////////////////////////////////////////////////
34
35#include <Collections/Collection.h>
36#include <Collections/CollectionOptimisation.h>
37#include <LocalRegions/HexExp.h>
38#include <SpatialDomains/MeshGraph.h>
39#include <boost/test/tools/floating_point_comparison.hpp>
40#include <boost/test/unit_test.hpp>
41namespace Nektar::HexCollectionTests
42{
43
44SpatialDomains::SegGeomUniquePtr CreateSegGeom(unsigned int id,
45 SpatialDomains::PointGeom *v0,
46 SpatialDomains::PointGeom *v1)
47{
48 std::array<SpatialDomains::PointGeom *, 2> vertices = {v0, v1};
49 SpatialDomains::SegGeomUniquePtr result(
50 new SpatialDomains::SegGeom(id, 3, vertices));
51 return result;
52}
53
54SpatialDomains::HexGeomUniquePtr CreateHex(
55 std::array<SpatialDomains::PointGeom *, 8> v,
56 std::array<SpatialDomains::SegGeomUniquePtr, 12> &segVec,
57 std::array<SpatialDomains::QuadGeomUniquePtr, 6> &faceVec)
58{
59 std::array<std::array<int, 2>, 12> edgeVerts = {{{{0, 1}},
60 {{1, 2}},
61 {{2, 3}},
62 {{3, 0}},
63 {{0, 4}},
64 {{1, 5}},
65 {{2, 6}},
66 {{3, 7}},
67 {{4, 5}},
68 {{5, 6}},
69 {{6, 7}},
70 {{7, 4}}}};
71 std::array<std::array<int, 4>, 6> faceEdges = {{{{0, 1, 2, 3}},
72 {{0, 5, 8, 4}},
73 {{1, 6, 9, 5}},
74 {{2, 6, 10, 7}},
75 {{3, 7, 11, 4}},
76 {{8, 9, 10, 11}}}};
77
78 // Create segments from vertices
79 for (int i = 0; i < 12; ++i)
80 {
81 segVec[i] = CreateSegGeom(i, v[edgeVerts[i][0]], v[edgeVerts[i][1]]);
82 }
83
84 // Create faces from edges
85 std::array<SpatialDomains::QuadGeom *, 6> faces;
86 for (int i = 0; i < 6; ++i)
87 {
88 std::array<SpatialDomains::SegGeom *, 4> face;
89 for (int j = 0; j < 4; ++j)
90 {
91 face[j] = segVec[faceEdges[i][j]].get();
92 }
93 faceVec[i] = SpatialDomains::QuadGeomUniquePtr(
94 new SpatialDomains::QuadGeom(i, face));
95 faces[i] = faceVec[i].get();
96 }
97
98 SpatialDomains::HexGeomUniquePtr hexGeom(
99 new SpatialDomains::HexGeom(0, faces));
100 return hexGeom;
101}
102
103BOOST_AUTO_TEST_CASE(TestHexBwdTrans_StdMat_UniformP)
104{
105 SpatialDomains::PointGeomUniquePtr v0(
106 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
107 SpatialDomains::PointGeomUniquePtr v1(
108 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
109 SpatialDomains::PointGeomUniquePtr v2(
110 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
111 SpatialDomains::PointGeomUniquePtr v3(
112 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
113 SpatialDomains::PointGeomUniquePtr v4(
114 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
115 SpatialDomains::PointGeomUniquePtr v5(
116 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
117 SpatialDomains::PointGeomUniquePtr v6(
118 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
119 SpatialDomains::PointGeomUniquePtr v7(
120 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
121
122 std::array<SpatialDomains::PointGeom *, 8> v = {
123 v0.get(), v1.get(), v2.get(), v3.get(),
124 v4.get(), v5.get(), v6.get(), v7.get()};
125 std::array<SpatialDomains::SegGeomUniquePtr, 12> segVec;
126 std::array<SpatialDomains::QuadGeomUniquePtr, 6> faceVec;
127 SpatialDomains::HexGeomUniquePtr hexGeom = CreateHex(v, segVec, faceVec);
128
129 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
130 Nektar::LibUtilities::eGaussLobattoLegendre;
131 Nektar::LibUtilities::BasisType basisTypeDir1 =
132 Nektar::LibUtilities::eModified_A;
133 unsigned int numQuadPoints = 6;
134 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
135 quadPointsTypeDir1);
136 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
137 quadPointsKeyDir1);
138
139 Nektar::LocalRegions::HexExpSharedPtr Exp =
140 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
141 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
142
143 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
144 CollExp.push_back(Exp);
145
146 LibUtilities::SessionReaderSharedPtr dummySession;
147 Collections::CollectionOptimisation colOpt(dummySession, 3,
148 Collections::eStdMat);
149 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(Exp);
150 Collections::Collection c(CollExp, impTypes);
151 c.Initialise(Collections::eBwdTrans);
152
153 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
154 for (int i = 0; i < coeffs.size(); ++i)
155 {
156 coeffs[i] = i + 1;
157 }
158 Array<OneD, NekDouble> phys1(Exp->GetTotPoints());
159 Array<OneD, NekDouble> phys2(Exp->GetTotPoints());
160
161 Exp->BwdTrans(coeffs, phys1);
162 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
163
164 double epsilon = 1.0e-8;
165 for (int i = 0; i < phys1.size(); ++i)
166 {
167 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
168 }
169}
170#if 0
171BOOST_AUTO_TEST_CASE(TestHexBwdTrans_StdMat_VariableP)
172{
173 SpatialDomains::PointGeomUniquePtr v0(
174 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
175 SpatialDomains::PointGeomUniquePtr v1(
176 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
177 SpatialDomains::PointGeomUniquePtr v2(
178 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
179 SpatialDomains::PointGeomUniquePtr v3(
180 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
181 SpatialDomains::PointGeomUniquePtr v4(
182 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
183 SpatialDomains::PointGeomUniquePtr v5(
184 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
185 SpatialDomains::PointGeomUniquePtr v6(
186 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
187 SpatialDomains::PointGeomUniquePtr v7(
188 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
189
190 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
191 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
192 SpatialDomains::HexGeomUniquePtr hexGeom =
193 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
194
195 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
196 Nektar::LibUtilities::eGaussLobattoLegendre;
197 Nektar::LibUtilities::BasisType basisTypeDir1 =
198 Nektar::LibUtilities::eModified_A;
199 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
200 quadPointsTypeDir1);
201 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
202 quadPointsTypeDir1);
203 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
204 quadPointsTypeDir1);
205 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
206 quadPointsKeyDir1);
207 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
208 quadPointsKeyDir2);
209 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
210 quadPointsKeyDir3);
211
212 Nektar::LocalRegions::HexExpSharedPtr Exp =
213 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
214 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
215
216 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
217 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
218 basisKeyDir1, basisKeyDir2, basisKeyDir3);
219
220 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
221 CollExp.push_back(Exp);
222
223 LibUtilities::SessionReaderSharedPtr dummySession;
224 Collections::CollectionOptimisation colOpt(dummySession, 3,
225 Collections::eStdMat);
226 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
227 Collections::Collection c(CollExp, impTypes);
228 c.Initialise(Collections::eBwdTrans);
229
230 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
231 for (int i = 0; i < coeffs.size(); ++i)
232 {
233 coeffs[i] = i + 1;
234 }
235 Array<OneD, NekDouble> phys1(Exp->GetTotPoints());
236 Array<OneD, NekDouble> phys2(Exp->GetTotPoints());
237
238 Exp->BwdTrans(coeffs, phys1);
239 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
240
241 double epsilon = 1.0e-8;
242 for (int i = 0; i < phys1.size(); ++i)
243 {
244 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
245 }
246}
247
248BOOST_AUTO_TEST_CASE(TestHexBwdTrans_IterPerExp_UniformP)
249{
250 SpatialDomains::PointGeomUniquePtr v0(
251 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
252 SpatialDomains::PointGeomUniquePtr v1(
253 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
254 SpatialDomains::PointGeomUniquePtr v2(
255 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
256 SpatialDomains::PointGeomUniquePtr v3(
257 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
258 SpatialDomains::PointGeomUniquePtr v4(
259 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
260 SpatialDomains::PointGeomUniquePtr v5(
261 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
262 SpatialDomains::PointGeomUniquePtr v6(
263 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
264 SpatialDomains::PointGeomUniquePtr v7(
265 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
266
267 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
268 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
269 SpatialDomains::HexGeomUniquePtr hexGeom =
270 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
271
272 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
273 Nektar::LibUtilities::eGaussLobattoLegendre;
274 Nektar::LibUtilities::BasisType basisTypeDir1 =
275 Nektar::LibUtilities::eModified_A;
276 unsigned int numQuadPoints = 6;
277 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
278 quadPointsTypeDir1);
279 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
280 quadPointsKeyDir1);
281
282 Nektar::LocalRegions::HexExpSharedPtr Exp =
283 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
284 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
285
286 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
287 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
288 basisKeyDir1, basisKeyDir1, basisKeyDir1);
289
290 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
291 CollExp.push_back(Exp);
292
293 LibUtilities::SessionReaderSharedPtr dummySession;
294 Collections::CollectionOptimisation colOpt(dummySession, 3,
295 Collections::eIterPerExp);
296 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
297 Collections::Collection c(CollExp, impTypes);
298 c.Initialise(Collections::eBwdTrans);
299
300 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
301 for (int i = 0; i < coeffs.size(); ++i)
302 {
303 coeffs[i] = i + 1;
304 }
305 Array<OneD, NekDouble> phys1(Exp->GetTotPoints());
306 Array<OneD, NekDouble> phys2(Exp->GetTotPoints());
307
308 Exp->BwdTrans(coeffs, phys1);
309 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
310
311 double epsilon = 1.0e-8;
312 for (int i = 0; i < phys1.size(); ++i)
313 {
314 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
315 }
316}
317
318BOOST_AUTO_TEST_CASE(TestHexBwdTrans_IterPerExp_VariableP)
319{
320 SpatialDomains::PointGeomUniquePtr v0(
321 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
322 SpatialDomains::PointGeomUniquePtr v1(
323 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
324 SpatialDomains::PointGeomUniquePtr v2(
325 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
326 SpatialDomains::PointGeomUniquePtr v3(
327 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
328 SpatialDomains::PointGeomUniquePtr v4(
329 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
330 SpatialDomains::PointGeomUniquePtr v5(
331 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
332 SpatialDomains::PointGeomUniquePtr v6(
333 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
334 SpatialDomains::PointGeomUniquePtr v7(
335 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
336
337 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
338 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
339 SpatialDomains::HexGeomUniquePtr hexGeom =
340 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
341
342 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
343 Nektar::LibUtilities::eGaussLobattoLegendre;
344 Nektar::LibUtilities::BasisType basisTypeDir1 =
345 Nektar::LibUtilities::eModified_A;
346 unsigned int numQuadPoints = 6;
347 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
348 quadPointsTypeDir1);
349 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
350 quadPointsKeyDir1);
351 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
352 quadPointsKeyDir1);
353 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
354 quadPointsKeyDir1);
355
356 Nektar::LocalRegions::HexExpSharedPtr Exp =
357 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
358 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
359
360 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
361 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
362 basisKeyDir1, basisKeyDir2, basisKeyDir3);
363
364 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
365 CollExp.push_back(Exp);
366
367 LibUtilities::SessionReaderSharedPtr dummySession;
368 Collections::CollectionOptimisation colOpt(dummySession, 3,
369 Collections::eIterPerExp);
370 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
371 Collections::Collection c(CollExp, impTypes);
372 c.Initialise(Collections::eBwdTrans);
373
374 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
375 for (int i = 0; i < coeffs.size(); ++i)
376 {
377 coeffs[i] = i + 1;
378 }
379 Array<OneD, NekDouble> phys1(Exp->GetTotPoints());
380 Array<OneD, NekDouble> phys2(Exp->GetTotPoints());
381
382 Exp->BwdTrans(coeffs, phys1);
383 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
384 c.Initialise(Collections::eBwdTrans);
385
386 double epsilon = 1.0e-8;
387 for (int i = 0; i < phys1.size(); ++i)
388 {
389 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
390 }
391}
392
393BOOST_AUTO_TEST_CASE(TestHexBwdTrans_IterPerExp_VariableP_MultiElmt)
394{
395 SpatialDomains::PointGeomUniquePtr v0(
396 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
397 SpatialDomains::PointGeomUniquePtr v1(
398 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
399 SpatialDomains::PointGeomUniquePtr v2(
400 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
401 SpatialDomains::PointGeomUniquePtr v3(
402 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
403 SpatialDomains::PointGeomUniquePtr v4(
404 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
405 SpatialDomains::PointGeomUniquePtr v5(
406 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
407 SpatialDomains::PointGeomUniquePtr v6(
408 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
409 SpatialDomains::PointGeomUniquePtr v7(
410 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
411
412 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
413 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
414 SpatialDomains::HexGeomUniquePtr hexGeom =
415 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
416
417 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
418 Nektar::LibUtilities::eGaussLobattoLegendre;
419 Nektar::LibUtilities::BasisType basisTypeDir1 =
420 Nektar::LibUtilities::eModified_A;
421 unsigned int numQuadPoints = 6;
422 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
423 quadPointsTypeDir1);
424 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
425 quadPointsKeyDir1);
426 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
427 quadPointsKeyDir1);
428 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
429 quadPointsKeyDir1);
430
431 Nektar::LocalRegions::HexExpSharedPtr Exp =
432 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
433 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
434
435 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
436 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
437 basisKeyDir1, basisKeyDir2, basisKeyDir3);
438
439 int nelmts = 10;
440
441 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
442 for (int i = 0; i < nelmts; ++i)
443 {
444 CollExp.push_back(Exp);
445 }
446
447 LibUtilities::SessionReaderSharedPtr dummySession;
448 Collections::CollectionOptimisation colOpt(dummySession, 3,
449 Collections::eIterPerExp);
450 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
451 Collections::Collection c(CollExp, impTypes);
452 c.Initialise(Collections::eBwdTrans);
453
454 Array<OneD, NekDouble> coeffs(nelmts * Exp->GetNcoeffs(), 1.0), tmp;
455 for (int i = 0; i < coeffs.size(); ++i)
456 {
457 coeffs[i] = i + 1;
458 }
459 Array<OneD, NekDouble> phys1(nelmts * Exp->GetTotPoints());
460 Array<OneD, NekDouble> phys2(nelmts * Exp->GetTotPoints());
461
462 for (int i = 0; i < nelmts; ++i)
463 {
464 Exp->BwdTrans(coeffs + i * Exp->GetNcoeffs(),
465 tmp = phys1 + i * Exp->GetTotPoints());
466 }
467
468 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
469
470 double epsilon = 1.0e-8;
471 for (int i = 0; i < phys1.size(); ++i)
472 {
473 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
474 }
475}
476
477BOOST_AUTO_TEST_CASE(TestHexBwdTrans_NoCollection_VariableP)
478{
479 SpatialDomains::PointGeomUniquePtr v0(
480 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
481 SpatialDomains::PointGeomUniquePtr v1(
482 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
483 SpatialDomains::PointGeomUniquePtr v2(
484 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
485 SpatialDomains::PointGeomUniquePtr v3(
486 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
487 SpatialDomains::PointGeomUniquePtr v4(
488 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
489 SpatialDomains::PointGeomUniquePtr v5(
490 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
491 SpatialDomains::PointGeomUniquePtr v6(
492 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
493 SpatialDomains::PointGeomUniquePtr v7(
494 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
495
496 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
497 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
498 SpatialDomains::HexGeomUniquePtr hexGeom =
499 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
500
501 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
502 Nektar::LibUtilities::eGaussLobattoLegendre;
503 Nektar::LibUtilities::BasisType basisTypeDir1 =
504 Nektar::LibUtilities::eModified_A;
505 unsigned int numQuadPoints = 6;
506 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
507 quadPointsTypeDir1);
508 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
509 quadPointsKeyDir1);
510 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
511 quadPointsKeyDir1);
512 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
513 quadPointsKeyDir1);
514
515 Nektar::LocalRegions::HexExpSharedPtr Exp =
516 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
517 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
518
519 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
520 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
521 basisKeyDir1, basisKeyDir2, basisKeyDir3);
522
523 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
524 CollExp.push_back(Exp);
525
526 LibUtilities::SessionReaderSharedPtr dummySession;
527 Collections::CollectionOptimisation colOpt(dummySession, 3,
528 Collections::eNoCollection);
529 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
530 Collections::Collection c(CollExp, impTypes);
531 c.Initialise(Collections::eBwdTrans);
532
533 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
534 for (int i = 0; i < coeffs.size(); ++i)
535 {
536 coeffs[i] = i + 1;
537 }
538 Array<OneD, NekDouble> phys1(Exp->GetTotPoints());
539 Array<OneD, NekDouble> phys2(Exp->GetTotPoints());
540
541 Exp->BwdTrans(coeffs, phys1);
542
543 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
544
545 double epsilon = 1.0e-8;
546 for (int i = 0; i < phys1.size(); ++i)
547 {
548 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
549 }
550}
551
552BOOST_AUTO_TEST_CASE(TestHexBwdTrans_SumFac_UniformP)
553{
554 SpatialDomains::PointGeomUniquePtr v0(
555 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
556 SpatialDomains::PointGeomUniquePtr v1(
557 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
558 SpatialDomains::PointGeomUniquePtr v2(
559 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
560 SpatialDomains::PointGeomUniquePtr v3(
561 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
562 SpatialDomains::PointGeomUniquePtr v4(
563 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
564 SpatialDomains::PointGeomUniquePtr v5(
565 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
566 SpatialDomains::PointGeomUniquePtr v6(
567 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
568 SpatialDomains::PointGeomUniquePtr v7(
569 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
570
571 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
572 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
573 SpatialDomains::HexGeomUniquePtr hexGeom =
574 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
575
576 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
577 Nektar::LibUtilities::eGaussLobattoLegendre;
578 Nektar::LibUtilities::BasisType basisTypeDir1 =
579 Nektar::LibUtilities::eModified_A;
580 unsigned int numQuadPoints = 6;
581 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
582 quadPointsTypeDir1);
583 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
584 quadPointsKeyDir1);
585
586 Nektar::LocalRegions::HexExpSharedPtr Exp =
587 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
588 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
589
590 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
591 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
592 basisKeyDir1, basisKeyDir1, basisKeyDir1);
593
594 int nelmts = 1;
595
596 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
597 for (int i = 0; i < nelmts; ++i)
598 {
599 CollExp.push_back(Exp);
600 }
601
602 LibUtilities::SessionReaderSharedPtr dummySession;
603 Collections::CollectionOptimisation colOpt(dummySession, 3,
604 Collections::eSumFac);
605 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
606 Collections::Collection c(CollExp, impTypes);
607 c.Initialise(Collections::eBwdTrans);
608
609 Array<OneD, NekDouble> coeffs(nelmts * Exp->GetNcoeffs(), 1.0), tmp;
610 for (int i = 0; i < coeffs.size(); ++i)
611 {
612 coeffs[i] = i + 1;
613 }
614 Array<OneD, NekDouble> phys1(nelmts * Exp->GetTotPoints());
615 Array<OneD, NekDouble> phys2(nelmts * Exp->GetTotPoints());
616
617 for (int i = 0; i < nelmts; ++i)
618 {
619 Exp->BwdTrans(coeffs + i * Exp->GetNcoeffs(),
620 tmp = phys1 + i * Exp->GetTotPoints());
621 }
622 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
623
624 double epsilon = 1.0e-8;
625 for (int i = 0; i < phys1.size(); ++i)
626 {
627 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
628 }
629}
630
631BOOST_AUTO_TEST_CASE(TestHexBwdTrans_SumFac_UniformP_MultiElmt)
632{
633 SpatialDomains::PointGeomUniquePtr v0(
634 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
635 SpatialDomains::PointGeomUniquePtr v1(
636 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
637 SpatialDomains::PointGeomUniquePtr v2(
638 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
639 SpatialDomains::PointGeomUniquePtr v3(
640 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
641 SpatialDomains::PointGeomUniquePtr v4(
642 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
643 SpatialDomains::PointGeomUniquePtr v5(
644 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
645 SpatialDomains::PointGeomUniquePtr v6(
646 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
647 SpatialDomains::PointGeomUniquePtr v7(
648 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
649
650 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
651 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
652 SpatialDomains::HexGeomUniquePtr hexGeom =
653 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
654
655 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
656 Nektar::LibUtilities::eGaussLobattoLegendre;
657 Nektar::LibUtilities::BasisType basisTypeDir1 =
658 Nektar::LibUtilities::eModified_A;
659 unsigned int numQuadPoints = 6;
660 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
661 quadPointsTypeDir1);
662 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
663 quadPointsKeyDir1);
664
665 Nektar::LocalRegions::HexExpSharedPtr Exp =
666 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
667 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
668
669 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
670 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
671 basisKeyDir1, basisKeyDir1, basisKeyDir1);
672
673 int nelmts = 10;
674
675 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
676 for (int i = 0; i < nelmts; ++i)
677 {
678 CollExp.push_back(Exp);
679 }
680
681 LibUtilities::SessionReaderSharedPtr dummySession;
682 Collections::CollectionOptimisation colOpt(dummySession, 3,
683 Collections::eSumFac);
684 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
685 Collections::Collection c(CollExp, impTypes);
686 c.Initialise(Collections::eBwdTrans);
687
688 Array<OneD, NekDouble> coeffs(nelmts * Exp->GetNcoeffs(), 1.0), tmp;
689 for (int i = 0; i < coeffs.size(); ++i)
690 {
691 coeffs[i] = i + 1;
692 }
693 Array<OneD, NekDouble> phys1(nelmts * Exp->GetTotPoints());
694 Array<OneD, NekDouble> phys2(nelmts * Exp->GetTotPoints());
695
696 for (int i = 0; i < nelmts; ++i)
697 {
698 Exp->BwdTrans(coeffs + i * Exp->GetNcoeffs(),
699 tmp = phys1 + i * Exp->GetTotPoints());
700 }
701 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
702
703 double epsilon = 1.0e-8;
704 for (int i = 0; i < phys1.size(); ++i)
705 {
706 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
707 }
708}
709
710BOOST_AUTO_TEST_CASE(TestHexBwdTrans_SumFac_VariableP)
711{
712 SpatialDomains::PointGeomUniquePtr v0(
713 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
714 SpatialDomains::PointGeomUniquePtr v1(
715 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
716 SpatialDomains::PointGeomUniquePtr v2(
717 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
718 SpatialDomains::PointGeomUniquePtr v3(
719 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
720 SpatialDomains::PointGeomUniquePtr v4(
721 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
722 SpatialDomains::PointGeomUniquePtr v5(
723 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
724 SpatialDomains::PointGeomUniquePtr v6(
725 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
726 SpatialDomains::PointGeomUniquePtr v7(
727 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
728
729 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
730 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
731 SpatialDomains::HexGeomUniquePtr hexGeom =
732 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
733
734 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
735 Nektar::LibUtilities::eGaussLobattoLegendre;
736 Nektar::LibUtilities::BasisType basisTypeDir1 =
737 Nektar::LibUtilities::eModified_A;
738 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
739 quadPointsTypeDir1);
740 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
741 quadPointsTypeDir1);
742 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
743 quadPointsTypeDir1);
744 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
745 quadPointsKeyDir1);
746 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
747 quadPointsKeyDir2);
748 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
749 quadPointsKeyDir3);
750
751 Nektar::LocalRegions::HexExpSharedPtr Exp =
752 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
753 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
754
755 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
756 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
757 basisKeyDir1, basisKeyDir2, basisKeyDir3);
758
759 int nelmts = 1;
760
761 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
762 for (int i = 0; i < nelmts; ++i)
763 {
764 CollExp.push_back(Exp);
765 }
766
767 LibUtilities::SessionReaderSharedPtr dummySession;
768 Collections::CollectionOptimisation colOpt(dummySession, 3,
769 Collections::eSumFac);
770 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
771 Collections::Collection c(CollExp, impTypes);
772 c.Initialise(Collections::eBwdTrans);
773
774 Array<OneD, NekDouble> coeffs(nelmts * Exp->GetNcoeffs(), 1.0), tmp;
775 for (int i = 0; i < coeffs.size(); ++i)
776 {
777 coeffs[i] = i + 1;
778 }
779 Array<OneD, NekDouble> phys1(nelmts * Exp->GetTotPoints());
780 Array<OneD, NekDouble> phys2(nelmts * Exp->GetTotPoints());
781
782 for (int i = 0; i < nelmts; ++i)
783 {
784 Exp->BwdTrans(coeffs + i * Exp->GetNcoeffs(),
785 tmp = phys1 + i * Exp->GetTotPoints());
786 }
787 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
788
789 double epsilon = 1.0e-8;
790 for (int i = 0; i < phys1.size(); ++i)
791 {
792 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
793 }
794}
795
796BOOST_AUTO_TEST_CASE(TestHexBwdTrans_SumFac_VariableP_MultiElmt)
797{
798 SpatialDomains::PointGeomUniquePtr v0(
799 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
800 SpatialDomains::PointGeomUniquePtr v1(
801 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
802 SpatialDomains::PointGeomUniquePtr v2(
803 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
804 SpatialDomains::PointGeomUniquePtr v3(
805 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
806 SpatialDomains::PointGeomUniquePtr v4(
807 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
808 SpatialDomains::PointGeomUniquePtr v5(
809 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
810 SpatialDomains::PointGeomUniquePtr v6(
811 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
812 SpatialDomains::PointGeomUniquePtr v7(
813 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
814
815 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
816 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
817 SpatialDomains::HexGeomUniquePtr hexGeom =
818 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
819
820 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
821 Nektar::LibUtilities::eGaussLobattoLegendre;
822 Nektar::LibUtilities::BasisType basisTypeDir1 =
823 Nektar::LibUtilities::eModified_A;
824 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
825 quadPointsTypeDir1);
826 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
827 quadPointsTypeDir1);
828 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
829 quadPointsTypeDir1);
830 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
831 quadPointsKeyDir1);
832 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
833 quadPointsKeyDir2);
834 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
835 quadPointsKeyDir3);
836
837 Nektar::LocalRegions::HexExpSharedPtr Exp =
838 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
839 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
840
841 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
842 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
843 basisKeyDir1, basisKeyDir2, basisKeyDir3);
844
845 int nelmts = 10;
846
847 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
848 for (int i = 0; i < nelmts; ++i)
849 {
850 CollExp.push_back(Exp);
851 }
852
853 LibUtilities::SessionReaderSharedPtr dummySession;
854 Collections::CollectionOptimisation colOpt(dummySession, 3,
855 Collections::eSumFac);
856 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
857 Collections::Collection c(CollExp, impTypes);
858 c.Initialise(Collections::eBwdTrans);
859
860 Array<OneD, NekDouble> coeffs(nelmts * Exp->GetNcoeffs(), 1.0), tmp;
861 for (int i = 0; i < coeffs.size(); ++i)
862 {
863 coeffs[i] = i + 1;
864 }
865 Array<OneD, NekDouble> phys1(nelmts * Exp->GetTotPoints());
866 Array<OneD, NekDouble> phys2(nelmts * Exp->GetTotPoints());
867
868 for (int i = 0; i < nelmts; ++i)
869 {
870 Exp->BwdTrans(coeffs + i * Exp->GetNcoeffs(),
871 tmp = phys1 + i * Exp->GetTotPoints());
872 }
873 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys2);
874
875 double epsilon = 1.0e-8;
876 for (int i = 0; i < phys1.size(); ++i)
877 {
878 BOOST_CHECK_CLOSE(phys1[i], phys2[i], epsilon);
879 }
880}
881
882BOOST_AUTO_TEST_CASE(TestHexBwdTrans_MatrixFree_UniformP)
883{
884 SpatialDomains::PointGeomUniquePtr v0(
885 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
886 SpatialDomains::PointGeomUniquePtr v1(
887 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
888 SpatialDomains::PointGeomUniquePtr v2(
889 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
890 SpatialDomains::PointGeomUniquePtr v3(
891 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
892 SpatialDomains::PointGeomUniquePtr v4(
893 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
894 SpatialDomains::PointGeomUniquePtr v5(
895 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
896 SpatialDomains::PointGeomUniquePtr v6(
897 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
898 SpatialDomains::PointGeomUniquePtr v7(
899 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
900
901 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
902 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
903 SpatialDomains::HexGeomUniquePtr hexGeom =
904 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
905
906 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
907 Nektar::LibUtilities::eGaussLobattoLegendre;
908 Nektar::LibUtilities::BasisType basisTypeDir1 =
909 Nektar::LibUtilities::eModified_A;
910 unsigned int numQuadPoints = 6;
911 unsigned int numModes = 4;
912 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
913 quadPointsTypeDir1);
914 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
915 quadPointsKeyDir1);
916
917 Nektar::LocalRegions::HexExpSharedPtr Exp =
918 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
919 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
920
921 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
922 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
923 basisKeyDir1, basisKeyDir1, basisKeyDir1);
924
925 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
926 CollExp.push_back(Exp);
927
928 LibUtilities::SessionReaderSharedPtr dummySession;
929 Collections::CollectionOptimisation colOpt(dummySession, 3,
930 Collections::eMatrixFree);
931 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
932 Collections::Collection c(CollExp, impTypes);
933 c.Initialise(Collections::eBwdTrans);
934
935 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs(), 1.0), tmp;
936 for (int i = 0; i < coeffs.size(); ++i)
937 {
938 coeffs[i] = i + 1;
939 }
940 Array<OneD, NekDouble> physRef(Exp->GetTotPoints());
941 Array<OneD, NekDouble> phys(Exp->GetTotPoints());
942
943 Exp->BwdTrans(coeffs, physRef);
944 c.ApplyOperator(Collections::eBwdTrans, coeffs, phys);
945
946 double epsilon = 1.0e-8;
947 for (int i = 0; i < physRef.size(); ++i)
948 {
949 BOOST_CHECK_CLOSE(physRef[i], phys[i], epsilon);
950 }
951}
952
953BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_StdMat_UniformP)
954{
955 SpatialDomains::PointGeomUniquePtr v0(
956 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
957 SpatialDomains::PointGeomUniquePtr v1(
958 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
959 SpatialDomains::PointGeomUniquePtr v2(
960 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
961 SpatialDomains::PointGeomUniquePtr v3(
962 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
963 SpatialDomains::PointGeomUniquePtr v4(
964 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
965 SpatialDomains::PointGeomUniquePtr v5(
966 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
967 SpatialDomains::PointGeomUniquePtr v6(
968 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
969 SpatialDomains::PointGeomUniquePtr v7(
970 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
971
972 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
973 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
974 SpatialDomains::HexGeomUniquePtr hexGeom =
975 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
976
977 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
978 Nektar::LibUtilities::eGaussLobattoLegendre;
979 Nektar::LibUtilities::BasisType basisTypeDir1 =
980 Nektar::LibUtilities::eModified_A;
981 unsigned int numQuadPoints = 6;
982 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
983 quadPointsTypeDir1);
984 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
985 quadPointsKeyDir1);
986
987 Nektar::LocalRegions::HexExpSharedPtr Exp =
988 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
989 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
990
991 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
992 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
993 basisKeyDir1, basisKeyDir1, basisKeyDir1);
994
995 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
996 CollExp.push_back(Exp);
997
998 LibUtilities::SessionReaderSharedPtr dummySession;
999 Collections::CollectionOptimisation colOpt(dummySession, 3,
1000 Collections::eStdMat);
1001 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1002 Collections::Collection c(CollExp, impTypes);
1003 c.Initialise(Collections::eIProductWRTBase);
1004
1005 const int nq = Exp->GetTotPoints();
1006 Array<OneD, NekDouble> phys(nq);
1007 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1008 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1009
1010 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1011
1012 Exp->GetCoords(xc, yc, zc);
1013
1014 for (int i = 0; i < nq; ++i)
1015 {
1016 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1017 }
1018
1019 Exp->IProductWRTBase(phys, coeffs1);
1020 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1021
1022 double epsilon = 1.0e-8;
1023 for (int i = 0; i < coeffs1.size(); ++i)
1024 {
1025 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1026 }
1027}
1028
1029BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_MatrixFree_UniformP_Undeformed)
1030{
1031 SpatialDomains::PointGeomUniquePtr v0(
1032 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
1033 SpatialDomains::PointGeomUniquePtr v1(
1034 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1035 SpatialDomains::PointGeomUniquePtr v2(
1036 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1037 SpatialDomains::PointGeomUniquePtr v3(
1038 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1039 SpatialDomains::PointGeomUniquePtr v4(
1040 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1041 SpatialDomains::PointGeomUniquePtr v5(
1042 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1043 SpatialDomains::PointGeomUniquePtr v6(
1044 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1045 SpatialDomains::PointGeomUniquePtr v7(
1046 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1047
1048 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1049 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1050 SpatialDomains::HexGeomUniquePtr hexGeom =
1051 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1052
1053 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1054 Nektar::LibUtilities::eGaussLobattoLegendre;
1055 Nektar::LibUtilities::BasisType basisTypeDir1 =
1056 Nektar::LibUtilities::eModified_A;
1057 unsigned int numQuadPoints = 5;
1058 unsigned int numModes = 4;
1059 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
1060 quadPointsTypeDir1);
1061 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
1062 quadPointsKeyDir1);
1063
1064 Nektar::LocalRegions::HexExpSharedPtr Exp =
1065 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1066 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
1067
1068 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1069 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1070 basisKeyDir1, basisKeyDir1, basisKeyDir1);
1071
1072 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1073 CollExp.push_back(Exp);
1074
1075 LibUtilities::SessionReaderSharedPtr dummySession;
1076 Collections::CollectionOptimisation colOpt(dummySession, 3,
1077 Collections::eMatrixFree);
1078 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1079 Collections::Collection c(CollExp, impTypes);
1080 c.Initialise(Collections::eIProductWRTBase);
1081
1082 const int nq = Exp->GetTotPoints();
1083 Array<OneD, NekDouble> phys(nq);
1084 Array<OneD, NekDouble> coeffsRef(Exp->GetNcoeffs());
1085 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs());
1086
1087 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1088
1089 Exp->GetCoords(xc, yc, zc);
1090
1091 for (int i = 0; i < nq; ++i)
1092 {
1093 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1094 }
1095
1096 Exp->IProductWRTBase(phys, coeffsRef);
1097 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs);
1098
1099 double epsilon = 1.0e-8;
1100 for (int i = 0; i < coeffsRef.size(); ++i)
1101 {
1102 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
1103 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
1104 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
1105 }
1106}
1107
1108BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_MatrixFree_UniformP_Deformed)
1109{
1110 SpatialDomains::PointGeomUniquePtr v0(
1111 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
1112 SpatialDomains::PointGeomUniquePtr v1(
1113 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1114 SpatialDomains::PointGeomUniquePtr v2(
1115 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1116 SpatialDomains::PointGeomUniquePtr v3(
1117 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1118 SpatialDomains::PointGeomUniquePtr v4(
1119 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1120 SpatialDomains::PointGeomUniquePtr v5(
1121 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1122 SpatialDomains::PointGeomUniquePtr v6(
1123 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
1124 SpatialDomains::PointGeomUniquePtr v7(
1125 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1126
1127 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1128 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1129 SpatialDomains::HexGeomUniquePtr hexGeom =
1130 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1131
1132 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1133 Nektar::LibUtilities::eGaussLobattoLegendre;
1134 Nektar::LibUtilities::BasisType basisTypeDir1 =
1135 Nektar::LibUtilities::eModified_A;
1136 unsigned int numQuadPoints = 5;
1137 unsigned int numModes = 4;
1138 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
1139 quadPointsTypeDir1);
1140 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
1141 quadPointsKeyDir1);
1142
1143 Nektar::LocalRegions::HexExpSharedPtr Exp =
1144 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1145 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
1146
1147 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1148 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1149 basisKeyDir1, basisKeyDir1, basisKeyDir1);
1150
1151 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1152 CollExp.push_back(Exp);
1153
1154 LibUtilities::SessionReaderSharedPtr dummySession;
1155 Collections::CollectionOptimisation colOpt(dummySession, 3,
1156 Collections::eMatrixFree);
1157 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1158 Collections::Collection c(CollExp, impTypes);
1159 c.Initialise(Collections::eIProductWRTBase);
1160
1161 const int nq = Exp->GetTotPoints();
1162 Array<OneD, NekDouble> phys(nq);
1163 Array<OneD, NekDouble> coeffsRef(Exp->GetNcoeffs());
1164 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs());
1165
1166 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1167
1168 Exp->GetCoords(xc, yc, zc);
1169
1170 for (int i = 0; i < nq; ++i)
1171 {
1172 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1173 }
1174
1175 Exp->IProductWRTBase(phys, coeffsRef);
1176 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs);
1177
1178 double epsilon = 1.0e-8;
1179 for (int i = 0; i < coeffsRef.size(); ++i)
1180 {
1181 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
1182 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
1183 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
1184 }
1185}
1186
1187BOOST_AUTO_TEST_CASE(
1188 TestHexIProductWRTBase_MatrixFree_UniformP_Deformed_OverInt)
1189{
1190 SpatialDomains::PointGeomUniquePtr v0(
1191 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
1192 SpatialDomains::PointGeomUniquePtr v1(
1193 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1194 SpatialDomains::PointGeomUniquePtr v2(
1195 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1196 SpatialDomains::PointGeomUniquePtr v3(
1197 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1198 SpatialDomains::PointGeomUniquePtr v4(
1199 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1200 SpatialDomains::PointGeomUniquePtr v5(
1201 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1202 SpatialDomains::PointGeomUniquePtr v6(
1203 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
1204 SpatialDomains::PointGeomUniquePtr v7(
1205 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1206
1207 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1208 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1209 SpatialDomains::HexGeomUniquePtr hexGeom =
1210 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1211
1212 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1213 Nektar::LibUtilities::eGaussLobattoLegendre;
1214 Nektar::LibUtilities::BasisType basisTypeDir1 =
1215 Nektar::LibUtilities::eModified_A;
1216 unsigned int numQuadPoints = 8;
1217 unsigned int numModes = 4;
1218 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
1219 quadPointsTypeDir1);
1220 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
1221 quadPointsKeyDir1);
1222
1223 Nektar::LocalRegions::HexExpSharedPtr Exp =
1224 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1225 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
1226
1227 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1228 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1229 basisKeyDir1, basisKeyDir1, basisKeyDir1);
1230
1231 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1232 CollExp.push_back(Exp);
1233
1234 LibUtilities::SessionReaderSharedPtr dummySession;
1235 Collections::CollectionOptimisation colOpt(dummySession, 3,
1236 Collections::eMatrixFree);
1237 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1238 Collections::Collection c(CollExp, impTypes);
1239 c.Initialise(Collections::eIProductWRTBase);
1240
1241 const int nq = Exp->GetTotPoints();
1242 Array<OneD, NekDouble> phys(nq);
1243 Array<OneD, NekDouble> coeffsRef(Exp->GetNcoeffs());
1244 Array<OneD, NekDouble> coeffs(Exp->GetNcoeffs());
1245
1246 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1247
1248 Exp->GetCoords(xc, yc, zc);
1249
1250 for (int i = 0; i < nq; ++i)
1251 {
1252 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1253 }
1254
1255 Exp->IProductWRTBase(phys, coeffsRef);
1256 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs);
1257
1258 double epsilon = 1.0e-8;
1259 for (int i = 0; i < coeffsRef.size(); ++i)
1260 {
1261 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
1262 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
1263 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
1264 }
1265}
1266
1267BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_StdMat_VariableP)
1268{
1269 SpatialDomains::PointGeomUniquePtr v0(
1270 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1271 SpatialDomains::PointGeomUniquePtr v1(
1272 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1273 SpatialDomains::PointGeomUniquePtr v2(
1274 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1275 SpatialDomains::PointGeomUniquePtr v3(
1276 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1277 SpatialDomains::PointGeomUniquePtr v4(
1278 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1279 SpatialDomains::PointGeomUniquePtr v5(
1280 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1281 SpatialDomains::PointGeomUniquePtr v6(
1282 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1283 SpatialDomains::PointGeomUniquePtr v7(
1284 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1285
1286 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1287 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1288 SpatialDomains::HexGeomUniquePtr hexGeom =
1289 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1290
1291 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1292 Nektar::LibUtilities::eGaussLobattoLegendre;
1293 Nektar::LibUtilities::BasisType basisTypeDir1 =
1294 Nektar::LibUtilities::eModified_A;
1295 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1296 quadPointsTypeDir1);
1297 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1298 quadPointsTypeDir1);
1299 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
1300 quadPointsTypeDir1);
1301 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1302 quadPointsKeyDir1);
1303 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1304 quadPointsKeyDir2);
1305 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1306 quadPointsKeyDir3);
1307
1308 Nektar::LocalRegions::HexExpSharedPtr Exp =
1309 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1310 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1311
1312 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1313 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1314 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1315
1316 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1317 CollExp.push_back(Exp);
1318
1319 LibUtilities::SessionReaderSharedPtr dummySession;
1320 Collections::CollectionOptimisation colOpt(dummySession, 3,
1321 Collections::eStdMat);
1322 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1323 Collections::Collection c(CollExp, impTypes);
1324 c.Initialise(Collections::eIProductWRTBase);
1325
1326 const int nq = Exp->GetTotPoints();
1327 Array<OneD, NekDouble> phys(nq);
1328 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1329 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1330
1331 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1332
1333 Exp->GetCoords(xc, yc, zc);
1334
1335 for (int i = 0; i < nq; ++i)
1336 {
1337 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1338 }
1339
1340 Exp->IProductWRTBase(phys, coeffs1);
1341 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1342 c.Initialise(Collections::eIProductWRTBase);
1343
1344 double epsilon = 1.0e-8;
1345 for (int i = 0; i < coeffs1.size(); ++i)
1346 {
1347 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1348 }
1349}
1350
1351BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_NoCollection_VariableP)
1352{
1353 SpatialDomains::PointGeomUniquePtr v0(
1354 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1355 SpatialDomains::PointGeomUniquePtr v1(
1356 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1357 SpatialDomains::PointGeomUniquePtr v2(
1358 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1359 SpatialDomains::PointGeomUniquePtr v3(
1360 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1361 SpatialDomains::PointGeomUniquePtr v4(
1362 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1363 SpatialDomains::PointGeomUniquePtr v5(
1364 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1365 SpatialDomains::PointGeomUniquePtr v6(
1366 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1367 SpatialDomains::PointGeomUniquePtr v7(
1368 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1369
1370 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1371 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1372 SpatialDomains::HexGeomUniquePtr hexGeom =
1373 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1374
1375 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1376 Nektar::LibUtilities::eGaussLobattoLegendre;
1377 Nektar::LibUtilities::BasisType basisTypeDir1 =
1378 Nektar::LibUtilities::eModified_A;
1379 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1380 quadPointsTypeDir1);
1381 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1382 quadPointsTypeDir1);
1383 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
1384 quadPointsTypeDir1);
1385 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1386 quadPointsKeyDir1);
1387 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1388 quadPointsKeyDir2);
1389 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1390 quadPointsKeyDir3);
1391
1392 Nektar::LocalRegions::HexExpSharedPtr Exp =
1393 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1394 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1395
1396 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1397 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1398 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1399
1400 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1401 CollExp.push_back(Exp);
1402
1403 LibUtilities::SessionReaderSharedPtr dummySession;
1404 Collections::CollectionOptimisation colOpt(dummySession, 3,
1405 Collections::eNoCollection);
1406 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1407 Collections::Collection c(CollExp, impTypes);
1408 c.Initialise(Collections::eIProductWRTBase);
1409
1410 const int nq = Exp->GetTotPoints();
1411 Array<OneD, NekDouble> phys(nq);
1412 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1413 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1414
1415 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1416
1417 Exp->GetCoords(xc, yc, zc);
1418
1419 for (int i = 0; i < nq; ++i)
1420 {
1421 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1422 }
1423
1424 Exp->IProductWRTBase(phys, coeffs1);
1425 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1426 c.Initialise(Collections::eIProductWRTBase);
1427
1428 double epsilon = 1.0e-8;
1429 for (int i = 0; i < coeffs1.size(); ++i)
1430 {
1431 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1432 }
1433}
1434
1435BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_VariableP_CollAll)
1436{
1437 SpatialDomains::PointGeomUniquePtr v0(
1438 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1439 SpatialDomains::PointGeomUniquePtr v1(
1440 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1441 SpatialDomains::PointGeomUniquePtr v2(
1442 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1443 SpatialDomains::PointGeomUniquePtr v3(
1444 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1445 SpatialDomains::PointGeomUniquePtr v4(
1446 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1447 SpatialDomains::PointGeomUniquePtr v5(
1448 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1449 SpatialDomains::PointGeomUniquePtr v6(
1450 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1451 SpatialDomains::PointGeomUniquePtr v7(
1452 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1453
1454 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1455 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1456 SpatialDomains::HexGeomUniquePtr hexGeom =
1457 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1458
1459 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1460 Nektar::LibUtilities::eGaussLobattoLegendre;
1461 Nektar::LibUtilities::BasisType basisTypeDir1 =
1462 Nektar::LibUtilities::eGLL_Lagrange;
1463 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(4,
1464 quadPointsTypeDir1);
1465 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
1466 quadPointsTypeDir1);
1467 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
1468 quadPointsTypeDir1);
1469 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1470 quadPointsKeyDir1);
1471 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1472 quadPointsKeyDir2);
1473 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1474 quadPointsKeyDir3);
1475
1476 Nektar::LocalRegions::HexExpSharedPtr Exp =
1477 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1478 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1479
1480 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1481 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1482 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1483
1484 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1485 CollExp.push_back(Exp);
1486
1487 LibUtilities::SessionReaderSharedPtr dummySession;
1488 Collections::CollectionOptimisation colOpt(dummySession, 3,
1489 Collections::eSumFac);
1490 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1491 Collections::Collection c(CollExp, impTypes);
1492 c.Initialise(Collections::eIProductWRTBase);
1493
1494 const int nq = Exp->GetTotPoints();
1495 Array<OneD, NekDouble> phys(nq);
1496 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1497 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1498
1499 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1500
1501 Exp->GetCoords(xc, yc, zc);
1502
1503 for (int i = 0; i < nq; ++i)
1504 {
1505 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1506 }
1507
1508 Exp->IProductWRTBase(phys, coeffs1);
1509 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1510
1511 double epsilon = 1.0e-8;
1512 for (int i = 0; i < coeffs1.size(); ++i)
1513 {
1514 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1515 }
1516}
1517
1518BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_VariableP_CollDir02)
1519{
1520 SpatialDomains::PointGeomUniquePtr v0(
1521 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1522 SpatialDomains::PointGeomUniquePtr v1(
1523 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1524 SpatialDomains::PointGeomUniquePtr v2(
1525 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1526 SpatialDomains::PointGeomUniquePtr v3(
1527 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1528 SpatialDomains::PointGeomUniquePtr v4(
1529 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1530 SpatialDomains::PointGeomUniquePtr v5(
1531 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1532 SpatialDomains::PointGeomUniquePtr v6(
1533 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1534 SpatialDomains::PointGeomUniquePtr v7(
1535 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1536
1537 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1538 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1539 SpatialDomains::HexGeomUniquePtr hexGeom =
1540 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1541
1542 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1543 Nektar::LibUtilities::eGaussLobattoLegendre;
1544 Nektar::LibUtilities::BasisType basisTypeDir1 =
1545 Nektar::LibUtilities::eGLL_Lagrange;
1546 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(4,
1547 quadPointsTypeDir1);
1548 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1549 quadPointsTypeDir1);
1550 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
1551 quadPointsTypeDir1);
1552 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1553 quadPointsKeyDir1);
1554 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1555 quadPointsKeyDir2);
1556 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1557 quadPointsKeyDir3);
1558
1559 Nektar::LocalRegions::HexExpSharedPtr Exp =
1560 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1561 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1562
1563 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1564 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1565 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1566
1567 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1568 CollExp.push_back(Exp);
1569
1570 LibUtilities::SessionReaderSharedPtr dummySession;
1571 Collections::CollectionOptimisation colOpt(dummySession, 3,
1572 Collections::eSumFac);
1573 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1574 Collections::Collection c(CollExp, impTypes);
1575 c.Initialise(Collections::eIProductWRTBase);
1576
1577 const int nq = Exp->GetTotPoints();
1578 Array<OneD, NekDouble> phys(nq);
1579 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1580 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1581
1582 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1583
1584 Exp->GetCoords(xc, yc, zc);
1585
1586 for (int i = 0; i < nq; ++i)
1587 {
1588 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1589 }
1590
1591 Exp->IProductWRTBase(phys, coeffs1);
1592 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1593
1594 double epsilon = 1.0e-8;
1595 for (int i = 0; i < coeffs1.size(); ++i)
1596 {
1597 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1598 }
1599}
1600
1601BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_VariableP_CollDir12)
1602{
1603 SpatialDomains::PointGeomUniquePtr v0(
1604 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1605 SpatialDomains::PointGeomUniquePtr v1(
1606 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1607 SpatialDomains::PointGeomUniquePtr v2(
1608 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1609 SpatialDomains::PointGeomUniquePtr v3(
1610 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1611 SpatialDomains::PointGeomUniquePtr v4(
1612 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1613 SpatialDomains::PointGeomUniquePtr v5(
1614 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1615 SpatialDomains::PointGeomUniquePtr v6(
1616 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1617 SpatialDomains::PointGeomUniquePtr v7(
1618 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1619
1620 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1621 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1622 SpatialDomains::HexGeomUniquePtr hexGeom =
1623 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1624
1625 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1626 Nektar::LibUtilities::eGaussLobattoLegendre;
1627 Nektar::LibUtilities::BasisType basisTypeDir1 =
1628 Nektar::LibUtilities::eGLL_Lagrange;
1629 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1630 quadPointsTypeDir1);
1631 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
1632 quadPointsTypeDir1);
1633 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
1634 quadPointsTypeDir1);
1635 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1636 quadPointsKeyDir1);
1637 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1638 quadPointsKeyDir2);
1639 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1640 quadPointsKeyDir3);
1641
1642 Nektar::LocalRegions::HexExpSharedPtr Exp =
1643 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1644 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1645
1646 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1647 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1648 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1649
1650 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1651 CollExp.push_back(Exp);
1652
1653 LibUtilities::SessionReaderSharedPtr dummySession;
1654 Collections::CollectionOptimisation colOpt(dummySession, 3,
1655 Collections::eSumFac);
1656 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1657 Collections::Collection c(CollExp, impTypes);
1658 c.Initialise(Collections::eIProductWRTBase);
1659
1660 const int nq = Exp->GetTotPoints();
1661 Array<OneD, NekDouble> phys(nq);
1662 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1663 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1664
1665 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1666
1667 Exp->GetCoords(xc, yc, zc);
1668
1669 for (int i = 0; i < nq; ++i)
1670 {
1671 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1672 }
1673
1674 Exp->IProductWRTBase(phys, coeffs1);
1675 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1676
1677 double epsilon = 1.0e-8;
1678 for (int i = 0; i < coeffs1.size(); ++i)
1679 {
1680 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1681 }
1682}
1683
1684BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_StdMat_VariableP_MultiElmt)
1685{
1686 SpatialDomains::PointGeomUniquePtr v0(
1687 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1688 SpatialDomains::PointGeomUniquePtr v1(
1689 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1690 SpatialDomains::PointGeomUniquePtr v2(
1691 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1692 SpatialDomains::PointGeomUniquePtr v3(
1693 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1694 SpatialDomains::PointGeomUniquePtr v4(
1695 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1696 SpatialDomains::PointGeomUniquePtr v5(
1697 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1698 SpatialDomains::PointGeomUniquePtr v6(
1699 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1700 SpatialDomains::PointGeomUniquePtr v7(
1701 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1702
1703 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1704 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1705 SpatialDomains::HexGeomUniquePtr hexGeom =
1706 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1707
1708 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1709 Nektar::LibUtilities::eGaussLobattoLegendre;
1710 Nektar::LibUtilities::BasisType basisTypeDir1 =
1711 Nektar::LibUtilities::eModified_A;
1712 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1713 quadPointsTypeDir1);
1714 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1715 quadPointsTypeDir1);
1716 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
1717 quadPointsTypeDir1);
1718 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1719 quadPointsKeyDir1);
1720 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1721 quadPointsKeyDir2);
1722 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1723 quadPointsKeyDir3);
1724
1725 Nektar::LocalRegions::HexExpSharedPtr Exp =
1726 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1727 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1728
1729 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1730 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1731 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1732
1733 int nelmts = 10;
1734
1735 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1736 for (int i = 0; i < nelmts; ++i)
1737 {
1738 CollExp.push_back(Exp);
1739 }
1740
1741 LibUtilities::SessionReaderSharedPtr dummySession;
1742 Collections::CollectionOptimisation colOpt(dummySession, 3,
1743 Collections::eStdMat);
1744 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1745 Collections::Collection c(CollExp, impTypes);
1746 c.Initialise(Collections::eIProductWRTBase);
1747
1748 const int nq = Exp->GetTotPoints();
1749 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
1750 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
1751 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
1752
1753 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1754
1755 Exp->GetCoords(xc, yc, zc);
1756
1757 for (int i = 0; i < nq; ++i)
1758 {
1759 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1760 }
1761 Exp->IProductWRTBase(phys, coeffs1);
1762
1763 for (int i = 1; i < nelmts; ++i)
1764 {
1765 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
1766 Exp->IProductWRTBase(phys + i * nq,
1767 tmp = coeffs1 + i * Exp->GetNcoeffs());
1768 }
1769
1770 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1771
1772 double epsilon = 1.0e-8;
1773 for (int i = 0; i < coeffs1.size(); ++i)
1774 {
1775 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1776 }
1777}
1778
1779BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_IterPerExp_VariableP_MultiElmt)
1780{
1781 SpatialDomains::PointGeomUniquePtr v0(
1782 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1783 SpatialDomains::PointGeomUniquePtr v1(
1784 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1785 SpatialDomains::PointGeomUniquePtr v2(
1786 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1787 SpatialDomains::PointGeomUniquePtr v3(
1788 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1789 SpatialDomains::PointGeomUniquePtr v4(
1790 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1791 SpatialDomains::PointGeomUniquePtr v5(
1792 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1793 SpatialDomains::PointGeomUniquePtr v6(
1794 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1795 SpatialDomains::PointGeomUniquePtr v7(
1796 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1797
1798 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1799 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1800 SpatialDomains::HexGeomUniquePtr hexGeom =
1801 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1802
1803 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1804 Nektar::LibUtilities::eGaussLobattoLegendre;
1805 Nektar::LibUtilities::BasisType basisTypeDir1 =
1806 Nektar::LibUtilities::eModified_A;
1807 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1808 quadPointsTypeDir1);
1809 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1810 quadPointsTypeDir1);
1811 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
1812 quadPointsTypeDir1);
1813 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1814 quadPointsKeyDir1);
1815 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1816 quadPointsKeyDir2);
1817 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1818 quadPointsKeyDir3);
1819
1820 Nektar::LocalRegions::HexExpSharedPtr Exp =
1821 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1822 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1823
1824 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1825 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1826 basisKeyDir1, basisKeyDir2, basisKeyDir3);
1827
1828 int nelmts = 10;
1829
1830 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1831 for (int i = 0; i < nelmts; ++i)
1832 {
1833 CollExp.push_back(Exp);
1834 }
1835
1836 LibUtilities::SessionReaderSharedPtr dummySession;
1837 Collections::CollectionOptimisation colOpt(dummySession, 3,
1838 Collections::eIterPerExp);
1839 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1840 Collections::Collection c(CollExp, impTypes);
1841 c.Initialise(Collections::eIProductWRTBase);
1842
1843 const int nq = Exp->GetTotPoints();
1844 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
1845 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
1846 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
1847
1848 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1849
1850 Exp->GetCoords(xc, yc, zc);
1851
1852 for (int i = 0; i < nq; ++i)
1853 {
1854 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1855 }
1856 Exp->IProductWRTBase(phys, coeffs1);
1857
1858 for (int i = 1; i < nelmts; ++i)
1859 {
1860 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
1861 Exp->IProductWRTBase(phys + i * nq,
1862 tmp = coeffs1 + i * Exp->GetNcoeffs());
1863 }
1864
1865 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1866
1867 double epsilon = 1.0e-8;
1868 for (int i = 0; i < coeffs1.size(); ++i)
1869 {
1870 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1871 }
1872}
1873
1874BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_UniformP)
1875{
1876 SpatialDomains::PointGeomUniquePtr v0(
1877 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
1878 SpatialDomains::PointGeomUniquePtr v1(
1879 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1880 SpatialDomains::PointGeomUniquePtr v2(
1881 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1882 SpatialDomains::PointGeomUniquePtr v3(
1883 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1884 SpatialDomains::PointGeomUniquePtr v4(
1885 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1886 SpatialDomains::PointGeomUniquePtr v5(
1887 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1888 SpatialDomains::PointGeomUniquePtr v6(
1889 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1890 SpatialDomains::PointGeomUniquePtr v7(
1891 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1892
1893 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1894 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1895 SpatialDomains::HexGeomUniquePtr hexGeom =
1896 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1897
1898 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1899 Nektar::LibUtilities::eGaussLobattoLegendre;
1900 Nektar::LibUtilities::BasisType basisTypeDir1 =
1901 Nektar::LibUtilities::eModified_A;
1902 unsigned int numQuadPoints = 6;
1903 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
1904 quadPointsTypeDir1);
1905 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1906 quadPointsKeyDir1);
1907
1908 Nektar::LocalRegions::HexExpSharedPtr Exp =
1909 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1910 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
1911
1912 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1913 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
1914 basisKeyDir1, basisKeyDir1, basisKeyDir1);
1915
1916 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
1917 CollExp.push_back(Exp);
1918
1919 LibUtilities::SessionReaderSharedPtr dummySession;
1920 Collections::CollectionOptimisation colOpt(dummySession, 3,
1921 Collections::eSumFac);
1922 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
1923 Collections::Collection c(CollExp, impTypes);
1924 c.Initialise(Collections::eIProductWRTBase);
1925
1926 const int nq = Exp->GetTotPoints();
1927 Array<OneD, NekDouble> phys(nq);
1928 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
1929 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
1930
1931 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
1932
1933 Exp->GetCoords(xc, yc, zc);
1934
1935 for (int i = 0; i < nq; ++i)
1936 {
1937 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
1938 }
1939
1940 Exp->IProductWRTBase(phys, coeffs1);
1941 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
1942
1943 double epsilon = 1.0e-6;
1944 for (int i = 0; i < coeffs1.size(); ++i)
1945 {
1946 // clamp values below 1e-16 to zero
1947 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-16) ? 0.0 : coeffs1[i];
1948 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-16) ? 0.0 : coeffs2[i];
1949 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
1950 }
1951}
1952
1953BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_VariableP)
1954{
1955 SpatialDomains::PointGeomUniquePtr v0(
1956 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
1957 SpatialDomains::PointGeomUniquePtr v1(
1958 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
1959 SpatialDomains::PointGeomUniquePtr v2(
1960 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
1961 SpatialDomains::PointGeomUniquePtr v3(
1962 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
1963 SpatialDomains::PointGeomUniquePtr v4(
1964 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
1965 SpatialDomains::PointGeomUniquePtr v5(
1966 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
1967 SpatialDomains::PointGeomUniquePtr v6(
1968 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
1969 SpatialDomains::PointGeomUniquePtr v7(
1970 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
1971
1972 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
1973 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
1974 SpatialDomains::HexGeomUniquePtr hexGeom =
1975 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
1976
1977 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
1978 Nektar::LibUtilities::eGaussLobattoLegendre;
1979 Nektar::LibUtilities::BasisType basisTypeDir1 =
1980 Nektar::LibUtilities::eModified_A;
1981 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
1982 quadPointsTypeDir1);
1983 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
1984 quadPointsTypeDir1);
1985 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
1986 quadPointsTypeDir1);
1987 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
1988 quadPointsKeyDir1);
1989 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
1990 quadPointsKeyDir2);
1991 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
1992 quadPointsKeyDir3);
1993
1994 Nektar::LocalRegions::HexExpSharedPtr Exp =
1995 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
1996 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
1997
1998 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
1999 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2000 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2001
2002 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2003 CollExp.push_back(Exp);
2004
2005 LibUtilities::SessionReaderSharedPtr dummySession;
2006 Collections::CollectionOptimisation colOpt(dummySession, 3,
2007 Collections::eSumFac);
2008 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2009 Collections::Collection c(CollExp, impTypes);
2010 c.Initialise(Collections::eIProductWRTBase);
2011
2012 const int nq = Exp->GetTotPoints();
2013 Array<OneD, NekDouble> phys(nq);
2014 Array<OneD, NekDouble> coeffs1(Exp->GetNcoeffs());
2015 Array<OneD, NekDouble> coeffs2(Exp->GetNcoeffs());
2016
2017 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2018
2019 Exp->GetCoords(xc, yc, zc);
2020
2021 for (int i = 0; i < nq; ++i)
2022 {
2023 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2024 }
2025
2026 Exp->IProductWRTBase(phys, coeffs1);
2027 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
2028
2029 double epsilon = 1.0e-6;
2030 for (int i = 0; i < coeffs1.size(); ++i)
2031 {
2032 // clamp values below 1e-16 to zero
2033 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-16) ? 0.0 : coeffs1[i];
2034 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-16) ? 0.0 : coeffs2[i];
2035 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
2036 }
2037}
2038
2039BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_UniformP_MultiElmt)
2040{
2041 SpatialDomains::PointGeomUniquePtr v0(
2042 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2043 SpatialDomains::PointGeomUniquePtr v1(
2044 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2045 SpatialDomains::PointGeomUniquePtr v2(
2046 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2047 SpatialDomains::PointGeomUniquePtr v3(
2048 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2049 SpatialDomains::PointGeomUniquePtr v4(
2050 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2051 SpatialDomains::PointGeomUniquePtr v5(
2052 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2053 SpatialDomains::PointGeomUniquePtr v6(
2054 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2055 SpatialDomains::PointGeomUniquePtr v7(
2056 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2057
2058 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2059 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2060 SpatialDomains::HexGeomUniquePtr hexGeom =
2061 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2062
2063 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2064 Nektar::LibUtilities::eGaussLobattoLegendre;
2065 Nektar::LibUtilities::BasisType basisTypeDir1 =
2066 Nektar::LibUtilities::eModified_A;
2067 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2068 quadPointsTypeDir1);
2069 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2070 quadPointsKeyDir1);
2071
2072 Nektar::LocalRegions::HexExpSharedPtr Exp =
2073 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2074 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
2075
2076 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2077 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2078 basisKeyDir1, basisKeyDir1, basisKeyDir1);
2079
2080 int nelmts = 10;
2081
2082 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2083 for (int i = 0; i < nelmts; ++i)
2084 {
2085 CollExp.push_back(Exp);
2086 }
2087
2088 LibUtilities::SessionReaderSharedPtr dummySession;
2089 Collections::CollectionOptimisation colOpt(dummySession, 3,
2090 Collections::eSumFac);
2091 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2092 Collections::Collection c(CollExp, impTypes);
2093 c.Initialise(Collections::eIProductWRTBase);
2094
2095 const int nq = Exp->GetTotPoints();
2096 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
2097 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
2098 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
2099
2100 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2101
2102 Exp->GetCoords(xc, yc, zc);
2103
2104 for (int i = 0; i < nq; ++i)
2105 {
2106 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2107 }
2108 Exp->IProductWRTBase(phys, coeffs1);
2109
2110 for (int i = 1; i < nelmts; ++i)
2111 {
2112 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
2113 Exp->IProductWRTBase(phys + i * nq,
2114 tmp = coeffs1 + i * Exp->GetNcoeffs());
2115 }
2116 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
2117
2118 double epsilon = 1.0e-6;
2119 for (int i = 0; i < coeffs1.size(); ++i)
2120 {
2121 // clamp values below 1e-16 to zero
2122 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-16) ? 0.0 : coeffs1[i];
2123 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-16) ? 0.0 : coeffs2[i];
2124 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
2125 }
2126}
2127
2128BOOST_AUTO_TEST_CASE(TestHexIProductWRTBase_SumFac_VariableP_MultiElmt)
2129{
2130 SpatialDomains::PointGeomUniquePtr v0(
2131 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2132 SpatialDomains::PointGeomUniquePtr v1(
2133 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2134 SpatialDomains::PointGeomUniquePtr v2(
2135 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2136 SpatialDomains::PointGeomUniquePtr v3(
2137 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2138 SpatialDomains::PointGeomUniquePtr v4(
2139 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2140 SpatialDomains::PointGeomUniquePtr v5(
2141 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2142 SpatialDomains::PointGeomUniquePtr v6(
2143 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2144 SpatialDomains::PointGeomUniquePtr v7(
2145 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2146
2147 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2148 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2149 SpatialDomains::HexGeomUniquePtr hexGeom =
2150 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2151
2152 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2153 Nektar::LibUtilities::eGaussLobattoLegendre;
2154 Nektar::LibUtilities::BasisType basisTypeDir1 =
2155 Nektar::LibUtilities::eModified_A;
2156 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2157 quadPointsTypeDir1);
2158 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(7,
2159 quadPointsTypeDir1);
2160 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(9,
2161 quadPointsTypeDir1);
2162 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2163 quadPointsKeyDir1);
2164 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2165 quadPointsKeyDir2);
2166 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2167 quadPointsKeyDir3);
2168
2169 Nektar::LocalRegions::HexExpSharedPtr Exp =
2170 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2171 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2172
2173 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2174 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2175 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2176
2177 int nelmts = 10;
2178
2179 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2180 for (int i = 0; i < nelmts; ++i)
2181 {
2182 CollExp.push_back(Exp);
2183 }
2184
2185 LibUtilities::SessionReaderSharedPtr dummySession;
2186 Collections::CollectionOptimisation colOpt(dummySession, 3,
2187 Collections::eSumFac);
2188 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2189 Collections::Collection c(CollExp, impTypes);
2190 c.Initialise(Collections::eIProductWRTBase);
2191
2192 const int nq = Exp->GetTotPoints();
2193 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
2194 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
2195 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
2196
2197 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2198
2199 Exp->GetCoords(xc, yc, zc);
2200
2201 for (int i = 0; i < nq; ++i)
2202 {
2203 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2204 }
2205 Exp->IProductWRTBase(phys, coeffs1);
2206
2207 for (int i = 1; i < nelmts; ++i)
2208 {
2209 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
2210 Exp->IProductWRTBase(phys + i * nq,
2211 tmp = coeffs1 + i * Exp->GetNcoeffs());
2212 }
2213 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
2214
2215 double epsilon = 1.0e-4;
2216 for (int i = 0; i < coeffs1.size(); ++i)
2217 {
2218 // clamp values below 1e-14 to zero
2219 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
2220 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
2221 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
2222 }
2223}
2224
2225BOOST_AUTO_TEST_CASE(
2226 TestHexIProductWRTBase_SumFac_VariableP_MultiElmt_CollDir02)
2227{
2228 SpatialDomains::PointGeomUniquePtr v0(
2229 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2230 SpatialDomains::PointGeomUniquePtr v1(
2231 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2232 SpatialDomains::PointGeomUniquePtr v2(
2233 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2234 SpatialDomains::PointGeomUniquePtr v3(
2235 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2236 SpatialDomains::PointGeomUniquePtr v4(
2237 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2238 SpatialDomains::PointGeomUniquePtr v5(
2239 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2240 SpatialDomains::PointGeomUniquePtr v6(
2241 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2242 SpatialDomains::PointGeomUniquePtr v7(
2243 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2244
2245 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2246 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2247 SpatialDomains::HexGeomUniquePtr hexGeom =
2248 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2249
2250 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2251 Nektar::LibUtilities::eGaussLobattoLegendre;
2252 Nektar::LibUtilities::BasisType basisTypeDir1 =
2253 Nektar::LibUtilities::eGLL_Lagrange;
2254 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(4,
2255 quadPointsTypeDir1);
2256 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
2257 quadPointsTypeDir1);
2258 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
2259 quadPointsTypeDir1);
2260 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2261 quadPointsKeyDir1);
2262 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2263 quadPointsKeyDir2);
2264 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2265 quadPointsKeyDir3);
2266
2267 Nektar::LocalRegions::HexExpSharedPtr Exp =
2268 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2269 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2270
2271 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2272 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2273 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2274
2275 int nelmts = 10;
2276
2277 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2278 for (int i = 0; i < nelmts; ++i)
2279 {
2280 CollExp.push_back(Exp);
2281 }
2282
2283 LibUtilities::SessionReaderSharedPtr dummySession;
2284 Collections::CollectionOptimisation colOpt(dummySession, 3,
2285 Collections::eSumFac);
2286 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2287 Collections::Collection c(CollExp, impTypes);
2288 c.Initialise(Collections::eIProductWRTBase);
2289
2290 const int nq = Exp->GetTotPoints();
2291 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
2292 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
2293 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
2294
2295 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2296
2297 Exp->GetCoords(xc, yc, zc);
2298
2299 for (int i = 0; i < nq; ++i)
2300 {
2301 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2302 }
2303 Exp->IProductWRTBase(phys, coeffs1);
2304
2305 for (int i = 1; i < nelmts; ++i)
2306 {
2307 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
2308 Exp->IProductWRTBase(phys + i * nq,
2309 tmp = coeffs1 + i * Exp->GetNcoeffs());
2310 }
2311 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
2312
2313 double epsilon = 1.0e-4;
2314 for (int i = 0; i < coeffs1.size(); ++i)
2315 {
2316 // clamp values below 1e-14 to zero
2317 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
2318 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
2319 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
2320 }
2321}
2322
2323BOOST_AUTO_TEST_CASE(
2324 TestHexIProductWRTBase_SumFac_VariableP_MultiElmt_CollDir12)
2325{
2326 SpatialDomains::PointGeomUniquePtr v0(
2327 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2328 SpatialDomains::PointGeomUniquePtr v1(
2329 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2330 SpatialDomains::PointGeomUniquePtr v2(
2331 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2332 SpatialDomains::PointGeomUniquePtr v3(
2333 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2334 SpatialDomains::PointGeomUniquePtr v4(
2335 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2336 SpatialDomains::PointGeomUniquePtr v5(
2337 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2338 SpatialDomains::PointGeomUniquePtr v6(
2339 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2340 SpatialDomains::PointGeomUniquePtr v7(
2341 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2342
2343 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2344 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2345 SpatialDomains::HexGeomUniquePtr hexGeom =
2346 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2347
2348 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2349 Nektar::LibUtilities::eGaussLobattoLegendre;
2350 Nektar::LibUtilities::BasisType basisTypeDir1 =
2351 Nektar::LibUtilities::eGLL_Lagrange;
2352 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2353 quadPointsTypeDir1);
2354 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
2355 quadPointsTypeDir1);
2356 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
2357 quadPointsTypeDir1);
2358 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2359 quadPointsKeyDir1);
2360 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2361 quadPointsKeyDir2);
2362 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2363 quadPointsKeyDir3);
2364
2365 Nektar::LocalRegions::HexExpSharedPtr Exp =
2366 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2367 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2368
2369 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2370 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2371 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2372
2373 int nelmts = 10;
2374
2375 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2376 for (int i = 0; i < nelmts; ++i)
2377 {
2378 CollExp.push_back(Exp);
2379 }
2380
2381 LibUtilities::SessionReaderSharedPtr dummySession;
2382 Collections::CollectionOptimisation colOpt(dummySession, 3,
2383 Collections::eSumFac);
2384 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2385 Collections::Collection c(CollExp, impTypes);
2386 c.Initialise(Collections::eIProductWRTBase);
2387
2388 const int nq = Exp->GetTotPoints();
2389 Array<OneD, NekDouble> phys(nelmts * nq), tmp;
2390 Array<OneD, NekDouble> coeffs1(nelmts * Exp->GetNcoeffs());
2391 Array<OneD, NekDouble> coeffs2(nelmts * Exp->GetNcoeffs());
2392
2393 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2394
2395 Exp->GetCoords(xc, yc, zc);
2396
2397 for (int i = 0; i < nq; ++i)
2398 {
2399 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2400 }
2401 Exp->IProductWRTBase(phys, coeffs1);
2402
2403 for (int i = 1; i < nelmts; ++i)
2404 {
2405 Vmath::Vcopy(nq, &phys[0], 1, &phys[i * nq], 1);
2406 Exp->IProductWRTBase(phys + i * nq,
2407 tmp = coeffs1 + i * Exp->GetNcoeffs());
2408 }
2409 c.ApplyOperator(Collections::eIProductWRTBase, phys, coeffs2);
2410
2411 double epsilon = 1.0e-4;
2412 for (int i = 0; i < coeffs1.size(); ++i)
2413 {
2414 // clamp values below 1e-14 to zero
2415 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
2416 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
2417 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
2418 }
2419}
2420
2421BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_IterPerExp_UniformP)
2422{
2423 SpatialDomains::PointGeomUniquePtr v0(
2424 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2425 SpatialDomains::PointGeomUniquePtr v1(
2426 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2427 SpatialDomains::PointGeomUniquePtr v2(
2428 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2429 SpatialDomains::PointGeomUniquePtr v3(
2430 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2431 SpatialDomains::PointGeomUniquePtr v4(
2432 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2433 SpatialDomains::PointGeomUniquePtr v5(
2434 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2435 SpatialDomains::PointGeomUniquePtr v6(
2436 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2437 SpatialDomains::PointGeomUniquePtr v7(
2438 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2439
2440 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2441 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2442 SpatialDomains::HexGeomUniquePtr hexGeom =
2443 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2444
2445 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2446 Nektar::LibUtilities::eGaussLobattoLegendre;
2447 Nektar::LibUtilities::BasisType basisTypeDir1 =
2448 Nektar::LibUtilities::eModified_A;
2449 unsigned int numQuadPoints = 6;
2450 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
2451 quadPointsTypeDir1);
2452 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2453 quadPointsKeyDir1);
2454
2455 Nektar::LocalRegions::HexExpSharedPtr Exp =
2456 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2457 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
2458
2459 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2460 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2461 basisKeyDir1, basisKeyDir1, basisKeyDir1);
2462
2463 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2464 CollExp.push_back(Exp);
2465
2466 LibUtilities::SessionReaderSharedPtr dummySession;
2467 Collections::CollectionOptimisation colOpt(dummySession, 3,
2468 Collections::eIterPerExp);
2469 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2470 Collections::Collection c(CollExp, impTypes);
2471 c.Initialise(Collections::ePhysDeriv);
2472
2473 const int nq = Exp->GetTotPoints();
2474 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2475 Array<OneD, NekDouble> phys(nq), tmp, tmp1, tmp2;
2476 Array<OneD, NekDouble> diff1(3 * nq);
2477 Array<OneD, NekDouble> diff2(3 * nq);
2478
2479 Exp->GetCoords(xc, yc, zc);
2480
2481 for (int i = 0; i < nq; ++i)
2482 {
2483 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2484 }
2485
2486 Exp->PhysDeriv(phys, diff1, tmp = diff1 + nq, tmp1 = diff1 + 2 * nq);
2487 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2, tmp = diff2 + nq,
2488 tmp2 = diff2 + 2 * nq);
2489
2490 double epsilon = 1.0e-8;
2491 for (int i = 0; i < diff1.size(); ++i)
2492 {
2493 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
2494 }
2495}
2496
2497BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_MatrixFree_UniformP_Undeformed)
2498{
2499 SpatialDomains::PointGeomUniquePtr v0(
2500 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
2501 SpatialDomains::PointGeomUniquePtr v1(
2502 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2503 SpatialDomains::PointGeomUniquePtr v2(
2504 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2505 SpatialDomains::PointGeomUniquePtr v3(
2506 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2507 SpatialDomains::PointGeomUniquePtr v4(
2508 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2509 SpatialDomains::PointGeomUniquePtr v5(
2510 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2511 SpatialDomains::PointGeomUniquePtr v6(
2512 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2513 SpatialDomains::PointGeomUniquePtr v7(
2514 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2515
2516 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2517 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2518 SpatialDomains::HexGeomUniquePtr hexGeom =
2519 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2520
2521 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2522 Nektar::LibUtilities::eGaussLobattoLegendre;
2523 Nektar::LibUtilities::BasisType basisTypeDir1 =
2524 Nektar::LibUtilities::eModified_A;
2525 unsigned int numQuadPoints = 6;
2526 unsigned int numModes = 2;
2527 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
2528 quadPointsTypeDir1);
2529 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
2530 quadPointsKeyDir1);
2531
2532 Nektar::LocalRegions::HexExpSharedPtr Exp =
2533 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2534 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
2535
2536 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2537 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2538 basisKeyDir1, basisKeyDir1, basisKeyDir1);
2539
2540 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2541 CollExp.push_back(Exp);
2542
2543 LibUtilities::SessionReaderSharedPtr dummySession;
2544 Collections::CollectionOptimisation colOpt(dummySession, 2,
2545 Collections::eMatrixFree);
2546 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2547 Collections::Collection c(CollExp, impTypes);
2548 c.Initialise(Collections::ePhysDeriv);
2549
2550 const int nq = Exp->GetTotPoints();
2551 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2552 Array<OneD, NekDouble> phys(nq), tmp, tmp1, tmp2;
2553 Array<OneD, NekDouble> diffRef(3 * nq);
2554 Array<OneD, NekDouble> diff(3 * nq);
2555
2556 Exp->GetCoords(xc, yc, zc);
2557
2558 for (int i = 0; i < nq; ++i)
2559 {
2560 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2561 }
2562
2563 Exp->PhysDeriv(phys, diffRef, tmp = diffRef + nq, tmp1 = diffRef + 2 * nq);
2564 c.ApplyOperator(Collections::ePhysDeriv, phys, diff, tmp = diff + nq,
2565 tmp2 = diff + 2 * nq);
2566
2567 double epsilon = 1.0e-8;
2568 for (int i = 0; i < diffRef.size(); ++i)
2569 {
2570 diffRef[i] = (std::abs(diffRef[i]) < 1e-14) ? 0.0 : diffRef[i];
2571 diff[i] = (std::abs(diff[i]) < 1e-14) ? 0.0 : diff[i];
2572 BOOST_CHECK_CLOSE(diffRef[i], diff[i], epsilon);
2573 }
2574}
2575
2576BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_MatrixFree_UniformP_Deformed)
2577{
2578 SpatialDomains::PointGeomUniquePtr v0(
2579 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
2580 SpatialDomains::PointGeomUniquePtr v1(
2581 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2582 SpatialDomains::PointGeomUniquePtr v2(
2583 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2584 SpatialDomains::PointGeomUniquePtr v3(
2585 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2586 SpatialDomains::PointGeomUniquePtr v4(
2587 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2588 SpatialDomains::PointGeomUniquePtr v5(
2589 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2590 SpatialDomains::PointGeomUniquePtr v6(
2591 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
2592 SpatialDomains::PointGeomUniquePtr v7(
2593 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2594
2595 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2596 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2597 SpatialDomains::HexGeomUniquePtr hexGeom =
2598 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2599
2600 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2601 Nektar::LibUtilities::eGaussLobattoLegendre;
2602 Nektar::LibUtilities::BasisType basisTypeDir1 =
2603 Nektar::LibUtilities::eModified_A;
2604 unsigned int numQuadPoints = 5;
2605 unsigned int numModes = 2;
2606 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
2607 quadPointsTypeDir1);
2608 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
2609 quadPointsKeyDir1);
2610
2611 Nektar::LocalRegions::HexExpSharedPtr Exp =
2612 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2613 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
2614
2615 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2616 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2617 basisKeyDir1, basisKeyDir1, basisKeyDir1);
2618
2619 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2620 CollExp.push_back(Exp);
2621
2622 LibUtilities::SessionReaderSharedPtr dummySession;
2623 Collections::CollectionOptimisation colOpt(dummySession, 2,
2624 Collections::eMatrixFree);
2625 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2626 Collections::Collection c(CollExp, impTypes);
2627 c.Initialise(Collections::ePhysDeriv);
2628
2629 const int nq = Exp->GetTotPoints();
2630 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2631 Array<OneD, NekDouble> phys(nq), tmp, tmp1, tmp2;
2632 Array<OneD, NekDouble> diffRef(3 * nq);
2633 Array<OneD, NekDouble> diff(3 * nq);
2634
2635 Exp->GetCoords(xc, yc, zc);
2636
2637 for (int i = 0; i < nq; ++i)
2638 {
2639 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2640 }
2641
2642 Exp->PhysDeriv(phys, diffRef, tmp = diffRef + nq, tmp1 = diffRef + 2 * nq);
2643 c.ApplyOperator(Collections::ePhysDeriv, phys, diff, tmp = diff + nq,
2644 tmp2 = diff + 2 * nq);
2645
2646 double epsilon = 1.0e-8;
2647 for (int i = 0; i < diffRef.size(); ++i)
2648 {
2649 diffRef[i] = (std::abs(diffRef[i]) < 1e-14) ? 0.0 : diffRef[i];
2650 diff[i] = (std::abs(diff[i]) < 1e-14) ? 0.0 : diff[i];
2651 BOOST_CHECK_CLOSE(diffRef[i], diff[i], epsilon);
2652 }
2653}
2654
2655BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_IterPerExp_VariableP_MultiElmt)
2656{
2657 SpatialDomains::PointGeomUniquePtr v0(
2658 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2659 SpatialDomains::PointGeomUniquePtr v1(
2660 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2661 SpatialDomains::PointGeomUniquePtr v2(
2662 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2663 SpatialDomains::PointGeomUniquePtr v3(
2664 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2665 SpatialDomains::PointGeomUniquePtr v4(
2666 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2667 SpatialDomains::PointGeomUniquePtr v5(
2668 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2669 SpatialDomains::PointGeomUniquePtr v6(
2670 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2671 SpatialDomains::PointGeomUniquePtr v7(
2672 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2673
2674 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2675 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2676 SpatialDomains::HexGeomUniquePtr hexGeom =
2677 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2678
2679 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2680 Nektar::LibUtilities::eGaussLobattoLegendre;
2681 Nektar::LibUtilities::BasisType basisTypeDir1 =
2682 Nektar::LibUtilities::eModified_A;
2683 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2684 quadPointsTypeDir1);
2685 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
2686 quadPointsTypeDir1);
2687 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
2688 quadPointsTypeDir1);
2689 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2690 quadPointsKeyDir1);
2691 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2692 quadPointsKeyDir2);
2693 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2694 quadPointsKeyDir3);
2695
2696 Nektar::LocalRegions::HexExpSharedPtr Exp =
2697 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2698 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2699
2700 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2701 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2702 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2703
2704 int nelmts = 10;
2705
2706 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2707 for (int i = 0; i < nelmts; ++i)
2708 {
2709 CollExp.push_back(Exp);
2710 }
2711
2712 LibUtilities::SessionReaderSharedPtr dummySession;
2713 Collections::CollectionOptimisation colOpt(dummySession, 3,
2714 Collections::eIterPerExp);
2715 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2716 Collections::Collection c(CollExp, impTypes);
2717 c.Initialise(Collections::ePhysDeriv);
2718
2719 const int nq = Exp->GetTotPoints();
2720 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2721 Array<OneD, NekDouble> phys(nelmts * nq), tmp, tmp1, tmp2;
2722 Array<OneD, NekDouble> diff1(3 * nelmts * nq);
2723 Array<OneD, NekDouble> diff2(3 * nelmts * nq);
2724
2725 Exp->GetCoords(xc, yc, zc);
2726
2727 for (int i = 0; i < nq; ++i)
2728 {
2729 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2730 }
2731 Exp->PhysDeriv(phys, tmp = diff1, tmp1 = diff1 + (nelmts)*nq,
2732 tmp2 = diff1 + (2 * nelmts) * nq);
2733 for (int i = 1; i < nelmts; ++i)
2734 {
2735 Vmath::Vcopy(nq, phys, 1, tmp = phys + i * nq, 1);
2736 Exp->PhysDeriv(phys, tmp = diff1 + i * nq,
2737 tmp1 = diff1 + (nelmts + i) * nq,
2738 tmp2 = diff1 + (2 * nelmts + i) * nq);
2739 }
2740
2741 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2,
2742 tmp = diff2 + nelmts * nq, tmp2 = diff2 + 2 * nelmts * nq);
2743
2744 double epsilon = 1.0e-8;
2745 for (int i = 0; i < diff1.size(); ++i)
2746 {
2747 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
2748 }
2749}
2750
2751BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_NoCollection_VariableP_MultiElmt)
2752{
2753 SpatialDomains::PointGeomUniquePtr v0(
2754 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2755 SpatialDomains::PointGeomUniquePtr v1(
2756 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2757 SpatialDomains::PointGeomUniquePtr v2(
2758 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2759 SpatialDomains::PointGeomUniquePtr v3(
2760 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2761 SpatialDomains::PointGeomUniquePtr v4(
2762 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2763 SpatialDomains::PointGeomUniquePtr v5(
2764 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2765 SpatialDomains::PointGeomUniquePtr v6(
2766 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2767 SpatialDomains::PointGeomUniquePtr v7(
2768 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2769
2770 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2771 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2772 SpatialDomains::HexGeomUniquePtr hexGeom =
2773 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2774
2775 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2776 Nektar::LibUtilities::eGaussLobattoLegendre;
2777 Nektar::LibUtilities::BasisType basisTypeDir1 =
2778 Nektar::LibUtilities::eModified_A;
2779 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2780 quadPointsTypeDir1);
2781 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
2782 quadPointsTypeDir1);
2783 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
2784 quadPointsTypeDir1);
2785 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2786 quadPointsKeyDir1);
2787 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2788 quadPointsKeyDir2);
2789 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2790 quadPointsKeyDir3);
2791
2792 Nektar::LocalRegions::HexExpSharedPtr Exp =
2793 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2794 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2795
2796 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2797 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2798 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2799
2800 int nelmts = 10;
2801
2802 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2803 for (int i = 0; i < nelmts; ++i)
2804 {
2805 CollExp.push_back(Exp);
2806 }
2807
2808 LibUtilities::SessionReaderSharedPtr dummySession;
2809 Collections::CollectionOptimisation colOpt(dummySession, 3,
2810 Collections::eNoCollection);
2811 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2812 Collections::Collection c(CollExp, impTypes);
2813 c.Initialise(Collections::ePhysDeriv);
2814
2815 const int nq = Exp->GetTotPoints();
2816 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2817 Array<OneD, NekDouble> phys(nelmts * nq), tmp, tmp1, tmp2;
2818 Array<OneD, NekDouble> diff1(3 * nelmts * nq);
2819 Array<OneD, NekDouble> diff2(3 * nelmts * nq);
2820
2821 Exp->GetCoords(xc, yc, zc);
2822
2823 for (int i = 0; i < nq; ++i)
2824 {
2825 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2826 }
2827 Exp->PhysDeriv(phys, tmp = diff1, tmp1 = diff1 + (nelmts)*nq,
2828 tmp2 = diff1 + (2 * nelmts) * nq);
2829
2830 for (int i = 1; i < nelmts; ++i)
2831 {
2832 Vmath::Vcopy(nq, phys, 1, tmp = phys + i * nq, 1);
2833 Exp->PhysDeriv(phys, tmp = diff1 + i * nq,
2834 tmp1 = diff1 + (nelmts + i) * nq,
2835 tmp2 = diff1 + (2 * nelmts + i) * nq);
2836 }
2837
2838 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2,
2839 tmp = diff2 + nelmts * nq, tmp2 = diff2 + 2 * nelmts * nq);
2840
2841 double epsilon = 1.0e-8;
2842 for (int i = 0; i < diff1.size(); ++i)
2843 {
2844 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
2845 }
2846}
2847
2848BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_StdMat_UniformP)
2849{
2850 SpatialDomains::PointGeomUniquePtr v0(
2851 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2852 SpatialDomains::PointGeomUniquePtr v1(
2853 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2854 SpatialDomains::PointGeomUniquePtr v2(
2855 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2856 SpatialDomains::PointGeomUniquePtr v3(
2857 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2858 SpatialDomains::PointGeomUniquePtr v4(
2859 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2860 SpatialDomains::PointGeomUniquePtr v5(
2861 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2862 SpatialDomains::PointGeomUniquePtr v6(
2863 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2864 SpatialDomains::PointGeomUniquePtr v7(
2865 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2866
2867 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2868 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2869 SpatialDomains::HexGeomUniquePtr hexGeom =
2870 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2871
2872 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2873 Nektar::LibUtilities::eGaussLobattoLegendre;
2874 Nektar::LibUtilities::BasisType basisTypeDir1 =
2875 Nektar::LibUtilities::eModified_A;
2876 unsigned int numQuadPoints = 4;
2877 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
2878 quadPointsTypeDir1);
2879 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 3,
2880 quadPointsKeyDir1);
2881
2882 Nektar::LocalRegions::HexExpSharedPtr Exp =
2883 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2884 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
2885
2886 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2887 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2888 basisKeyDir1, basisKeyDir1, basisKeyDir1);
2889
2890 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2891 CollExp.push_back(Exp);
2892
2893 LibUtilities::SessionReaderSharedPtr dummySession;
2894 Collections::CollectionOptimisation colOpt(dummySession, 3,
2895 Collections::eStdMat);
2896 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2897 Collections::Collection c(CollExp, impTypes);
2898 c.Initialise(Collections::ePhysDeriv);
2899
2900 const int nq = Exp->GetTotPoints();
2901 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2902 Array<OneD, NekDouble> phys(nq), tmp, tmp1, tmp2;
2903 Array<OneD, NekDouble> diff1(3 * nq);
2904 Array<OneD, NekDouble> diff2(3 * nq);
2905
2906 Exp->GetCoords(xc, yc, zc);
2907
2908 for (int i = 0; i < nq; ++i)
2909 {
2910 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2911 }
2912
2913 Exp->PhysDeriv(phys, diff1, tmp = diff1 + nq, tmp1 = diff1 + 2 * nq);
2914 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2, tmp = diff2 + nq,
2915 tmp2 = diff2 + 2 * nq);
2916
2917 double epsilon = 1.0e-8;
2918 for (int i = 0; i < diff1.size(); ++i)
2919 {
2920 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
2921 }
2922}
2923
2924BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_StdMat_VariableP_MultiElmt)
2925{
2926 SpatialDomains::PointGeomUniquePtr v0(
2927 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
2928 SpatialDomains::PointGeomUniquePtr v1(
2929 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
2930 SpatialDomains::PointGeomUniquePtr v2(
2931 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
2932 SpatialDomains::PointGeomUniquePtr v3(
2933 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
2934 SpatialDomains::PointGeomUniquePtr v4(
2935 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
2936 SpatialDomains::PointGeomUniquePtr v5(
2937 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
2938 SpatialDomains::PointGeomUniquePtr v6(
2939 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
2940 SpatialDomains::PointGeomUniquePtr v7(
2941 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
2942
2943 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
2944 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
2945 SpatialDomains::HexGeomUniquePtr hexGeom =
2946 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
2947
2948 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
2949 Nektar::LibUtilities::eGaussLobattoLegendre;
2950 Nektar::LibUtilities::BasisType basisTypeDir1 =
2951 Nektar::LibUtilities::eModified_A;
2952 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
2953 quadPointsTypeDir1);
2954 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
2955 quadPointsTypeDir1);
2956 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
2957 quadPointsTypeDir1);
2958 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
2959 quadPointsKeyDir1);
2960 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
2961 quadPointsKeyDir2);
2962 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
2963 quadPointsKeyDir3);
2964
2965 Nektar::LocalRegions::HexExpSharedPtr Exp =
2966 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
2967 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
2968
2969 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
2970 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
2971 basisKeyDir1, basisKeyDir2, basisKeyDir3);
2972
2973 int nelmts = 10;
2974
2975 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
2976 for (int i = 0; i < nelmts; ++i)
2977 {
2978 CollExp.push_back(Exp);
2979 }
2980
2981 LibUtilities::SessionReaderSharedPtr dummySession;
2982 Collections::CollectionOptimisation colOpt(dummySession, 3,
2983 Collections::eStdMat);
2984 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
2985 Collections::Collection c(CollExp, impTypes);
2986 c.Initialise(Collections::ePhysDeriv);
2987
2988 const int nq = Exp->GetTotPoints();
2989 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
2990 Array<OneD, NekDouble> phys(nelmts * nq), tmp, tmp1, tmp2;
2991 Array<OneD, NekDouble> diff1(3 * nelmts * nq);
2992 Array<OneD, NekDouble> diff2(3 * nelmts * nq);
2993
2994 Exp->GetCoords(xc, yc, zc);
2995
2996 for (int i = 0; i < nq; ++i)
2997 {
2998 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
2999 }
3000 Exp->PhysDeriv(phys, tmp = diff1, tmp1 = diff1 + (nelmts)*nq,
3001 tmp2 = diff1 + (2 * nelmts) * nq);
3002 for (int i = 1; i < nelmts; ++i)
3003 {
3004 Vmath::Vcopy(nq, phys, 1, tmp = phys + i * nq, 1);
3005 Exp->PhysDeriv(phys, tmp = diff1 + i * nq,
3006 tmp1 = diff1 + (nelmts + i) * nq,
3007 tmp2 = diff1 + (2 * nelmts + i) * nq);
3008 }
3009
3010 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2,
3011 tmp = diff2 + nelmts * nq, tmp2 = diff2 + 2 * nelmts * nq);
3012
3013 double epsilon = 1.0e-8;
3014 for (int i = 0; i < diff1.size(); ++i)
3015 {
3016 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
3017 }
3018}
3019
3020BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_SumFac_UniformP)
3021{
3022 SpatialDomains::PointGeomUniquePtr v0(
3023 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3024 SpatialDomains::PointGeomUniquePtr v1(
3025 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3026 SpatialDomains::PointGeomUniquePtr v2(
3027 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3028 SpatialDomains::PointGeomUniquePtr v3(
3029 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3030 SpatialDomains::PointGeomUniquePtr v4(
3031 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3032 SpatialDomains::PointGeomUniquePtr v5(
3033 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3034 SpatialDomains::PointGeomUniquePtr v6(
3035 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3036 SpatialDomains::PointGeomUniquePtr v7(
3037 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3038
3039 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3040 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3041 SpatialDomains::HexGeomUniquePtr hexGeom =
3042 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3043
3044 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3045 Nektar::LibUtilities::eGaussLobattoLegendre;
3046 Nektar::LibUtilities::BasisType basisTypeDir1 =
3047 Nektar::LibUtilities::eModified_A;
3048 unsigned int numQuadPoints = 6;
3049 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3050 quadPointsTypeDir1);
3051 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3052 quadPointsKeyDir1);
3053
3054 Nektar::LocalRegions::HexExpSharedPtr Exp =
3055 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3056 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3057
3058 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3059 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3060 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3061
3062 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3063 CollExp.push_back(Exp);
3064
3065 LibUtilities::SessionReaderSharedPtr dummySession;
3066 Collections::CollectionOptimisation colOpt(dummySession, 3,
3067 Collections::eSumFac);
3068 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3069 Collections::Collection c(CollExp, impTypes);
3070 c.Initialise(Collections::ePhysDeriv);
3071
3072 const int nq = Exp->GetTotPoints();
3073 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3074 Array<OneD, NekDouble> phys(nq), tmp, tmp1, tmp2;
3075 Array<OneD, NekDouble> diff1(3 * nq);
3076 Array<OneD, NekDouble> diff2(3 * nq);
3077
3078 Exp->GetCoords(xc, yc, zc);
3079
3080 for (int i = 0; i < nq; ++i)
3081 {
3082 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3083 }
3084
3085 Exp->PhysDeriv(phys, diff1, tmp = diff1 + nq, tmp1 = diff1 + 2 * nq);
3086 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2, tmp = diff2 + nq,
3087 tmp2 = diff2 + 2 * nq);
3088
3089 double epsilon = 1.0e-8;
3090 for (int i = 0; i < diff1.size(); ++i)
3091 {
3092 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
3093 }
3094}
3095
3096BOOST_AUTO_TEST_CASE(TestHexPhysDeriv_SumFac_VariableP_MultiElmt)
3097{
3098 SpatialDomains::PointGeomUniquePtr v0(
3099 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3100 SpatialDomains::PointGeomUniquePtr v1(
3101 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3102 SpatialDomains::PointGeomUniquePtr v2(
3103 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3104 SpatialDomains::PointGeomUniquePtr v3(
3105 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3106 SpatialDomains::PointGeomUniquePtr v4(
3107 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3108 SpatialDomains::PointGeomUniquePtr v5(
3109 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3110 SpatialDomains::PointGeomUniquePtr v6(
3111 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3112 SpatialDomains::PointGeomUniquePtr v7(
3113 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3114
3115 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3116 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3117 SpatialDomains::HexGeomUniquePtr hexGeom =
3118 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3119
3120 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3121 Nektar::LibUtilities::eGaussLobattoLegendre;
3122 Nektar::LibUtilities::BasisType basisTypeDir1 =
3123 Nektar::LibUtilities::eModified_A;
3124 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
3125 quadPointsTypeDir1);
3126 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
3127 quadPointsTypeDir1);
3128 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
3129 quadPointsTypeDir1);
3130 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3131 quadPointsKeyDir1);
3132 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
3133 quadPointsKeyDir2);
3134 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
3135 quadPointsKeyDir3);
3136
3137 Nektar::LocalRegions::HexExpSharedPtr Exp =
3138 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3139 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
3140
3141 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3142 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3143 basisKeyDir1, basisKeyDir2, basisKeyDir3);
3144
3145 int nelmts = 10;
3146
3147 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3148 for (int i = 0; i < nelmts; ++i)
3149 {
3150 CollExp.push_back(Exp);
3151 }
3152
3153 LibUtilities::SessionReaderSharedPtr dummySession;
3154 Collections::CollectionOptimisation colOpt(dummySession, 3,
3155 Collections::eSumFac);
3156 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3157 Collections::Collection c(CollExp, impTypes);
3158 c.Initialise(Collections::ePhysDeriv);
3159
3160 const int nq = Exp->GetTotPoints();
3161 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3162 Array<OneD, NekDouble> phys(nelmts * nq), tmp, tmp1, tmp2;
3163 Array<OneD, NekDouble> diff1(3 * nelmts * nq);
3164 Array<OneD, NekDouble> diff2(3 * nelmts * nq);
3165
3166 Exp->GetCoords(xc, yc, zc);
3167
3168 for (int i = 0; i < nq; ++i)
3169 {
3170 phys[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3171 }
3172 Exp->PhysDeriv(phys, tmp = diff1, tmp1 = diff1 + (nelmts)*nq,
3173 tmp2 = diff1 + (2 * nelmts) * nq);
3174 for (int i = 1; i < nelmts; ++i)
3175 {
3176 Vmath::Vcopy(nq, phys, 1, tmp = phys + i * nq, 1);
3177 Exp->PhysDeriv(phys, tmp = diff1 + i * nq,
3178 tmp1 = diff1 + (nelmts + i) * nq,
3179 tmp2 = diff1 + (2 * nelmts + i) * nq);
3180 }
3181
3182 c.ApplyOperator(Collections::ePhysDeriv, phys, diff2,
3183 tmp = diff2 + nelmts * nq, tmp2 = diff2 + 2 * nelmts * nq);
3184
3185 double epsilon = 1.0e-8;
3186 for (int i = 0; i < diff1.size(); ++i)
3187 {
3188 BOOST_CHECK_CLOSE(diff1[i], diff2[i], epsilon);
3189 }
3190}
3191
3192BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_Iterperexp_UniformP)
3193{
3194 SpatialDomains::PointGeomUniquePtr v0(
3195 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3196 SpatialDomains::PointGeomUniquePtr v1(
3197 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3198 SpatialDomains::PointGeomUniquePtr v2(
3199 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3200 SpatialDomains::PointGeomUniquePtr v3(
3201 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3202 SpatialDomains::PointGeomUniquePtr v4(
3203 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3204 SpatialDomains::PointGeomUniquePtr v5(
3205 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3206 SpatialDomains::PointGeomUniquePtr v6(
3207 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3208 SpatialDomains::PointGeomUniquePtr v7(
3209 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3210
3211 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3212 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3213 SpatialDomains::HexGeomUniquePtr hexGeom =
3214 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3215
3216 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3217 Nektar::LibUtilities::eGaussLobattoLegendre;
3218 Nektar::LibUtilities::BasisType basisTypeDir1 =
3219 Nektar::LibUtilities::eModified_A;
3220 unsigned int numQuadPoints = 6;
3221 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3222 quadPointsTypeDir1);
3223 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3224 quadPointsKeyDir1);
3225
3226 Nektar::LocalRegions::HexExpSharedPtr Exp =
3227 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3228 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3229
3230 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3231 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3232 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3233
3234 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3235 CollExp.push_back(Exp);
3236
3237 LibUtilities::SessionReaderSharedPtr dummySession;
3238 Collections::CollectionOptimisation colOpt(dummySession, 3,
3239 Collections::eIterPerExp);
3240 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3241 Collections::Collection c(CollExp, impTypes);
3242 c.Initialise(Collections::eIProductWRTDerivBase);
3243
3244 const int nq = Exp->GetTotPoints();
3245 const int nm = Exp->GetNcoeffs();
3246 Array<OneD, NekDouble> phys1(nq, 0.0);
3247 Array<OneD, NekDouble> phys2(nq, 0.0);
3248 Array<OneD, NekDouble> phys3(nq, 0.0);
3249 Array<OneD, NekDouble> coeffs1(nm, 0.0);
3250 Array<OneD, NekDouble> coeffs2(nm, 0.0);
3251
3252 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3253
3254 Exp->GetCoords(xc, yc, zc);
3255
3256 for (int i = 0; i < nq; ++i)
3257 {
3258 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3259 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3260 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3261 }
3262
3263 // Standard routines
3264 Exp->IProductWRTDerivBase(0, phys1, coeffs1);
3265 Exp->IProductWRTDerivBase(1, phys2, coeffs2);
3266 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
3267 Exp->IProductWRTDerivBase(2, phys3, coeffs2);
3268 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
3269
3270 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3271 coeffs2);
3272
3273 double epsilon = 1.0e-8;
3274 for (int i = 0; i < nm; ++i)
3275 {
3276 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
3277 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
3278 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
3279 }
3280}
3281
3282BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_MatrixFree_UniformP_Undeformed)
3283{
3284 SpatialDomains::PointGeomUniquePtr v0(
3285 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
3286 SpatialDomains::PointGeomUniquePtr v1(
3287 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3288 SpatialDomains::PointGeomUniquePtr v2(
3289 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3290 SpatialDomains::PointGeomUniquePtr v3(
3291 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3292 SpatialDomains::PointGeomUniquePtr v4(
3293 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3294 SpatialDomains::PointGeomUniquePtr v5(
3295 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3296 SpatialDomains::PointGeomUniquePtr v6(
3297 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3298 SpatialDomains::PointGeomUniquePtr v7(
3299 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3300
3301 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3302 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3303 SpatialDomains::HexGeomUniquePtr hexGeom =
3304 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3305
3306 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3307 Nektar::LibUtilities::eGaussLobattoLegendre;
3308 Nektar::LibUtilities::BasisType basisTypeDir1 =
3309 Nektar::LibUtilities::eModified_A;
3310 unsigned int numQuadPoints = 5;
3311 unsigned int numModes = 4;
3312 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3313 quadPointsTypeDir1);
3314 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
3315 quadPointsKeyDir1);
3316
3317 Nektar::LocalRegions::HexExpSharedPtr Exp =
3318 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3319 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3320
3321 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3322 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3323 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3324
3325 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3326 CollExp.push_back(Exp);
3327
3328 LibUtilities::SessionReaderSharedPtr dummySession;
3329 Collections::CollectionOptimisation colOpt(dummySession, 2,
3330 Collections::eMatrixFree);
3331 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3332 Collections::Collection c(CollExp, impTypes);
3333 c.Initialise(Collections::eIProductWRTDerivBase);
3334
3335 const int nq = Exp->GetTotPoints();
3336 const int nm = Exp->GetNcoeffs();
3337 Array<OneD, NekDouble> phys1(nq, 0.0);
3338 Array<OneD, NekDouble> phys2(nq, 0.0);
3339 Array<OneD, NekDouble> phys3(nq, 0.0);
3340 Array<OneD, NekDouble> coeffsRef(nm, 0.0);
3341 Array<OneD, NekDouble> coeffs(nm, 0.0);
3342
3343 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3344
3345 Exp->GetCoords(xc, yc, zc);
3346
3347 for (int i = 0; i < nq; ++i)
3348 {
3349 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3350 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3351 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3352 }
3353
3354 // Standard routines
3355 Exp->IProductWRTDerivBase(0, phys1, coeffsRef);
3356 Exp->IProductWRTDerivBase(1, phys2, coeffs);
3357 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3358 Exp->IProductWRTDerivBase(2, phys3, coeffs);
3359 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3360
3361 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3362 coeffs);
3363
3364 double epsilon = 1.0e-8;
3365 for (int i = 0; i < nm; ++i)
3366 {
3367 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
3368 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
3369 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
3370 }
3371}
3372
3373BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_MatrixFree_UniformP_Deformed)
3374{
3375 SpatialDomains::PointGeomUniquePtr v0(
3376 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
3377 SpatialDomains::PointGeomUniquePtr v1(
3378 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3379 SpatialDomains::PointGeomUniquePtr v2(
3380 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3381 SpatialDomains::PointGeomUniquePtr v3(
3382 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3383 SpatialDomains::PointGeomUniquePtr v4(
3384 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3385 SpatialDomains::PointGeomUniquePtr v5(
3386 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3387 SpatialDomains::PointGeomUniquePtr v6(
3388 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
3389 SpatialDomains::PointGeomUniquePtr v7(
3390 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3391
3392 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3393 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3394 SpatialDomains::HexGeomUniquePtr hexGeom =
3395 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3396
3397 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3398 Nektar::LibUtilities::eGaussLobattoLegendre;
3399 Nektar::LibUtilities::BasisType basisTypeDir1 =
3400 Nektar::LibUtilities::eModified_A;
3401 unsigned int numQuadPoints = 5;
3402 unsigned int numModes = 4;
3403 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3404 quadPointsTypeDir1);
3405 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
3406 quadPointsKeyDir1);
3407
3408 Nektar::LocalRegions::HexExpSharedPtr Exp =
3409 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3410 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3411
3412 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3413 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3414 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3415
3416 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3417 CollExp.push_back(Exp);
3418
3419 LibUtilities::SessionReaderSharedPtr dummySession;
3420 Collections::CollectionOptimisation colOpt(dummySession, 2,
3421 Collections::eMatrixFree);
3422 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3423 Collections::Collection c(CollExp, impTypes);
3424 c.Initialise(Collections::eIProductWRTDerivBase);
3425
3426 const int nq = Exp->GetTotPoints();
3427 const int nm = Exp->GetNcoeffs();
3428 Array<OneD, NekDouble> phys1(nq, 0.0);
3429 Array<OneD, NekDouble> phys2(nq, 0.0);
3430 Array<OneD, NekDouble> phys3(nq, 0.0);
3431 Array<OneD, NekDouble> coeffsRef(nm, 0.0);
3432 Array<OneD, NekDouble> coeffs(nm, 0.0);
3433
3434 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3435
3436 Exp->GetCoords(xc, yc, zc);
3437
3438 for (int i = 0; i < nq; ++i)
3439 {
3440 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3441 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3442 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3443 }
3444
3445 // Standard routines
3446 Exp->IProductWRTDerivBase(0, phys1, coeffsRef);
3447 Exp->IProductWRTDerivBase(1, phys2, coeffs);
3448 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3449 Exp->IProductWRTDerivBase(2, phys3, coeffs);
3450 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3451
3452 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3453 coeffs);
3454
3455 double epsilon = 1.0e-8;
3456 for (int i = 0; i < nm; ++i)
3457 {
3458 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
3459 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
3460 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
3461 }
3462}
3463
3464BOOST_AUTO_TEST_CASE(
3465 TestHexIProductWRTDerivBase_MatrixFree_UniformP_Deformed_OverInt)
3466{
3467 SpatialDomains::PointGeomUniquePtr v0(
3468 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
3469 SpatialDomains::PointGeomUniquePtr v1(
3470 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3471 SpatialDomains::PointGeomUniquePtr v2(
3472 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3473 SpatialDomains::PointGeomUniquePtr v3(
3474 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3475 SpatialDomains::PointGeomUniquePtr v4(
3476 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3477 SpatialDomains::PointGeomUniquePtr v5(
3478 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3479 SpatialDomains::PointGeomUniquePtr v6(
3480 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
3481 SpatialDomains::PointGeomUniquePtr v7(
3482 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3483
3484 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3485 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3486 SpatialDomains::HexGeomUniquePtr hexGeom =
3487 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3488
3489 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3490 Nektar::LibUtilities::eGaussLobattoLegendre;
3491 Nektar::LibUtilities::BasisType basisTypeDir1 =
3492 Nektar::LibUtilities::eModified_A;
3493 unsigned int numQuadPoints = 8;
3494 unsigned int numModes = 4;
3495 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3496 quadPointsTypeDir1);
3497 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
3498 quadPointsKeyDir1);
3499
3500 Nektar::LocalRegions::HexExpSharedPtr Exp =
3501 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3502 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3503
3504 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3505 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3506 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3507
3508 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3509 CollExp.push_back(Exp);
3510
3511 LibUtilities::SessionReaderSharedPtr dummySession;
3512 Collections::CollectionOptimisation colOpt(dummySession, 2,
3513 Collections::eMatrixFree);
3514 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3515 Collections::Collection c(CollExp, impTypes);
3516 c.Initialise(Collections::eIProductWRTDerivBase);
3517
3518 const int nq = Exp->GetTotPoints();
3519 const int nm = Exp->GetNcoeffs();
3520 Array<OneD, NekDouble> phys1(nq, 0.0);
3521 Array<OneD, NekDouble> phys2(nq, 0.0);
3522 Array<OneD, NekDouble> phys3(nq, 0.0);
3523 Array<OneD, NekDouble> coeffsRef(nm, 0.0);
3524 Array<OneD, NekDouble> coeffs(nm, 0.0);
3525
3526 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3527
3528 Exp->GetCoords(xc, yc, zc);
3529
3530 for (int i = 0; i < nq; ++i)
3531 {
3532 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3533 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3534 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3535 }
3536
3537 // Standard routines
3538 Exp->IProductWRTDerivBase(0, phys1, coeffsRef);
3539 Exp->IProductWRTDerivBase(1, phys2, coeffs);
3540 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3541 Exp->IProductWRTDerivBase(2, phys3, coeffs);
3542 Vmath::Vadd(nm, coeffsRef, 1, coeffs, 1, coeffsRef, 1);
3543
3544 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3545 coeffs);
3546
3547 double epsilon = 1.0e-8;
3548 for (int i = 0; i < nm; ++i)
3549 {
3550 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
3551 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
3552 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
3553 }
3554}
3555
3556BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_IterPerExp_VariableP_MultiElmt)
3557{
3558 SpatialDomains::PointGeomUniquePtr v0(
3559 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3560 SpatialDomains::PointGeomUniquePtr v1(
3561 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3562 SpatialDomains::PointGeomUniquePtr v2(
3563 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3564 SpatialDomains::PointGeomUniquePtr v3(
3565 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3566 SpatialDomains::PointGeomUniquePtr v4(
3567 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3568 SpatialDomains::PointGeomUniquePtr v5(
3569 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3570 SpatialDomains::PointGeomUniquePtr v6(
3571 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3572 SpatialDomains::PointGeomUniquePtr v7(
3573 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3574
3575 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3576 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3577 SpatialDomains::HexGeomUniquePtr hexGeom =
3578 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3579
3580 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3581 Nektar::LibUtilities::eGaussLobattoLegendre;
3582 Nektar::LibUtilities::BasisType basisTypeDir1 =
3583 Nektar::LibUtilities::eModified_A;
3584 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
3585 quadPointsTypeDir1);
3586 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
3587 quadPointsTypeDir1);
3588 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
3589 quadPointsTypeDir1);
3590 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3591 quadPointsKeyDir1);
3592 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
3593 quadPointsKeyDir2);
3594 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
3595 quadPointsKeyDir3);
3596
3597 Nektar::LocalRegions::HexExpSharedPtr Exp =
3598 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3599 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
3600
3601 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3602 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3603 basisKeyDir1, basisKeyDir2, basisKeyDir3);
3604
3605 int nelmts = 10;
3606
3607 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3608 for (int i = 0; i < nelmts; ++i)
3609 {
3610 CollExp.push_back(Exp);
3611 }
3612
3613 LibUtilities::SessionReaderSharedPtr dummySession;
3614 Collections::CollectionOptimisation colOpt(dummySession, 3,
3615 Collections::eIterPerExp);
3616 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3617 Collections::Collection c(CollExp, impTypes);
3618 c.Initialise(Collections::eIProductWRTDerivBase);
3619
3620 const int nq = Exp->GetTotPoints();
3621 const int nm = Exp->GetNcoeffs();
3622 Array<OneD, NekDouble> phys1(nelmts * nq, 0.0);
3623 Array<OneD, NekDouble> phys2(nelmts * nq, 0.0);
3624 Array<OneD, NekDouble> phys3(nelmts * nq, 0.0);
3625 Array<OneD, NekDouble> coeffs1(nelmts * nm, 0.0);
3626 Array<OneD, NekDouble> coeffs2(nelmts * nm, 0.0);
3627 Array<OneD, NekDouble> tmp;
3628
3629 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3630
3631 Exp->GetCoords(xc, yc, zc);
3632
3633 for (int i = 0; i < nq; ++i)
3634 {
3635 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3636 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3637 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3638 }
3639
3640 for (int i = 1; i < nelmts; ++i)
3641 {
3642 Vmath::Vcopy(nq, phys1, 1, tmp = phys1 + i * nq, 1);
3643 Vmath::Vcopy(nq, phys2, 1, tmp = phys2 + i * nq, 1);
3644 Vmath::Vcopy(nq, phys3, 1, tmp = phys3 + i * nq, 1);
3645 }
3646
3647 // Standard routines
3648 for (int i = 0; i < nelmts; ++i)
3649 {
3650
3651 Exp->IProductWRTDerivBase(0, phys1 + i * nq, tmp = coeffs1 + i * nm);
3652 Exp->IProductWRTDerivBase(1, phys2 + i * nq, tmp = coeffs2 + i * nm);
3653 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3654 tmp = coeffs1 + i * nm, 1);
3655 Exp->IProductWRTDerivBase(2, phys3 + i * nq, tmp = coeffs2 + i * nm);
3656 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3657 tmp = coeffs1 + i * nm, 1);
3658 }
3659
3660 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3661 coeffs2);
3662
3663 double epsilon = 1.0e-6;
3664 for (int i = 0; i < coeffs1.size(); ++i)
3665 {
3666 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
3667 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
3668 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
3669 }
3670}
3671
3672BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_StdMat_UniformP)
3673{
3674 SpatialDomains::PointGeomUniquePtr v0(
3675 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3676 SpatialDomains::PointGeomUniquePtr v1(
3677 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3678 SpatialDomains::PointGeomUniquePtr v2(
3679 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3680 SpatialDomains::PointGeomUniquePtr v3(
3681 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3682 SpatialDomains::PointGeomUniquePtr v4(
3683 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3684 SpatialDomains::PointGeomUniquePtr v5(
3685 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3686 SpatialDomains::PointGeomUniquePtr v6(
3687 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3688 SpatialDomains::PointGeomUniquePtr v7(
3689 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3690
3691 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3692 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3693 SpatialDomains::HexGeomUniquePtr hexGeom =
3694 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3695
3696 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3697 Nektar::LibUtilities::eGaussLobattoLegendre;
3698 Nektar::LibUtilities::BasisType basisTypeDir1 =
3699 Nektar::LibUtilities::eModified_A;
3700 unsigned int numQuadPoints = 6;
3701 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
3702 quadPointsTypeDir1);
3703 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3704 quadPointsKeyDir1);
3705
3706 Nektar::LocalRegions::HexExpSharedPtr Exp =
3707 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3708 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
3709
3710 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3711 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3712 basisKeyDir1, basisKeyDir1, basisKeyDir1);
3713
3714 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3715 CollExp.push_back(Exp);
3716
3717 LibUtilities::SessionReaderSharedPtr dummySession;
3718 Collections::CollectionOptimisation colOpt(dummySession, 3,
3719 Collections::eStdMat);
3720 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3721 Collections::Collection c(CollExp, impTypes);
3722 c.Initialise(Collections::eIProductWRTDerivBase);
3723
3724 const int nq = Exp->GetTotPoints();
3725 const int nm = Exp->GetNcoeffs();
3726 Array<OneD, NekDouble> phys1(nq, 0.0);
3727 Array<OneD, NekDouble> phys2(nq, 0.0);
3728 Array<OneD, NekDouble> phys3(nq, 0.0);
3729 Array<OneD, NekDouble> coeffs1(nm, 0.0);
3730 Array<OneD, NekDouble> coeffs2(nm, 0.0);
3731
3732 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3733
3734 Exp->GetCoords(xc, yc, zc);
3735
3736 for (int i = 0; i < nq; ++i)
3737 {
3738 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3739 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3740 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3741 }
3742
3743 // Standard routines
3744 Exp->IProductWRTDerivBase(0, phys1, coeffs1);
3745 Exp->IProductWRTDerivBase(1, phys2, coeffs2);
3746 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
3747 Exp->IProductWRTDerivBase(2, phys3, coeffs2);
3748 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
3749
3750 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3751 coeffs2);
3752
3753 double epsilon = 1.0e-8;
3754 for (int i = 0; i < coeffs1.size(); ++i)
3755 {
3756 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
3757 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
3758 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
3759 }
3760}
3761
3762BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_StdMat_VariableP_MultiElmt)
3763{
3764 SpatialDomains::PointGeomUniquePtr v0(
3765 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3766 SpatialDomains::PointGeomUniquePtr v1(
3767 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3768 SpatialDomains::PointGeomUniquePtr v2(
3769 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3770 SpatialDomains::PointGeomUniquePtr v3(
3771 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3772 SpatialDomains::PointGeomUniquePtr v4(
3773 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3774 SpatialDomains::PointGeomUniquePtr v5(
3775 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3776 SpatialDomains::PointGeomUniquePtr v6(
3777 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3778 SpatialDomains::PointGeomUniquePtr v7(
3779 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3780
3781 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3782 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3783 SpatialDomains::HexGeomUniquePtr hexGeom =
3784 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3785
3786 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3787 Nektar::LibUtilities::eGaussLobattoLegendre;
3788 Nektar::LibUtilities::BasisType basisTypeDir1 =
3789 Nektar::LibUtilities::eModified_A;
3790 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
3791 quadPointsTypeDir1);
3792 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
3793 quadPointsTypeDir1);
3794 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
3795 quadPointsTypeDir1);
3796 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3797 quadPointsKeyDir1);
3798 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
3799 quadPointsKeyDir2);
3800 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
3801 quadPointsKeyDir3);
3802
3803 Nektar::LocalRegions::HexExpSharedPtr Exp =
3804 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3805 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
3806
3807 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3808 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3809 basisKeyDir1, basisKeyDir2, basisKeyDir3);
3810
3811 int nelmts = 10;
3812
3813 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3814 for (int i = 0; i < nelmts; ++i)
3815 {
3816 CollExp.push_back(Exp);
3817 }
3818
3819 LibUtilities::SessionReaderSharedPtr dummySession;
3820 Collections::CollectionOptimisation colOpt(dummySession, 3,
3821 Collections::eStdMat);
3822 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3823 Collections::Collection c(CollExp, impTypes);
3824 c.Initialise(Collections::eIProductWRTDerivBase);
3825
3826 const int nq = Exp->GetTotPoints();
3827 const int nm = Exp->GetNcoeffs();
3828 Array<OneD, NekDouble> phys1(nelmts * nq, 0.0);
3829 Array<OneD, NekDouble> phys2(nelmts * nq, 0.0);
3830 Array<OneD, NekDouble> phys3(nelmts * nq, 0.0);
3831 Array<OneD, NekDouble> coeffs1(nelmts * nm, 0.0);
3832 Array<OneD, NekDouble> coeffs2(nelmts * nm, 0.0);
3833 Array<OneD, NekDouble> tmp;
3834
3835 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3836
3837 Exp->GetCoords(xc, yc, zc);
3838
3839 for (int i = 0; i < nq; ++i)
3840 {
3841 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3842 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3843 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3844 }
3845
3846 for (int i = 1; i < nelmts; ++i)
3847 {
3848 Vmath::Vcopy(nq, phys1, 1, tmp = phys1 + i * nq, 1);
3849 Vmath::Vcopy(nq, phys2, 1, tmp = phys2 + i * nq, 1);
3850 Vmath::Vcopy(nq, phys3, 1, tmp = phys3 + i * nq, 1);
3851 }
3852
3853 // Standard routines
3854 for (int i = 0; i < nelmts; ++i)
3855 {
3856
3857 Exp->IProductWRTDerivBase(0, phys1 + i * nq, tmp = coeffs1 + i * nm);
3858 Exp->IProductWRTDerivBase(1, phys2 + i * nq, tmp = coeffs2 + i * nm);
3859 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3860 tmp = coeffs1 + i * nm, 1);
3861 Exp->IProductWRTDerivBase(2, phys3 + i * nq, tmp = coeffs2 + i * nm);
3862 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3863 tmp = coeffs1 + i * nm, 1);
3864 }
3865
3866 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3867 coeffs2);
3868
3869 double epsilon = 1.0e-8;
3870 for (int i = 0; i < coeffs1.size(); ++i)
3871 {
3872 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
3873 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
3874 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
3875 }
3876}
3877
3878BOOST_AUTO_TEST_CASE(
3879 TestHexIProductWRTDerivBase_NoCollection_VariableP_MultiElmt)
3880{
3881 SpatialDomains::PointGeomUniquePtr v0(
3882 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3883 SpatialDomains::PointGeomUniquePtr v1(
3884 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
3885 SpatialDomains::PointGeomUniquePtr v2(
3886 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
3887 SpatialDomains::PointGeomUniquePtr v3(
3888 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
3889 SpatialDomains::PointGeomUniquePtr v4(
3890 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
3891 SpatialDomains::PointGeomUniquePtr v5(
3892 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
3893 SpatialDomains::PointGeomUniquePtr v6(
3894 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
3895 SpatialDomains::PointGeomUniquePtr v7(
3896 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
3897
3898 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
3899 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
3900 SpatialDomains::HexGeomUniquePtr hexGeom =
3901 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
3902
3903 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
3904 Nektar::LibUtilities::eGaussLobattoLegendre;
3905 Nektar::LibUtilities::BasisType basisTypeDir1 =
3906 Nektar::LibUtilities::eModified_A;
3907 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
3908 quadPointsTypeDir1);
3909 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
3910 quadPointsTypeDir1);
3911 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
3912 quadPointsTypeDir1);
3913 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
3914 quadPointsKeyDir1);
3915 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
3916 quadPointsKeyDir2);
3917 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
3918 quadPointsKeyDir3);
3919
3920 Nektar::LocalRegions::HexExpSharedPtr Exp =
3921 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
3922 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
3923
3924 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
3925 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
3926 basisKeyDir1, basisKeyDir2, basisKeyDir3);
3927
3928 int nelmts = 10;
3929
3930 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
3931 for (int i = 0; i < nelmts; ++i)
3932 {
3933 CollExp.push_back(Exp);
3934 }
3935
3936 LibUtilities::SessionReaderSharedPtr dummySession;
3937 Collections::CollectionOptimisation colOpt(dummySession, 3,
3938 Collections::eNoCollection);
3939 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
3940 Collections::Collection c(CollExp, impTypes);
3941 c.Initialise(Collections::eIProductWRTDerivBase);
3942
3943 const int nq = Exp->GetTotPoints();
3944 const int nm = Exp->GetNcoeffs();
3945 Array<OneD, NekDouble> phys1(nelmts * nq, 0.0);
3946 Array<OneD, NekDouble> phys2(nelmts * nq, 0.0);
3947 Array<OneD, NekDouble> phys3(nelmts * nq, 0.0);
3948 Array<OneD, NekDouble> coeffs1(nelmts * nm, 0.0);
3949 Array<OneD, NekDouble> coeffs2(nelmts * nm, 0.0);
3950 Array<OneD, NekDouble> tmp;
3951 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
3952
3953 Exp->GetCoords(xc, yc, zc);
3954
3955 for (int i = 0; i < nq; ++i)
3956 {
3957 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
3958 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
3959 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
3960 }
3961
3962 for (int i = 1; i < nelmts; ++i)
3963 {
3964 Vmath::Vcopy(nq, phys1, 1, tmp = phys1 + i * nq, 1);
3965 Vmath::Vcopy(nq, phys2, 1, tmp = phys2 + i * nq, 1);
3966 Vmath::Vcopy(nq, phys3, 1, tmp = phys3 + i * nq, 1);
3967 }
3968
3969 // Standard routines
3970 for (int i = 0; i < nelmts; ++i)
3971 {
3972
3973 Exp->IProductWRTDerivBase(0, phys1 + i * nq, tmp = coeffs1 + i * nm);
3974 Exp->IProductWRTDerivBase(1, phys2 + i * nq, tmp = coeffs2 + i * nm);
3975 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3976 tmp = coeffs1 + i * nm, 1);
3977 Exp->IProductWRTDerivBase(2, phys3 + i * nq, tmp = coeffs2 + i * nm);
3978 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
3979 tmp = coeffs1 + i * nm, 1);
3980 }
3981
3982 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
3983 coeffs2);
3984
3985 double epsilon = 1.0e-6;
3986 for (int i = 0; i < coeffs1.size(); ++i)
3987 {
3988 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
3989 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
3990 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
3991 }
3992}
3993
3994BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_SumFac_UniformP)
3995{
3996 SpatialDomains::PointGeomUniquePtr v0(
3997 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
3998 SpatialDomains::PointGeomUniquePtr v1(
3999 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4000 SpatialDomains::PointGeomUniquePtr v2(
4001 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4002 SpatialDomains::PointGeomUniquePtr v3(
4003 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4004 SpatialDomains::PointGeomUniquePtr v4(
4005 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4006 SpatialDomains::PointGeomUniquePtr v5(
4007 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4008 SpatialDomains::PointGeomUniquePtr v6(
4009 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
4010 SpatialDomains::PointGeomUniquePtr v7(
4011 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4012
4013 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4014 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4015 SpatialDomains::HexGeomUniquePtr hexGeom =
4016 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4017
4018 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4019 Nektar::LibUtilities::eGaussLobattoLegendre;
4020 Nektar::LibUtilities::BasisType basisTypeDir1 =
4021 Nektar::LibUtilities::eModified_A;
4022 unsigned int numQuadPoints = 6;
4023 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4024 quadPointsTypeDir1);
4025 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
4026 quadPointsKeyDir1);
4027
4028 Nektar::LocalRegions::HexExpSharedPtr Exp =
4029 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4030 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4031
4032 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4033 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4034 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4035
4036 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4037 CollExp.push_back(Exp);
4038
4039 LibUtilities::SessionReaderSharedPtr dummySession;
4040 Collections::CollectionOptimisation colOpt(dummySession, 3,
4041 Collections::eSumFac);
4042 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4043 Collections::Collection c(CollExp, impTypes);
4044 c.Initialise(Collections::eIProductWRTDerivBase);
4045
4046 const int nq = Exp->GetTotPoints();
4047 const int nm = Exp->GetNcoeffs();
4048 Array<OneD, NekDouble> phys1(nq, 0.0);
4049 Array<OneD, NekDouble> phys2(nq, 0.0);
4050 Array<OneD, NekDouble> phys3(nq, 0.0);
4051 Array<OneD, NekDouble> coeffs1(nm, 0.0);
4052 Array<OneD, NekDouble> coeffs2(nm, 0.0);
4053
4054 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
4055
4056 Exp->GetCoords(xc, yc, zc);
4057
4058 for (int i = 0; i < nq; ++i)
4059 {
4060 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
4061 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
4062 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
4063 }
4064
4065 // Standard routines
4066 Exp->IProductWRTDerivBase(0, phys1, coeffs1);
4067 Exp->IProductWRTDerivBase(1, phys2, coeffs2);
4068 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
4069 Exp->IProductWRTDerivBase(2, phys3, coeffs2);
4070 Vmath::Vadd(nm, coeffs1, 1, coeffs2, 1, coeffs1, 1);
4071
4072 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
4073 coeffs2);
4074
4075 double epsilon = 1.0e-8;
4076 for (int i = 0; i < coeffs1.size(); ++i)
4077 {
4078 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
4079 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
4080 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
4081 }
4082}
4083
4084BOOST_AUTO_TEST_CASE(TestHexIProductWRTDerivBase_SumFac_VariableP_MultiElmt)
4085{
4086 SpatialDomains::PointGeomUniquePtr v0(
4087 new SpatialDomains::PointGeom(3u, 0u, -1.5, -1.5, -1.5));
4088 SpatialDomains::PointGeomUniquePtr v1(
4089 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4090 SpatialDomains::PointGeomUniquePtr v2(
4091 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4092 SpatialDomains::PointGeomUniquePtr v3(
4093 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4094 SpatialDomains::PointGeomUniquePtr v4(
4095 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4096 SpatialDomains::PointGeomUniquePtr v5(
4097 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4098 SpatialDomains::PointGeomUniquePtr v6(
4099 new SpatialDomains::PointGeom(3u, 6u, 1.0, 1.0, 1.0));
4100 SpatialDomains::PointGeomUniquePtr v7(
4101 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4102
4103 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4104 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4105 SpatialDomains::HexGeomUniquePtr hexGeom =
4106 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4107
4108 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4109 Nektar::LibUtilities::eGaussLobattoLegendre;
4110 Nektar::LibUtilities::BasisType basisTypeDir1 =
4111 Nektar::LibUtilities::eModified_A;
4112 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(5,
4113 quadPointsTypeDir1);
4114 const Nektar::LibUtilities::PointsKey quadPointsKeyDir2(6,
4115 quadPointsTypeDir1);
4116 const Nektar::LibUtilities::PointsKey quadPointsKeyDir3(8,
4117 quadPointsTypeDir1);
4118 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, 4,
4119 quadPointsKeyDir1);
4120 const Nektar::LibUtilities::BasisKey basisKeyDir2(basisTypeDir1, 6,
4121 quadPointsKeyDir2);
4122 const Nektar::LibUtilities::BasisKey basisKeyDir3(basisTypeDir1, 8,
4123 quadPointsKeyDir3);
4124
4125 Nektar::LocalRegions::HexExpSharedPtr Exp =
4126 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4127 basisKeyDir1, basisKeyDir2, basisKeyDir3, hexGeom.get());
4128
4129 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4130 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4131 basisKeyDir1, basisKeyDir2, basisKeyDir3);
4132
4133 int nelmts = 10;
4134
4135 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4136 for (int i = 0; i < nelmts; ++i)
4137 {
4138 CollExp.push_back(Exp);
4139 }
4140
4141 LibUtilities::SessionReaderSharedPtr dummySession;
4142 Collections::CollectionOptimisation colOpt(dummySession, 3,
4143 Collections::eSumFac);
4144 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4145 Collections::Collection c(CollExp, impTypes);
4146 c.Initialise(Collections::eIProductWRTDerivBase);
4147
4148 const int nq = Exp->GetTotPoints();
4149 const int nm = Exp->GetNcoeffs();
4150 Array<OneD, NekDouble> phys1(nelmts * nq, 0.0);
4151 Array<OneD, NekDouble> phys2(nelmts * nq, 0.0);
4152 Array<OneD, NekDouble> phys3(nelmts * nq, 0.0);
4153 Array<OneD, NekDouble> coeffs1(nelmts * nm, 0.0);
4154 Array<OneD, NekDouble> coeffs2(nelmts * nm, 0.0);
4155 Array<OneD, NekDouble> tmp;
4156
4157 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
4158
4159 Exp->GetCoords(xc, yc, zc);
4160
4161 for (int i = 0; i < nq; ++i)
4162 {
4163 phys1[i] = sin(xc[i]) * cos(yc[i]) * sin(zc[i]);
4164 phys2[i] = cos(xc[i]) * sin(yc[i]) * cos(zc[i]);
4165 phys2[i] = cos(xc[i]) * sin(yc[i]) * sin(zc[i]);
4166 }
4167
4168 for (int i = 1; i < nelmts; ++i)
4169 {
4170 Vmath::Vcopy(nq, phys1, 1, tmp = phys1 + i * nq, 1);
4171 Vmath::Vcopy(nq, phys2, 1, tmp = phys2 + i * nq, 1);
4172 Vmath::Vcopy(nq, phys3, 1, tmp = phys3 + i * nq, 1);
4173 }
4174
4175 // Standard routines
4176 for (int i = 0; i < nelmts; ++i)
4177 {
4178 Exp->IProductWRTDerivBase(0, phys1 + i * nq, tmp = coeffs1 + i * nm);
4179 Exp->IProductWRTDerivBase(1, phys2 + i * nq, tmp = coeffs2 + i * nm);
4180 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
4181 tmp = coeffs1 + i * nm, 1);
4182 Exp->IProductWRTDerivBase(2, phys3 + i * nq, tmp = coeffs2 + i * nm);
4183 Vmath::Vadd(nm, coeffs1 + i * nm, 1, coeffs2 + i * nm, 1,
4184 tmp = coeffs1 + i * nm, 1);
4185 }
4186
4187 c.ApplyOperator(Collections::eIProductWRTDerivBase, phys1, phys2, phys3,
4188 coeffs2);
4189
4190 double epsilon = 1.0e-8;
4191 for (int i = 0; i < coeffs1.size(); ++i)
4192 {
4193 coeffs1[i] = (std::abs(coeffs1[i]) < 1e-14) ? 0.0 : coeffs1[i];
4194 coeffs2[i] = (std::abs(coeffs2[i]) < 1e-14) ? 0.0 : coeffs2[i];
4195 BOOST_CHECK_CLOSE(coeffs1[i], coeffs2[i], epsilon);
4196 }
4197}
4198
4199BOOST_AUTO_TEST_CASE(TestHexHelmholtz_NoCollection_UniformP)
4200{
4201 SpatialDomains::PointGeomUniquePtr v0(
4202 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4203 SpatialDomains::PointGeomUniquePtr v1(
4204 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4205 SpatialDomains::PointGeomUniquePtr v2(
4206 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4207 SpatialDomains::PointGeomUniquePtr v3(
4208 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4209 SpatialDomains::PointGeomUniquePtr v4(
4210 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4211 SpatialDomains::PointGeomUniquePtr v5(
4212 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4213 SpatialDomains::PointGeomUniquePtr v6(
4214 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4215 SpatialDomains::PointGeomUniquePtr v7(
4216 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4217
4218 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4219 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4220 SpatialDomains::HexGeomUniquePtr hexGeom =
4221 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4222
4223 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4224 Nektar::LibUtilities::eGaussLobattoLegendre;
4225 Nektar::LibUtilities::BasisType basisTypeDir1 =
4226 Nektar::LibUtilities::eModified_A;
4227 unsigned int numQuadPoints = 5;
4228 unsigned int numModes = 4;
4229
4230 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4231 quadPointsTypeDir1);
4232 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4233 quadPointsKeyDir1);
4234
4235 Nektar::LocalRegions::HexExpSharedPtr Exp =
4236 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4237 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4238
4239 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4240 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4241 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4242
4243 int nelmts = 10;
4244
4245 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4246 for (int i = 0; i < nelmts; ++i)
4247 {
4248 CollExp.push_back(Exp);
4249 }
4250
4251 LibUtilities::SessionReaderSharedPtr dummySession;
4252 Collections::CollectionOptimisation colOpt(dummySession, 2,
4253 Collections::eNoCollection);
4254 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4255 Collections::Collection c(CollExp, impTypes);
4256 StdRegions::ConstFactorMap factors;
4257 factors[StdRegions::eFactorLambda] = 1.5;
4258
4259 c.Initialise(Collections::eHelmholtz, factors);
4260
4261 const int nm = Exp->GetNcoeffs();
4262 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4263 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4264 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4265
4266 for (int i = 0; i < nm; ++i)
4267 {
4268 coeffsIn[i] = 1.0;
4269 }
4270
4271 for (int i = 1; i < nelmts; ++i)
4272 {
4273 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4274 }
4275
4276 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4277 *Exp, factors);
4278
4279 for (int i = 0; i < nelmts; ++i)
4280 {
4281 // Standard routines
4282 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4283 }
4284
4285 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4286
4287 double epsilon = 1.0e-8;
4288 for (int i = 0; i < coeffsRef.size(); ++i)
4289 {
4290 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4291 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4292 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4293 }
4294}
4295
4296BOOST_AUTO_TEST_CASE(TestHexHelmholtz_IterPerExp_UniformP)
4297{
4298 SpatialDomains::PointGeomUniquePtr v0(
4299 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4300 SpatialDomains::PointGeomUniquePtr v1(
4301 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4302 SpatialDomains::PointGeomUniquePtr v2(
4303 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4304 SpatialDomains::PointGeomUniquePtr v3(
4305 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4306 SpatialDomains::PointGeomUniquePtr v4(
4307 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4308 SpatialDomains::PointGeomUniquePtr v5(
4309 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4310 SpatialDomains::PointGeomUniquePtr v6(
4311 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4312 SpatialDomains::PointGeomUniquePtr v7(
4313 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4314
4315 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4316 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4317 SpatialDomains::HexGeomUniquePtr hexGeom =
4318 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4319
4320 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4321 Nektar::LibUtilities::eGaussLobattoLegendre;
4322 Nektar::LibUtilities::BasisType basisTypeDir1 =
4323 Nektar::LibUtilities::eModified_A;
4324 unsigned int numQuadPoints = 5;
4325 unsigned int numModes = 4;
4326
4327 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4328 quadPointsTypeDir1);
4329 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4330 quadPointsKeyDir1);
4331
4332 Nektar::LocalRegions::HexExpSharedPtr Exp =
4333 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4334 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4335
4336 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4337 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4338 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4339
4340 int nelmts = 10;
4341
4342 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4343 for (int i = 0; i < nelmts; ++i)
4344 {
4345 CollExp.push_back(Exp);
4346 }
4347
4348 LibUtilities::SessionReaderSharedPtr dummySession;
4349 Collections::CollectionOptimisation colOpt(dummySession, 2,
4350 Collections::eIterPerExp);
4351 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4352 Collections::Collection c(CollExp, impTypes);
4353 StdRegions::ConstFactorMap factors;
4354 factors[StdRegions::eFactorLambda] = 1.5;
4355
4356 c.Initialise(Collections::eHelmholtz, factors);
4357
4358 const int nm = Exp->GetNcoeffs();
4359 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4360 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4361 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4362
4363 for (int i = 0; i < nm; ++i)
4364 {
4365 coeffsIn[i] = 1.0;
4366 }
4367
4368 for (int i = 1; i < nelmts; ++i)
4369 {
4370 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4371 }
4372
4373 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4374 *Exp, factors);
4375
4376 for (int i = 0; i < nelmts; ++i)
4377 {
4378 // Standard routines
4379 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4380 }
4381
4382 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4383
4384 double epsilon = 1.0e-8;
4385 for (int i = 0; i < coeffsRef.size(); ++i)
4386 {
4387 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4388 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4389 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4390 }
4391}
4392
4393BOOST_AUTO_TEST_CASE(TestHexHelmholtz_IterPerExp_UniformP_ConstVarDiff)
4394{
4395 SpatialDomains::PointGeomUniquePtr v0(
4396 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4397 SpatialDomains::PointGeomUniquePtr v1(
4398 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4399 SpatialDomains::PointGeomUniquePtr v2(
4400 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4401 SpatialDomains::PointGeomUniquePtr v3(
4402 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4403 SpatialDomains::PointGeomUniquePtr v4(
4404 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4405 SpatialDomains::PointGeomUniquePtr v5(
4406 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4407 SpatialDomains::PointGeomUniquePtr v6(
4408 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4409 SpatialDomains::PointGeomUniquePtr v7(
4410 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4411
4412 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4413 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4414 SpatialDomains::HexGeomUniquePtr hexGeom =
4415 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4416
4417 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4418 Nektar::LibUtilities::eGaussLobattoLegendre;
4419 Nektar::LibUtilities::BasisType basisTypeDir1 =
4420 Nektar::LibUtilities::eModified_A;
4421 unsigned int numQuadPoints = 5;
4422 unsigned int numModes = 4;
4423
4424 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4425 quadPointsTypeDir1);
4426 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4427 quadPointsKeyDir1);
4428
4429 Nektar::LocalRegions::HexExpSharedPtr Exp =
4430 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4431 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4432
4433 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4434 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4435 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4436
4437 int nelmts = 10;
4438
4439 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4440 for (int i = 0; i < nelmts; ++i)
4441 {
4442 CollExp.push_back(Exp);
4443 }
4444
4445 LibUtilities::SessionReaderSharedPtr dummySession;
4446 Collections::CollectionOptimisation colOpt(dummySession, 2,
4447 Collections::eIterPerExp);
4448 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4449 Collections::Collection c(CollExp, impTypes);
4450 StdRegions::ConstFactorMap factors;
4451 factors[StdRegions::eFactorLambda] = 1.5;
4452 factors[StdRegions::eFactorCoeffD00] = 1.25;
4453 factors[StdRegions::eFactorCoeffD01] = 0.25;
4454 factors[StdRegions::eFactorCoeffD11] = 1.25;
4455 factors[StdRegions::eFactorCoeffD02] = 0.25;
4456 factors[StdRegions::eFactorCoeffD12] = 0.25;
4457 factors[StdRegions::eFactorCoeffD22] = 1.25;
4458
4459 c.Initialise(Collections::eHelmholtz, factors);
4460
4461 const int nm = Exp->GetNcoeffs();
4462 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4463 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4464 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4465
4466 for (int i = 0; i < nm; ++i)
4467 {
4468 coeffsIn[i] = 1.0;
4469 }
4470
4471 for (int i = 1; i < nelmts; ++i)
4472 {
4473 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4474 }
4475
4476 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4477 *Exp, factors);
4478
4479 for (int i = 0; i < nelmts; ++i)
4480 {
4481 // Standard routines
4482 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4483 }
4484
4485 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4486
4487 double epsilon = 1.0e-8;
4488 for (int i = 0; i < coeffsRef.size(); ++i)
4489 {
4490 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4491 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4492 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4493 }
4494}
4495
4496BOOST_AUTO_TEST_CASE(TestHexHelmholtz_MatrixFree_UniformP)
4497{
4498 SpatialDomains::PointGeomUniquePtr v0(
4499 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4500 SpatialDomains::PointGeomUniquePtr v1(
4501 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4502 SpatialDomains::PointGeomUniquePtr v2(
4503 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4504 SpatialDomains::PointGeomUniquePtr v3(
4505 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4506 SpatialDomains::PointGeomUniquePtr v4(
4507 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4508 SpatialDomains::PointGeomUniquePtr v5(
4509 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4510 SpatialDomains::PointGeomUniquePtr v6(
4511 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4512 SpatialDomains::PointGeomUniquePtr v7(
4513 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4514
4515 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4516 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4517 SpatialDomains::HexGeomUniquePtr hexGeom =
4518 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4519
4520 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4521 Nektar::LibUtilities::eGaussLobattoLegendre;
4522 Nektar::LibUtilities::BasisType basisTypeDir1 =
4523 Nektar::LibUtilities::eModified_A;
4524 unsigned int numQuadPoints = 5;
4525 unsigned int numModes = 4;
4526
4527 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4528 quadPointsTypeDir1);
4529 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4530 quadPointsKeyDir1);
4531
4532 Nektar::LocalRegions::HexExpSharedPtr Exp =
4533 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4534 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4535
4536 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4537 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4538 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4539
4540 int nelmts = 10;
4541
4542 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4543 for (int i = 0; i < nelmts; ++i)
4544 {
4545 CollExp.push_back(Exp);
4546 }
4547
4548 LibUtilities::SessionReaderSharedPtr dummySession;
4549 Collections::CollectionOptimisation colOpt(dummySession, 2,
4550 Collections::eMatrixFree);
4551 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4552 Collections::Collection c(CollExp, impTypes);
4553 StdRegions::ConstFactorMap factors;
4554 factors[StdRegions::eFactorLambda] = 0.0;
4555
4556 c.Initialise(Collections::eHelmholtz, factors);
4557
4558 const int nm = Exp->GetNcoeffs();
4559 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4560 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4561 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4562
4563 for (int i = 0; i < nm; ++i)
4564 {
4565 coeffsIn[i] = 1.0;
4566 }
4567
4568 for (int i = 1; i < nelmts; ++i)
4569 {
4570 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4571 }
4572
4573 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4574 *Exp, factors);
4575
4576 for (int i = 0; i < nelmts; ++i)
4577 {
4578 // Standard routines
4579 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4580 }
4581
4582 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4583
4584 double epsilon = 1.0e-8;
4585 for (int i = 0; i < coeffsRef.size(); ++i)
4586 {
4587 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4588 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4589 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4590 }
4591}
4592
4593BOOST_AUTO_TEST_CASE(TestHexHelmholtz_MatrixFree_UniformP_Deformed_OverInt)
4594{
4595 SpatialDomains::PointGeomUniquePtr v0(
4596 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4597 SpatialDomains::PointGeomUniquePtr v1(
4598 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4599 SpatialDomains::PointGeomUniquePtr v2(
4600 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4601 SpatialDomains::PointGeomUniquePtr v3(
4602 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4603 SpatialDomains::PointGeomUniquePtr v4(
4604 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4605 SpatialDomains::PointGeomUniquePtr v5(
4606 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4607 SpatialDomains::PointGeomUniquePtr v6(
4608 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4609 SpatialDomains::PointGeomUniquePtr v7(
4610 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4611
4612 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4613 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4614 SpatialDomains::HexGeomUniquePtr hexGeom =
4615 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4616
4617 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4618 Nektar::LibUtilities::eGaussLobattoLegendre;
4619 Nektar::LibUtilities::BasisType basisTypeDir1 =
4620 Nektar::LibUtilities::eModified_A;
4621 unsigned int numQuadPoints = 8;
4622 unsigned int numModes = 4;
4623
4624 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4625 quadPointsTypeDir1);
4626 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4627 quadPointsKeyDir1);
4628
4629 Nektar::LocalRegions::HexExpSharedPtr Exp =
4630 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4631 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4632
4633 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4634 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4635 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4636
4637 int nelmts = 10;
4638
4639 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4640 for (int i = 0; i < nelmts; ++i)
4641 {
4642 CollExp.push_back(Exp);
4643 }
4644
4645 LibUtilities::SessionReaderSharedPtr dummySession;
4646 Collections::CollectionOptimisation colOpt(dummySession, 2,
4647 Collections::eMatrixFree);
4648 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4649 Collections::Collection c(CollExp, impTypes);
4650 StdRegions::ConstFactorMap factors;
4651 factors[StdRegions::eFactorLambda] = 0.0;
4652
4653 c.Initialise(Collections::eHelmholtz, factors);
4654
4655 const int nm = Exp->GetNcoeffs();
4656 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4657 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4658 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4659
4660 for (int i = 0; i < nm; ++i)
4661 {
4662 coeffsIn[i] = 1.0;
4663 }
4664
4665 for (int i = 1; i < nelmts; ++i)
4666 {
4667 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4668 }
4669
4670 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4671 *Exp, factors);
4672
4673 for (int i = 0; i < nelmts; ++i)
4674 {
4675 // Standard routines
4676 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4677 }
4678
4679 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4680
4681 double epsilon = 1.0e-8;
4682 for (int i = 0; i < coeffsRef.size(); ++i)
4683 {
4684 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4685 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4686 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4687 }
4688}
4689
4690BOOST_AUTO_TEST_CASE(TestHexHelmholtz_MatrixFree_UniformP_ConstVarDiff)
4691{
4692 SpatialDomains::PointGeomUniquePtr v0(
4693 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4694 SpatialDomains::PointGeomUniquePtr v1(
4695 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4696 SpatialDomains::PointGeomUniquePtr v2(
4697 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4698 SpatialDomains::PointGeomUniquePtr v3(
4699 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4700 SpatialDomains::PointGeomUniquePtr v4(
4701 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4702 SpatialDomains::PointGeomUniquePtr v5(
4703 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4704 SpatialDomains::PointGeomUniquePtr v6(
4705 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4706 SpatialDomains::PointGeomUniquePtr v7(
4707 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4708
4709 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4710 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4711 SpatialDomains::HexGeomUniquePtr hexGeom =
4712 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4713
4714 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4715 Nektar::LibUtilities::eGaussLobattoLegendre;
4716 Nektar::LibUtilities::BasisType basisTypeDir1 =
4717 Nektar::LibUtilities::eModified_A;
4718 unsigned int numQuadPoints = 5;
4719 unsigned int numModes = 4;
4720
4721 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4722 quadPointsTypeDir1);
4723 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4724 quadPointsKeyDir1);
4725
4726 Nektar::LocalRegions::HexExpSharedPtr Exp =
4727 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4728 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4729
4730 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4731 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4732 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4733
4734 int nelmts = 10;
4735
4736 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4737 for (int i = 0; i < nelmts; ++i)
4738 {
4739 CollExp.push_back(Exp);
4740 }
4741
4742 LibUtilities::SessionReaderSharedPtr dummySession;
4743 Collections::CollectionOptimisation colOpt(dummySession, 2,
4744 Collections::eMatrixFree);
4745 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(stdExp);
4746 Collections::Collection c(CollExp, impTypes);
4747 StdRegions::ConstFactorMap factors;
4748 factors[StdRegions::eFactorLambda] = 1.0;
4749 factors[StdRegions::eFactorCoeffD00] = 1.25;
4750 factors[StdRegions::eFactorCoeffD01] = 0.25;
4751 factors[StdRegions::eFactorCoeffD11] = 1.25;
4752 factors[StdRegions::eFactorCoeffD02] = 0.25;
4753 factors[StdRegions::eFactorCoeffD12] = 0.25;
4754 factors[StdRegions::eFactorCoeffD22] = 1.25;
4755
4756 c.Initialise(Collections::eHelmholtz, factors);
4757
4758 const int nm = Exp->GetNcoeffs();
4759 Array<OneD, NekDouble> coeffsIn(nelmts * nm);
4760 Array<OneD, NekDouble> coeffsRef(nelmts * nm);
4761 Array<OneD, NekDouble> coeffs(nelmts * nm), tmp;
4762
4763 for (int i = 0; i < nm; ++i)
4764 {
4765 coeffsIn[i] = 1.0;
4766 }
4767
4768 for (int i = 1; i < nelmts; ++i)
4769 {
4770 Vmath::Vcopy(nm, coeffsIn, 1, tmp = coeffsIn + i * nm, 1);
4771 }
4772
4773 StdRegions::StdMatrixKey mkey(StdRegions::eHelmholtz, Exp->DetShapeType(),
4774 *Exp, factors);
4775
4776 for (int i = 0; i < nelmts; ++i)
4777 {
4778 // Standard routines
4779 Exp->GeneralMatrixOp(coeffsIn + i * nm, tmp = coeffsRef + i * nm, mkey);
4780 }
4781
4782 c.ApplyOperator(Collections::eHelmholtz, coeffsIn, coeffs);
4783
4784 double epsilon = 1.0e-8;
4785 for (int i = 0; i < coeffsRef.size(); ++i)
4786 {
4787 coeffsRef[i] = (std::abs(coeffsRef[i]) < 1e-14) ? 0.0 : coeffsRef[i];
4788 coeffs[i] = (std::abs(coeffs[i]) < 1e-14) ? 0.0 : coeffs[i];
4789 BOOST_CHECK_CLOSE(coeffsRef[i], coeffs[i], epsilon);
4790 }
4791}
4792
4793BOOST_AUTO_TEST_CASE(TestHexPhysInterp1D_NoCollection_UniformP)
4794{
4795 SpatialDomains::PointGeomUniquePtr v0(
4796 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4797 SpatialDomains::PointGeomUniquePtr v1(
4798 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4799 SpatialDomains::PointGeomUniquePtr v2(
4800 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4801 SpatialDomains::PointGeomUniquePtr v3(
4802 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4803 SpatialDomains::PointGeomUniquePtr v4(
4804 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4805 SpatialDomains::PointGeomUniquePtr v5(
4806 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4807 SpatialDomains::PointGeomUniquePtr v6(
4808 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4809 SpatialDomains::PointGeomUniquePtr v7(
4810 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4811
4812 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4813 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4814 SpatialDomains::HexGeomUniquePtr hexGeom =
4815 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4816
4817 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4818 Nektar::LibUtilities::eGaussLobattoLegendre;
4819 Nektar::LibUtilities::BasisType basisTypeDir1 =
4820 Nektar::LibUtilities::eModified_A;
4821 unsigned int numQuadPoints = 5;
4822 unsigned int numModes = 4;
4823
4824 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4825 quadPointsTypeDir1);
4826 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4827 quadPointsKeyDir1);
4828
4829 Nektar::LocalRegions::HexExpSharedPtr Exp =
4830 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4831 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4832
4833 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4834 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4835 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4836
4837 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4838 CollExp.push_back(Exp);
4839
4840 LibUtilities::SessionReaderSharedPtr dummySession;
4841 Collections::CollectionOptimisation colOpt(dummySession, 3,
4842 Collections::eNoCollection);
4843 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(Exp);
4844 Collections::Collection c(CollExp, impTypes);
4845
4846 StdRegions::ConstFactorMap factors;
4847 factors[StdRegions::eFactorConst] = 1.5;
4848 c.Initialise(Collections::ePhysInterp1DScaled, factors);
4849
4850 const int nq = Exp->GetTotPoints();
4851
4852 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
4853 Array<OneD, NekDouble> phys(nq), tmp;
4854
4855 Exp->GetCoords(xc, yc, zc);
4856
4857 for (int i = 0; i < nq; ++i)
4858 {
4859 phys[i] = pow(xc[i], 3) + pow(yc[i], 3) + pow(zc[i], 3);
4860 }
4861
4862 const int nq1 = c.GetOutputSize(Collections::ePhysInterp1DScaled);
4863 Array<OneD, NekDouble> xc1(nq1);
4864 Array<OneD, NekDouble> yc1(nq1);
4865 Array<OneD, NekDouble> zc1(nq1);
4866 Array<OneD, NekDouble> phys1(nq1);
4867
4868 c.ApplyOperator(Collections::ePhysInterp1DScaled, xc, xc1);
4869 c.ApplyOperator(Collections::ePhysInterp1DScaled, yc, yc1);
4870 c.ApplyOperator(Collections::ePhysInterp1DScaled, zc, zc1);
4871 c.ApplyOperator(Collections::ePhysInterp1DScaled, phys, phys1);
4872
4873 double epsilon = 1.0e-8;
4874 // since solution is a polynomial should be able to compare soln directly
4875 for (int i = 0; i < nq1; ++i)
4876 {
4877 NekDouble exact = pow(xc1[i], 3) + pow(yc1[i], 3) + pow(zc1[i], 3);
4878 phys1[i] = (fabs(phys1[i]) < 1e-14) ? 0.0 : phys1[i];
4879 exact = (fabs(exact) < 1e-14) ? 0.0 : exact;
4880 BOOST_CHECK_CLOSE(phys1[i], exact, epsilon);
4881 }
4882}
4883
4884BOOST_AUTO_TEST_CASE(TestHexPhysInterp1D_MatrixFree_UniformP)
4885{
4886 SpatialDomains::PointGeomUniquePtr v0(
4887 new SpatialDomains::PointGeom(3u, 0u, -1.0, -1.0, -1.0));
4888 SpatialDomains::PointGeomUniquePtr v1(
4889 new SpatialDomains::PointGeom(3u, 1u, 1.0, -1.0, -1.0));
4890 SpatialDomains::PointGeomUniquePtr v2(
4891 new SpatialDomains::PointGeom(3u, 2u, 1.0, 1.0, -1.0));
4892 SpatialDomains::PointGeomUniquePtr v3(
4893 new SpatialDomains::PointGeom(3u, 3u, -1.0, 1.0, -1.0));
4894 SpatialDomains::PointGeomUniquePtr v4(
4895 new SpatialDomains::PointGeom(3u, 4u, -1.0, -1.0, 1.0));
4896 SpatialDomains::PointGeomUniquePtr v5(
4897 new SpatialDomains::PointGeom(3u, 5u, 1.0, -1.0, 1.0));
4898 SpatialDomains::PointGeomUniquePtr v6(
4899 new SpatialDomains::PointGeom(3u, 6u, 2.0, 3.0, 4.0));
4900 SpatialDomains::PointGeomUniquePtr v7(
4901 new SpatialDomains::PointGeom(3u, 7u, -1.0, 1.0, 1.0));
4902
4903 std::vector<SpatialDomains::SegGeomUniquePtr> segVec;
4904 std::vector<SpatialDomains::QuadGeomUniquePtr> faceVec;
4905 SpatialDomains::HexGeomUniquePtr hexGeom =
4906 CreateHex(v0.get(), v1.get(), v2.get(), v3.get(), v4.get(), v5.get(), v6.get(), v7.get(), segVec, faceVec);
4907
4908 Nektar::LibUtilities::PointsType quadPointsTypeDir1 =
4909 Nektar::LibUtilities::eGaussLobattoLegendre;
4910 Nektar::LibUtilities::BasisType basisTypeDir1 =
4911 Nektar::LibUtilities::eModified_A;
4912 unsigned int numQuadPoints = 5;
4913 unsigned int numModes = 4;
4914
4915 const Nektar::LibUtilities::PointsKey quadPointsKeyDir1(numQuadPoints,
4916 quadPointsTypeDir1);
4917 const Nektar::LibUtilities::BasisKey basisKeyDir1(basisTypeDir1, numModes,
4918 quadPointsKeyDir1);
4919
4920 Nektar::LocalRegions::HexExpSharedPtr Exp =
4921 MemoryManager<Nektar::LocalRegions::HexExp>::AllocateSharedPtr(
4922 basisKeyDir1, basisKeyDir1, basisKeyDir1, hexGeom.get());
4923
4924 Nektar::StdRegions::StdHexExpSharedPtr stdExp =
4925 MemoryManager<Nektar::StdRegions::StdHexExp>::AllocateSharedPtr(
4926 basisKeyDir1, basisKeyDir1, basisKeyDir1);
4927
4928 std::vector<StdRegions::StdExpansionSharedPtr> CollExp;
4929 CollExp.push_back(Exp);
4930
4931 LibUtilities::SessionReaderSharedPtr dummySession;
4932 Collections::CollectionOptimisation colOpt(dummySession, 3,
4933 Collections::eMatrixFree);
4934 Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(Exp);
4935 Collections::Collection c(CollExp, impTypes);
4936
4937 StdRegions::ConstFactorMap factors;
4938 factors[StdRegions::eFactorConst] = 1.5;
4939 c.Initialise(Collections::ePhysInterp1DScaled, factors);
4940
4941 const int nq = Exp->GetTotPoints();
4942
4943 Array<OneD, NekDouble> xc(nq), yc(nq), zc(nq);
4944 Array<OneD, NekDouble> phys(nq), tmp;
4945
4946 Exp->GetCoords(xc, yc, zc);
4947
4948 for (int i = 0; i < nq; ++i)
4949 {
4950 phys[i] = pow(xc[i], 3) + pow(yc[i], 3) + pow(zc[i], 3);
4951 }
4952
4953 const int nq1 = c.GetOutputSize(Collections::ePhysInterp1DScaled);
4954 Array<OneD, NekDouble> xc1(nq1);
4955 Array<OneD, NekDouble> yc1(nq1);
4956 Array<OneD, NekDouble> zc1(nq1);
4957 Array<OneD, NekDouble> phys1(nq1);
4958
4959 c.ApplyOperator(Collections::ePhysInterp1DScaled, xc, xc1);
4960 c.ApplyOperator(Collections::ePhysInterp1DScaled, yc, yc1);
4961 c.ApplyOperator(Collections::ePhysInterp1DScaled, zc, zc1);
4962 c.ApplyOperator(Collections::ePhysInterp1DScaled, phys, phys1);
4963
4964 double epsilon = 1.0e-8;
4965 // since solution is a polynomial should be able to compare soln directly
4966 for (int i = 0; i < nq1; ++i)
4967 {
4968 NekDouble exact = pow(xc1[i], 3) + pow(yc1[i], 3) + pow(zc1[i], 3);
4969 phys1[i] = (fabs(phys1[i]) < 1e-14) ? 0.0 : phys1[i];
4970 exact = (fabs(exact) < 1e-14) ? 0.0 : exact;
4971 BOOST_CHECK_CLOSE(phys1[i], exact, epsilon);
4972 }
4973}
4974#endif
4975} // namespace Nektar::HexCollectionTests
Nektar::Collections::Collection::Initialise
COLLECTIONS_EXPORT void Initialise(const OperatorType opType, StdRegions::FactorMap factors=StdRegions::NullFactorMap)
Definition Collection.cpp:66
Nektar::Collections::Collection::ApplyOperator
void ApplyOperator(const OperatorType &op, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &output)
Definition Collection.h:138
Nektar::Collections::CollectionOptimisation::GetOperatorImpMap
COLLECTIONS_EXPORT OperatorImpMap GetOperatorImpMap(StdRegions::StdExpansionSharedPtr pExp)
Get Operator Implementation Map from XMl or using default;.
Definition CollectionOptimisation.cpp:356
Nektar::LibUtilities::BasisKey
Describes the specification for a Basis.
Definition Basis.h:45
Nektar::LibUtilities::PointsKey
Defines a specification for a set of points.
Definition Points.h:50
Nektar::MemoryManager::AllocateSharedPtr
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Definition NekMemoryManager.hpp:166
Nektar::Collections::OperatorImpMap
std::map< OperatorType, ImplementationType > OperatorImpMap
Definition Operator.h:131
Nektar::HexCollectionTests
The above copyright notice and this permission notice shall be included.
Definition TestHexCollection.cpp:42
Nektar::HexCollectionTests::BOOST_AUTO_TEST_CASE
BOOST_AUTO_TEST_CASE(TestHexBwdTrans_StdMat_UniformP)
Definition TestHexCollection.cpp:103
Nektar::HexCollectionTests::CreateSegGeom
SpatialDomains::SegGeomUniquePtr CreateSegGeom(unsigned int id, SpatialDomains::PointGeom *v0, SpatialDomains::PointGeom *v1)
Definition TestHexCollection.cpp:44
Nektar::HexCollectionTests::CreateHex
SpatialDomains::HexGeomUniquePtr CreateHex(std::array< SpatialDomains::PointGeom *, 8 > v, std::array< SpatialDomains::SegGeomUniquePtr, 12 > &segVec, std::array< SpatialDomains::QuadGeomUniquePtr, 6 > &faceVec)
Definition TestHexCollection.cpp:54
Nektar::LibUtilities::SessionReaderSharedPtr
std::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition SessionReader.h:119
Nektar::LibUtilities::eGaussLobattoLegendre
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition PointsType.h:51
Nektar::LibUtilities::eGLL_Lagrange
@ eGLL_Lagrange
Lagrange for SEM basis .
Definition BasisType.h:56
Nektar::LibUtilities::eModified_A
@ eModified_A
Principle Modified Functions .
Definition BasisType.h:48
Nektar::LocalRegions::HexExpSharedPtr
std::shared_ptr< HexExp > HexExpSharedPtr
Definition HexExp.h:258
Nektar::SpatialDomains::HexGeomUniquePtr
unique_ptr_objpool< HexGeom > HexGeomUniquePtr
Definition MeshGraph.h:104
Nektar::SpatialDomains::QuadGeomUniquePtr
unique_ptr_objpool< QuadGeom > QuadGeomUniquePtr
Definition MeshGraph.h:100
Nektar::SpatialDomains::SegGeomUniquePtr
unique_ptr_objpool< SegGeom > SegGeomUniquePtr
Definition MeshGraph.h:98
Nektar::SpatialDomains::PointGeomUniquePtr
unique_ptr_objpool< PointGeom > PointGeomUniquePtr
Definition MeshGraph.h:93
Nektar::StdRegions::ConstFactorMap
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition StdRegions.hpp:430
Nektar::StdRegions::StdHexExpSharedPtr
std::shared_ptr< StdHexExp > StdHexExpSharedPtr
Definition StdHexExp.h:225
Nektar::VarcoeffHashingTest::factors
StdRegions::ConstFactorMap factors
Definition TestVarcoeffHashing.cpp:52
Vmath::Vadd
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition Vmath.hpp:180
Vmath::Vcopy
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition Vmath.hpp:825
Generated on Sun Nov 9 2025 20:17:00 for Nektar++ by doxygen 1.9.8