Nektar++
AssemblyMap.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: AssemblyMap.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: Assembly (e.g. local to global) base mapping routines
32//
33///////////////////////////////////////////////////////////////////////////////
34
35#include <boost/algorithm/string.hpp>
36
38
39using namespace std;
40
42{
43/**
44 * @class AssemblyMap
45 * This class acts as a base class for constructing mappings between
46 * local, global and boundary degrees of freedom. It holds the storage
47 * for the maps and provides the accessors needed to retrieve them.
48 *
49 * There are two derived classes: AssemblyMapCG and
50 * AssemblyMapDG. These perform the actual construction of the
51 * maps within their specific contexts.
52 *
53 */
54
55/// Rounds a double precision number to an integer.
57{
58 return int(x > 0.0 ? x + 0.5 : x - 0.5);
59}
60
61/// Rounds an array of double precision numbers to integers.
63 Array<OneD, int> outarray)
64{
65 int size = inarray.size();
66 ASSERTL1(outarray.size() >= size, "Array sizes not compatible");
67
68 NekDouble x;
69 for (int i = 0; i < size; i++)
70 {
71 x = inarray[i];
72 outarray[i] = int(x > 0.0 ? x + 0.5 : x - 0.5);
73 }
74}
75
76/**
77 * Initialises an empty mapping.
78 */
80 : m_session(), m_comm(), m_hash(0), m_numLocalBndCoeffs(0),
81 m_numGlobalBndCoeffs(0), m_numLocalDirBndCoeffs(0),
82 m_numGlobalDirBndCoeffs(0), m_solnType(eNoSolnType),
83 m_bndSystemBandWidth(0), m_successiveRHS(0),
84 m_linSysIterSolver("ConjugateGradient"), m_gsh(nullptr), m_bndGsh(nullptr)
85{
86}
87
90 const std::string variable)
91 : m_session(pSession), m_comm(comm), m_variable(variable), m_hash(0),
92 m_numLocalBndCoeffs(0), m_numGlobalBndCoeffs(0),
93 m_numLocalDirBndCoeffs(0), m_numGlobalDirBndCoeffs(0),
94 m_bndSystemBandWidth(0), m_successiveRHS(0),
95 m_linSysIterSolver("ConjugateGradient"), m_gsh(nullptr), m_bndGsh(nullptr)
96{
97 // Default value from Solver Info
99 pSession->GetSolverInfoAsEnum<GlobalSysSolnType>("GlobalSysSoln");
100
101 // Override values with data from GlobalSysSolnInfo section
102 if (pSession->DefinesGlobalSysSolnInfo(variable, "GlobalSysSoln"))
103 {
104 std::string sysSoln =
105 pSession->GetGlobalSysSolnInfo(variable, "GlobalSysSoln");
106 m_solnType = pSession->GetValueAsEnum<GlobalSysSolnType>(
107 "GlobalSysSoln", sysSoln);
108 }
109
110 // Set up preconditioner with SysSoln as override then solverinfo otherwise
111 // set a default as diagonal
112 if (pSession->DefinesGlobalSysSolnInfo(variable, "Preconditioner"))
113 {
115 pSession->GetGlobalSysSolnInfo(variable, "Preconditioner");
116 }
117 else if (pSession->DefinesSolverInfo("Preconditioner"))
118 {
119 m_preconType = pSession->GetSolverInfo("Preconditioner");
120 }
121 else
122 { // Possibly preconditioner is default registered in
123 // GlobLinSysIterative.cpp
124 m_preconType = "Diagonal";
125 }
126
127 if (pSession->DefinesGlobalSysSolnInfo(variable, "AbsoluteTolerance"))
128 {
129 std::string abstol =
130 pSession->GetGlobalSysSolnInfo(variable, "AbsoluteTolerance");
132 boost::iequals(boost::to_upper_copy(abstol), "TRUE");
133 }
134 else
135 {
136 m_isAbsoluteTolerance = false;
137 }
138
139 if (pSession->DefinesGlobalSysSolnInfo(variable, "SuccessiveRHS"))
140 {
141 m_successiveRHS = std::stoi(
142 pSession->GetGlobalSysSolnInfo(variable, "SuccessiveRHS").c_str());
143 }
144 else
145 {
146 pSession->LoadParameter("SuccessiveRHS", m_successiveRHS, 0);
147 }
148
149 if (pSession->DefinesGlobalSysSolnInfo(variable, "LinSysIterSolver"))
150 {
152 pSession->GetGlobalSysSolnInfo(variable, "LinSysIterSolver");
153 }
154 else if (pSession->DefinesSolverInfo("LinSysIterSolver"))
155 {
156 m_linSysIterSolver = pSession->GetSolverInfo("LinSysIterSolver");
157 }
158 else
159 {
160 m_linSysIterSolver = "ConjugateGradient";
161 }
162}
163
164/**
165 * Create a new level of mapping using the information in
166 * multiLevelGraph and performing the following steps:
167 */
169 AssemblyMap *oldLevelMap,
170 const BottomUpSubStructuredGraphSharedPtr &multiLevelGraph)
171 : m_session(oldLevelMap->m_session), m_comm(oldLevelMap->GetComm()),
172 m_hash(0), m_solnType(oldLevelMap->m_solnType),
173 m_preconType(oldLevelMap->m_preconType),
174 m_successiveRHS(oldLevelMap->m_successiveRHS),
175 m_linSysIterSolver(oldLevelMap->m_linSysIterSolver),
176 m_gsh(oldLevelMap->m_gsh), m_bndGsh(oldLevelMap->m_bndGsh),
177 m_lowestStaticCondLevel(oldLevelMap->m_lowestStaticCondLevel)
178{
179 int i;
180 int j;
181 int cnt;
182
183 //--------------------------------------------------------------
184 // -- Extract information from the input argument
185 int numGlobalBndCoeffsOld = oldLevelMap->GetNumGlobalBndCoeffs();
186 int numGlobalDirBndCoeffsOld = oldLevelMap->GetNumGlobalDirBndCoeffs();
187 int numLocalBndCoeffsOld = oldLevelMap->GetNumLocalBndCoeffs();
188 int numLocalDirBndCoeffsOld = oldLevelMap->GetNumLocalDirBndCoeffs();
189 bool signChangeOld = oldLevelMap->GetSignChange();
190
191 int staticCondLevelOld = oldLevelMap->GetStaticCondLevel();
192 int numPatchesOld = oldLevelMap->GetNumPatches();
193 GlobalSysSolnType solnTypeOld = oldLevelMap->GetGlobalSysSolnType();
194 const Array<OneD, const unsigned int> &numLocalBndCoeffsPerPatchOld =
195 oldLevelMap->GetNumLocalBndCoeffsPerPatch();
196 //--------------------------------------------------------------
197
198 //--------------------------------------------------------------
199 int newLevel = staticCondLevelOld + 1;
200 /** - STEP 1: setup a mask array to determine to which patch
201 * of the new level every patch of the current
202 * level belongs. To do so we make four arrays,
203 * #gloPatchMask, #globHomPatchMask,
204 * #locPatchMask_NekDouble and #locPatchMask.
205 * These arrays are then used to check which local
206 * dofs of the old level belong to which patch of
207 * the new level
208 */
209 Array<OneD, NekDouble> globPatchMask(numGlobalBndCoeffsOld, -1.0);
210 Array<OneD, NekDouble> globHomPatchMask(globPatchMask +
211 numGlobalDirBndCoeffsOld);
212 Array<OneD, NekDouble> locPatchMask_NekDouble(numLocalBndCoeffsOld, -3.0);
213 Array<OneD, int> locPatchMask(numLocalBndCoeffsOld);
214
215 // Fill the array globPatchMask as follows:
216 // - The first part (i.e. the glob bnd dofs) is filled with the
217 // value -1
218 // - The second part (i.e. the glob interior dofs) is numbered
219 // according to the patch it belongs to (i.e. dofs in first block
220 // all are numbered 0, the second block numbered are 1, etc...)
221 multiLevelGraph->MaskPatches(newLevel, globHomPatchMask);
222
223 // Map from Global Dofs to Local Dofs
224 // As a result, we know for each local dof whether
225 // it is mapped to the boundary of the next level, or to which
226 // patch. Based upon this, we can than later associate every patch
227 // of the current level with a patch in the next level.
228 oldLevelMap->GlobalToLocalBndWithoutSign(globPatchMask,
229 locPatchMask_NekDouble);
230
231 // Convert the result to an array of integers rather than doubles
232 RoundNekDoubleToInt(locPatchMask_NekDouble, locPatchMask);
233
234 /** - STEP 2: We calculate how many local bnd dofs of the
235 * old level belong to the boundaries of each patch at
236 * the new level. We need this information to set up the
237 * mapping between different levels.
238 */
239
240 // Retrieve the number of patches at the next level
241 int numPatchesWithIntNew =
242 multiLevelGraph->GetNpatchesWithInterior(newLevel);
243 int numPatchesNew = numPatchesWithIntNew;
244
245 // Allocate memory to store the number of local dofs associated to
246 // each of elemental boundaries of these patches
247 std::map<int, int> numLocalBndCoeffsPerPatchNew;
248 for (int i = 0; i < numPatchesNew; i++)
249 {
250 numLocalBndCoeffsPerPatchNew[i] = 0;
251 }
252
253 int minval;
254 int maxval;
255 int curPatch;
256 for (i = cnt = 0; i < numPatchesOld; i++)
257 {
258 // For every patch at the current level, the mask array
259 // locPatchMask should be filled with
260 // - the same (positive) number for each entry (which will
261 // correspond to the patch at the next level it belongs to)
262 // - the same (positive) number for each entry, except some
263 // entries that are -1 (the enties correspond to -1, will be
264 // mapped to the local boundary of the next level patch given
265 // by the positive number)
266 // - -1 for all entries. In this case, we will make an
267 // additional patch only consisting of boundaries at the next
268 // level
269 minval =
270 *min_element(&locPatchMask[cnt],
271 &locPatchMask[cnt] + numLocalBndCoeffsPerPatchOld[i]);
272 maxval =
273 *max_element(&locPatchMask[cnt],
274 &locPatchMask[cnt] + numLocalBndCoeffsPerPatchOld[i]);
275 ASSERTL0((minval == maxval) || (minval == -1),
276 "These values should never be the same");
277
278 if (maxval == -1)
279 {
280 curPatch = numPatchesNew;
281 numLocalBndCoeffsPerPatchNew[curPatch] = 0;
282 numPatchesNew++;
283 }
284 else
285 {
286 curPatch = maxval;
287 }
288
289 for (j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++)
290 {
291 ASSERTL0((locPatchMask[cnt] == maxval) ||
292 (locPatchMask[cnt] == minval),
293 "These values should never be the same");
294 if (locPatchMask[cnt] == -1)
295 {
296 ++numLocalBndCoeffsPerPatchNew[curPatch];
297 }
298 cnt++;
299 }
300 }
301
302 // Count how many local dofs of the old level are mapped
303 // to the local boundary dofs of the new level
306 m_numPatches = numLocalBndCoeffsPerPatchNew.size();
309 multiLevelGraph->GetNintDofsPerPatch(newLevel, m_numLocalIntCoeffsPerPatch);
310
311 for (int i = 0; i < m_numPatches; i++)
312 {
314 (unsigned int)numLocalBndCoeffsPerPatchNew[i];
318 }
319
320 // Also initialise some more data members
321 m_solnType = solnTypeOld;
326 "This method should only be called for in "
327 "case of multi-level static condensation.");
328 m_staticCondLevel = newLevel;
329 m_signChange = signChangeOld;
330 m_numLocalDirBndCoeffs = numLocalDirBndCoeffsOld;
331 m_numGlobalDirBndCoeffs = numGlobalDirBndCoeffsOld;
333 multiLevelGraph->GetInteriorOffset(newLevel) + m_numGlobalDirBndCoeffs;
335 multiLevelGraph->GetNumGlobalDofs(newLevel) + m_numGlobalDirBndCoeffs;
337
341
342 // local to bnd map is just a copy
343 for (int i = 0; i < m_numLocalBndCoeffs; ++i)
344 {
346 }
347
348 // local to int map is just a copy plus offset
349 for (int i = m_numLocalBndCoeffs; i < m_numLocalCoeffs; ++i)
350 {
352 }
353
354 if (m_signChange)
355 {
357 }
358
360 MemoryManager<PatchMap>::AllocateSharedPtr(numLocalBndCoeffsOld);
361
366
367 // Set up an offset array that denotes the offset of the local
368 // boundary degrees of freedom of the next level
369 Array<OneD, int> numLocalBndCoeffsPerPatchOffset(m_numPatches + 1, 0);
370 for (int i = 1; i < m_numPatches + 1; i++)
371 {
372 numLocalBndCoeffsPerPatchOffset[i] +=
373 numLocalBndCoeffsPerPatchOffset[i - 1] +
374 numLocalBndCoeffsPerPatchNew[i - 1];
375 }
376
377 int additionalPatchCnt = numPatchesWithIntNew;
378 int newid;
379 int blockid;
380 bool isBndDof;
382 Array<OneD, int> bndDofPerPatchCnt(m_numPatches, 0);
383
384 for (i = cnt = 0; i < numPatchesOld; i++)
385 {
386 minval =
387 *min_element(&locPatchMask[cnt],
388 &locPatchMask[cnt] + numLocalBndCoeffsPerPatchOld[i]);
389 maxval =
390 *max_element(&locPatchMask[cnt],
391 &locPatchMask[cnt] + numLocalBndCoeffsPerPatchOld[i]);
392 ASSERTL0((minval == maxval) || (minval == -1),
393 "These values should never be the same");
394
395 if (maxval == -1)
396 {
397 curPatch = additionalPatchCnt;
398 additionalPatchCnt++;
399 }
400 else
401 {
402 curPatch = maxval;
403 }
404
405 for (j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++)
406 {
407 ASSERTL0((locPatchMask[cnt] == maxval) ||
408 (locPatchMask[cnt] == minval),
409 "These values should never be the same");
410
411 sign = oldLevelMap->GetLocalToGlobalBndSign(cnt);
412
413 if (locPatchMask[cnt] == -1)
414 {
415 newid = numLocalBndCoeffsPerPatchOffset[curPatch];
416
417 m_localToGlobalBndMap[newid] =
418 oldLevelMap->GetLocalToGlobalBndMap(cnt);
419
420 if (m_signChange)
421 {
423 }
424
425 blockid = bndDofPerPatchCnt[curPatch];
426 isBndDof = true;
427
428 numLocalBndCoeffsPerPatchOffset[curPatch]++;
429 bndDofPerPatchCnt[curPatch]++;
430 }
431 else
432 {
433 newid = oldLevelMap->GetLocalToGlobalBndMap(cnt) -
435
436 blockid =
437 oldLevelMap->GetLocalToGlobalBndMap(cnt) -
439 multiLevelGraph->GetInteriorOffset(newLevel, curPatch);
440 isBndDof = false;
441 }
442
443 sign = isBndDof ? 1.0 : sign;
444
445 m_patchMapFromPrevLevel->SetPatchMap(cnt, curPatch, blockid,
446 isBndDof, sign);
447 cnt++;
448 }
449 }
450
451 // set up local to local mapping from previous to new level
454
456
457 // Postprocess the computed information - Update the old
458 // level with the mapping to new level
459 // oldLevelMap->SetLocalBndToNextLevelMap(oldLocalBndToNewLevelMap,oldLocalBndToNewLevelSign);
460
461 // - Construct the next level mapping object
462 if (m_staticCondLevel < (multiLevelGraph->GetNlevels() - 1))
463 {
466 multiLevelGraph);
467 }
468}
469
471{
472}
473
474/**
475 * The bandwidth calculated corresponds to what is referred to as
476 * half-bandwidth. If the elements of the matrix are designated as
477 * a_ij, it corresponds to the maximum value of |i-j| for non-zero
478 * a_ij. As a result, the value also corresponds to the number of
479 * sub- or super-diagonals.
480 *
481 * The bandwith can be calculated elementally as it corresponds to the
482 * maximal elemental bandwith (i.e. the maximal difference in global
483 * DOF index for every element).
484 *
485 * We here calculate the bandwith of the global boundary system (as
486 * used for static condensation).
487 */
489{
490 int i, j;
491 int cnt = 0;
492 int locSize;
493 int maxId;
494 int minId;
495 int bwidth = -1;
496 for (i = 0; i < m_numPatches; ++i)
497 {
498 locSize = m_numLocalBndCoeffsPerPatch[i];
499 maxId = -1;
500 minId = m_numLocalBndCoeffs + 1;
501 for (j = 0; j < locSize; j++)
502 {
504 {
505 if (m_localToGlobalBndMap[cnt + j] > maxId)
506 {
507 maxId = m_localToGlobalBndMap[cnt + j];
508 }
509
510 if (m_localToGlobalBndMap[cnt + j] < minId)
511 {
512 minId = m_localToGlobalBndMap[cnt + j];
513 }
514 }
515 }
516 bwidth = (bwidth > (maxId - minId)) ? bwidth : (maxId - minId);
517
518 cnt += locSize;
519 }
520
521 m_bndSystemBandWidth = bwidth;
522}
523
524int AssemblyMap::v_GetLocalToGlobalMap([[maybe_unused]] const int i) const
525{
526 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
527 return 0;
528}
529
530int AssemblyMap::v_GetGlobalToUniversalMap([[maybe_unused]] const int i) const
531{
532 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
533 return 0;
534}
535
537 [[maybe_unused]] const int i) const
538{
539 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
540 return 0;
541}
542
544{
545 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
546 static Array<OneD, const int> result;
547 return result;
548}
549
551{
552 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
553 static Array<OneD, const int> result;
554 return result;
555}
556
558{
559 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
560 static Array<OneD, const int> result;
561 return result;
562}
563
565 [[maybe_unused]] const int i) const
566{
567 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
568 return 0.0;
569}
570
572{
573 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
574 static Array<OneD, NekDouble> result;
575 return result;
576}
577
579 [[maybe_unused]] const Array<OneD, const NekDouble> &loc,
580 [[maybe_unused]] Array<OneD, NekDouble> &global,
581 [[maybe_unused]] bool useComm) const
582{
583 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
584}
585
587 [[maybe_unused]] const Array<OneD, const NekDouble> &global,
588 [[maybe_unused]] Array<OneD, NekDouble> &loc) const
589{
590 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
591}
592
594 [[maybe_unused]] const NekVector<NekDouble> &global,
595 [[maybe_unused]] NekVector<NekDouble> &loc) const
596{
597 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
598}
599
601 [[maybe_unused]] const Array<OneD, const NekDouble> &loc,
602 [[maybe_unused]] Array<OneD, NekDouble> &global) const
603{
604 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
605}
606
608 [[maybe_unused]] const NekVector<NekDouble> &loc,
609 [[maybe_unused]] NekVector<NekDouble> &global) const
610{
611 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
612}
613
615 [[maybe_unused]] Array<OneD, NekDouble> &pGlobal) const
616{
617 // Do nothing here since multi-level static condensation uses a
618 // AssemblyMap and thus will call this routine in serial.
619}
620
622 [[maybe_unused]] Array<OneD, NekDouble> &pGlobal,
623 [[maybe_unused]] int offset) const
624{
625 // Do nothing here since multi-level static condensation uses a
626 // AssemblyMap and thus will call this routine in serial.
627}
628
630{
631 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
632 return 0;
633}
634
636{
637 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
638 return 0;
639}
640
642{
643 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
644 return 0;
645}
646
648{
649 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
650 return 0;
651}
652
654{
655 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
656 return 0;
657}
658
660{
661 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
662 return 0;
663}
664
666{
667 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
668 return 0;
669}
670
672{
673 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
674 return 0;
675}
676
678{
679 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
680 static Array<OneD, const int> result;
681 return result;
682}
683
684std::shared_ptr<AssemblyMap> AssemblyMap::v_LinearSpaceMap(
685 [[maybe_unused]] const ExpList &locexp,
686 [[maybe_unused]] GlobalSysSolnType solnType)
687{
688 NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
689 static std::shared_ptr<AssemblyMap> result;
690 return result;
691}
692
694{
695 return m_comm;
696}
697
699{
700 return m_variable;
701}
702
704{
705 return m_hash;
706}
707
709{
710 return v_GetLocalToGlobalMap(i);
711}
712
714{
716}
717
719{
721}
722
724{
725 return v_GetLocalToGlobalMap();
726}
727
729{
731}
732
734{
736}
737
739{
740 return v_GetLocalToGlobalSign(i);
741}
742
744{
745 return v_GetLocalToGlobalSign();
746}
747
750 bool useComm) const
751{
752 v_LocalToGlobal(loc, global, useComm);
753}
754
756 NekVector<NekDouble> &global,
757 bool useComm) const
758{
759 v_LocalToGlobal(loc.GetPtr(), global.GetPtr(), useComm);
760}
761
764{
765 v_GlobalToLocal(global, loc);
766}
767
770{
771 v_GlobalToLocal(global, loc);
772}
773
775 Array<OneD, NekDouble> &global) const
776{
777 v_Assemble(loc, global);
778}
779
781 NekVector<NekDouble> &global) const
782{
783 v_Assemble(loc, global);
784}
785
787{
788 v_UniversalAssemble(pGlobal);
789}
790
792{
793 v_UniversalAssemble(pGlobal.GetPtr());
794}
795
797 int offset) const
798{
799 v_UniversalAssemble(pGlobal, offset);
800}
801
803 Array<OneD, NekDouble> &global) const
804{
806
807 Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
809
810 if (global.data() == loc.data())
811 {
812 local = Array<OneD, NekDouble>(map.size(), loc.data());
813 }
814 else
815 {
816 local = loc; // create reference
817 }
818
819 Vmath::Scatr(map.size(), sign.data(), local.data(), map.data(),
820 global.data());
821}
822
825{
827
828 Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
830
831 if (global.data() == loc.data())
832 {
833 glo = Array<OneD, NekDouble>(global.size(), global.data());
834 }
835 else
836 {
837 glo = global; // create reference
838 }
839
840 Vmath::Gathr(map.size(), sign.data(), glo.data(), map.data(), loc.data());
841}
842
844 Array<OneD, NekDouble> &global) const
845{
847 Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
849
850 if (global.data() == loc.data())
851 {
852 local = Array<OneD, NekDouble>(map.size(), loc.data());
853 }
854 else
855 {
856 local = loc; // create reference
857 }
858
859 // since we are calling mapping from level down from array
860 // the m_numLocaBndCoeffs represents the size of the
861 // boundary elements we need to assemble into
862 Vmath::Zero(m_numLocalCoeffs, global.data(), 1);
863
864 Vmath::Assmb(map.size(), sign.data(), local.data(), map.data(),
865 global.data());
866}
867
869{
871}
872
874{
876}
877
879{
881}
882
884{
886}
887
889{
890 return v_GetNumDirEdges();
891}
892
894{
895 return v_GetNumDirFaces();
896}
897
899{
900 return v_GetNumNonDirEdges();
901}
902
904{
905 return v_GetNumNonDirFaces();
906}
907
909{
910 return v_GetExtraDirEdges();
911}
912
913std::shared_ptr<AssemblyMap> AssemblyMap::LinearSpaceMap(
914 const ExpList &locexp, GlobalSysSolnType solnType)
915{
916 return v_LinearSpaceMap(locexp, solnType);
917}
918
920{
921 return m_localToGlobalBndMap[i];
922}
923
925{
927}
928
930{
931 return m_signChange;
932}
933
935{
937}
938
940{
942}
943
945{
947}
948
950{
951 if (m_signChange)
952 {
953 return m_localToGlobalBndSign[i];
954 }
955 else
956 {
957 return 1.0;
958 }
959}
960
962{
964}
965
967{
969}
970
972{
974}
975
977{
978 ASSERTL1(i < m_bndCondIDToGlobalTraceID.size(), "Index out of range.");
980}
981
983{
985}
986
988{
990}
991
993{
995}
996
998{
999 return m_numLocalBndCoeffs;
1000}
1001
1003{
1004 return m_numGlobalBndCoeffs;
1005}
1006
1008{
1009 return m_numLocalCoeffs;
1010}
1011
1013{
1014 return m_numGlobalCoeffs;
1015}
1016
1018{
1019 return m_systemSingular;
1020}
1021
1023 NekVector<NekDouble> &loc, int offset) const
1024{
1025 GlobalToLocalBnd(global.GetPtr(), loc.GetPtr(), offset);
1026}
1027
1030{
1031 GlobalToLocalBnd(global.GetPtr(), loc.GetPtr());
1032}
1033
1036 int offset) const
1037{
1039 "Local vector is not of correct dimension");
1040 ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1041 "Global vector is not of correct dimension");
1042
1043 // offset input data by length "offset" for Dirichlet boundary conditions.
1045 Vmath::Vcopy(m_numGlobalBndCoeffs - offset, global.data(), 1,
1046 tmp.data() + offset, 1);
1047
1048 if (m_signChange)
1049 {
1051 tmp.data(), m_localToGlobalBndMap.data(), loc.data());
1052 }
1053 else
1054 {
1056 m_localToGlobalBndMap.data(), loc.data());
1057 }
1058}
1059
1062{
1064 "Local vector is not of correct dimension");
1065 ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1066 "Global vector is not of correct dimension");
1067
1069 if (global.data() == loc.data())
1070 {
1071 glo = Array<OneD, NekDouble>(m_numLocalBndCoeffs, global.data());
1072 }
1073 else
1074 {
1075 glo = global; // create reference
1076 }
1077
1078 if (m_signChange)
1079 {
1081 glo.data(), m_localToGlobalBndMap.data(), loc.data());
1082 }
1083 else
1084 {
1086 m_localToGlobalBndMap.data(), loc.data());
1087 }
1088}
1089
1091 Array<OneD, NekDouble> &global, int offset,
1092 bool UseComm) const
1093{
1095 "Local vector is not of correct dimension");
1096 ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1097 "Global vector is not of correct dimension");
1098
1099 // offset input data by length "offset" for Dirichlet boundary conditions.
1101
1102 if (m_signChange)
1103 {
1105 loc.data(), m_localToGlobalBndMap.data(), tmp.data());
1106 }
1107 else
1108 {
1110 m_localToGlobalBndMap.data(), tmp.data());
1111 }
1112
1113 // Ensure each processor has unique value with a max gather.
1114 if (UseComm)
1115 {
1117 }
1118 Vmath::Vcopy(m_numGlobalBndCoeffs - offset, tmp.data() + offset, 1,
1119 global.data(), 1);
1120}
1121
1123 Array<OneD, NekDouble> &global,
1124 bool UseComm) const
1125{
1127 "Local vector is not of correct dimension");
1128 ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1129 "Global vector is not of correct dimension");
1130
1131 if (m_signChange)
1132 {
1134 loc.data(), m_localToGlobalBndMap.data(), global.data());
1135 }
1136 else
1137 {
1139 m_localToGlobalBndMap.data(), global.data());
1140 }
1141 if (UseComm)
1142 {
1143 Gs::Gather(global, Gs::gs_max, m_bndGsh);
1144 }
1145}
1146
1148 Array<OneD, NekDouble> &locbnd) const
1149{
1150 ASSERTL1(locbnd.size() >= m_numLocalBndCoeffs,
1151 "LocBnd vector is not of correct dimension");
1152 ASSERTL1(local.size() >= m_numLocalCoeffs,
1153 "Local vector is not of correct dimension");
1154
1156 locbnd.data());
1157}
1158
1160 Array<OneD, NekDouble> &locint) const
1161{
1163 "Locint vector is not of correct dimension");
1164 ASSERTL1(local.size() >= m_numLocalCoeffs,
1165 "Local vector is not of correct dimension");
1166
1168 m_localToLocalIntMap.data(), locint.data());
1169}
1170
1172 Array<OneD, NekDouble> &local) const
1173{
1174 ASSERTL1(locbnd.size() >= m_numLocalBndCoeffs,
1175 "LocBnd vector is not of correct dimension");
1176 ASSERTL1(local.size() >= m_numLocalCoeffs,
1177 "Local vector is not of correct dimension");
1178
1179 Vmath::Scatr(m_numLocalBndCoeffs, locbnd.data(),
1180 m_localToLocalBndMap.data(), local.data());
1181}
1182
1184 Array<OneD, NekDouble> &local) const
1185{
1187 "LocBnd vector is not of correct dimension");
1188 ASSERTL1(local.size() >= m_numLocalCoeffs,
1189 "Local vector is not of correct dimension");
1190
1192 m_localToLocalIntMap.data(), local.data());
1193}
1194
1196 NekVector<NekDouble> &global, int offset) const
1197{
1198 AssembleBnd(loc.GetPtr(), global.GetPtr(), offset);
1199}
1200
1202 NekVector<NekDouble> &global) const
1203{
1204 AssembleBnd(loc.GetPtr(), global.GetPtr());
1205}
1206
1208 Array<OneD, NekDouble> &global, int offset) const
1209{
1211 "Local array is not of correct dimension");
1212 ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1213 "Global array is not of correct dimension");
1215
1216 if (m_signChange)
1217 {
1219 loc.data(), m_localToGlobalBndMap.data(), tmp.data());
1220 }
1221 else
1222 {
1224 m_localToGlobalBndMap.data(), tmp.data());
1225 }
1227 Vmath::Vcopy(m_numGlobalBndCoeffs - offset, tmp.data() + offset, 1,
1228 global.data(), 1);
1229}
1230
1232 Array<OneD, NekDouble> &global) const
1233{
1235 "Local vector is not of correct dimension");
1236 ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1237 "Global vector is not of correct dimension");
1238
1239 Vmath::Zero(m_numGlobalBndCoeffs, global.data(), 1);
1240
1241 if (m_signChange)
1242 {
1244 loc.data(), m_localToGlobalBndMap.data(), global.data());
1245 }
1246 else
1247 {
1249 m_localToGlobalBndMap.data(), global.data());
1250 }
1251 UniversalAssembleBnd(global);
1252}
1253
1255{
1256 ASSERTL1(pGlobal.size() >= m_numGlobalBndCoeffs, "Wrong size.");
1257 Gs::Gather(pGlobal, Gs::gs_add, m_bndGsh);
1258}
1259
1261{
1262 UniversalAssembleBnd(pGlobal.GetPtr());
1263}
1264
1266 int offset) const
1267{
1268 Array<OneD, NekDouble> tmp(offset);
1269 if (offset > 0)
1270 {
1271 Vmath::Vcopy(offset, pGlobal, 1, tmp, 1);
1272 }
1273 UniversalAssembleBnd(pGlobal);
1274 if (offset > 0)
1275 {
1276 Vmath::Vcopy(offset, tmp, 1, pGlobal, 1);
1277 }
1278}
1279
1281{
1283}
1284
1286{
1287 return m_bndSystemBandWidth;
1288}
1289
1291{
1292 return m_staticCondLevel;
1293}
1294
1296{
1297 return m_numPatches;
1298}
1299
1302{
1304}
1305
1308{
1310}
1311
1313{
1315}
1316
1318{
1320}
1321
1323{
1324 return !(m_nextLevelLocalToGlobalMap.get());
1325}
1326
1328{
1329 return m_solnType;
1330}
1331
1333{
1334 return m_preconType;
1335}
1336
1338{
1339 return m_isAbsoluteTolerance;
1340}
1341
1343{
1344 return m_successiveRHS;
1345}
1346
1348{
1349 return m_linSysIterSolver;
1350}
1351
1354{
1356 "Local vector is not of correct dimension");
1357 ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1358 "Global vector is not of correct dimension");
1359
1360 Vmath::Gathr(m_numLocalBndCoeffs, global.data(),
1361 m_localToGlobalBndMap.data(), loc.data());
1362}
1363
1364void AssemblyMap::PrintStats(std::ostream &out, std::string variable,
1365 bool printHeader) const
1366{
1367 LibUtilities::CommSharedPtr vRowComm = m_session->GetComm()->GetRowComm();
1368 bool isRoot = m_session->GetComm()->IsParallelInTime()
1369 ? m_session->GetComm()->GetRank() == 0
1370 : vRowComm->GetRank() == 0;
1371 int n = vRowComm->GetSize();
1372 int i;
1373
1374 // Determine number of global degrees of freedom.
1375 int globBndCnt = 0, globDirCnt = 0;
1376
1377 for (i = 0; i < m_numGlobalBndCoeffs; ++i)
1378 {
1380 {
1381 globBndCnt++;
1382
1384 {
1385 globDirCnt++;
1386 }
1387 }
1388 }
1389
1390 int globCnt = m_numGlobalCoeffs - m_numGlobalBndCoeffs + globBndCnt;
1391
1392 // Calculate maximum valency
1395
1396 Vmath::Assmb(m_numLocalBndCoeffs, tmpLoc.data(),
1397 m_localToGlobalBndMap.data(), tmpGlob.data());
1398 UniversalAssembleBnd(tmpGlob);
1399
1400 int totGlobDof = globCnt;
1401 int totGlobBndDof = globBndCnt;
1402 int totGlobDirDof = globDirCnt;
1403 int totLocalDof = m_numLocalCoeffs;
1404 int totLocalBndDof = m_numLocalBndCoeffs;
1405 int totLocalDirDof = m_numLocalDirBndCoeffs;
1406
1407 int meanValence = 0;
1408 int maxValence = 0;
1409 int minValence = 10000000;
1410 for (int i = 0; i < m_numGlobalBndCoeffs; ++i)
1411 {
1413 {
1414 continue;
1415 }
1416
1417 if (tmpGlob[i] > maxValence)
1418 {
1419 maxValence = tmpGlob[i];
1420 }
1421 if (tmpGlob[i] < minValence)
1422 {
1423 minValence = tmpGlob[i];
1424 }
1425 meanValence += tmpGlob[i];
1426 }
1427
1428 vRowComm->AllReduce(maxValence, LibUtilities::ReduceMax);
1429 vRowComm->AllReduce(minValence, LibUtilities::ReduceMin);
1430 vRowComm->AllReduce(meanValence, LibUtilities::ReduceSum);
1431 vRowComm->AllReduce(totGlobDof, LibUtilities::ReduceSum);
1432 vRowComm->AllReduce(totGlobBndDof, LibUtilities::ReduceSum);
1433 vRowComm->AllReduce(totGlobDirDof, LibUtilities::ReduceSum);
1434 vRowComm->AllReduce(totLocalDof, LibUtilities::ReduceSum);
1435 vRowComm->AllReduce(totLocalBndDof, LibUtilities::ReduceSum);
1436 vRowComm->AllReduce(totLocalDirDof, LibUtilities::ReduceSum);
1437
1438 meanValence /= totGlobBndDof;
1439
1440 if (isRoot)
1441 {
1442 if (printHeader)
1443 {
1444 out << "Assembly map statistics for field " << variable << ":"
1445 << endl;
1446 }
1447
1448 out << " - Number of local/global dof : " << totLocalDof
1449 << " " << totGlobDof << endl;
1450 out << " - Number of local/global boundary dof : " << totLocalBndDof
1451 << " " << totGlobBndDof << endl;
1452 out << " - Number of local/global Dirichlet dof : " << totLocalDirDof
1453 << " " << totGlobDirDof << endl;
1454 out << " - dof valency (min/max/mean) : " << minValence
1455 << " " << maxValence << " " << meanValence << endl;
1456
1457 if (n > 1)
1458 {
1459 NekDouble mean = m_numLocalCoeffs, mean2 = mean * mean;
1460 NekDouble minval = mean, maxval = mean;
1462
1463 for (i = 1; i < n; ++i)
1464 {
1465 vRowComm->Recv(i, tmp);
1466 mean += tmp[0];
1467 mean2 += tmp[0] * tmp[0];
1468
1469 if (tmp[0] > maxval)
1470 {
1471 maxval = tmp[0];
1472 }
1473 if (tmp[0] < minval)
1474 {
1475 minval = tmp[0];
1476 }
1477 }
1478
1479 if (maxval > 0.1)
1480 {
1481 out << " - Local dof dist. (min/max/mean/dev) : " << minval
1482 << " " << maxval << " " << (mean / n) << " "
1483 << sqrt(mean2 / n - mean * mean / n / n) << endl;
1484 }
1485
1486 vRowComm->Block();
1487
1488 mean = minval = maxval = m_numLocalBndCoeffs;
1489 mean2 = mean * mean;
1490
1491 for (i = 1; i < n; ++i)
1492 {
1493 vRowComm->Recv(i, tmp);
1494 mean += tmp[0];
1495 mean2 += tmp[0] * tmp[0];
1496
1497 if (tmp[0] > maxval)
1498 {
1499 maxval = tmp[0];
1500 }
1501 if (tmp[0] < minval)
1502 {
1503 minval = tmp[0];
1504 }
1505 }
1506
1507 out << " - Local bnd dof dist. (min/max/mean/dev) : " << minval
1508 << " " << maxval << " " << (mean / n) << " "
1509 << sqrt(mean2 / n - mean * mean / n / n) << endl;
1510 }
1511 }
1512 else
1513 {
1515 tmp[0] = m_numLocalCoeffs;
1516 vRowComm->Send(0, tmp);
1517 vRowComm->Block();
1518 tmp[0] = m_numLocalBndCoeffs;
1519 vRowComm->Send(0, tmp);
1520 }
1521
1522 // Either we have no more levels in the static condensation, or we
1523 // are not multi-level.
1525 {
1526 return;
1527 }
1528
1529 int level = 2;
1531 while (tmp->m_nextLevelLocalToGlobalMap)
1532 {
1533 tmp = tmp->m_nextLevelLocalToGlobalMap;
1534 ++level;
1535 }
1536
1537 // Print out multi-level static condensation information.
1538 if (n > 1)
1539 {
1540 if (isRoot)
1541 {
1542 NekDouble mean = level, mean2 = mean * mean;
1543 int minval = level, maxval = level;
1544
1545 Array<OneD, NekDouble> tmpRecv(1);
1546 for (i = 1; i < n; ++i)
1547 {
1548 vRowComm->Recv(i, tmpRecv);
1549 mean += tmpRecv[0];
1550 mean2 += tmpRecv[0] * tmpRecv[0];
1551
1552 if (tmpRecv[0] > maxval)
1553 {
1554 maxval = (int)(tmpRecv[0] + 0.5);
1555 }
1556 if (tmpRecv[0] < minval)
1557 {
1558 minval = (int)(tmpRecv[0] + 0.5);
1559 }
1560 }
1561
1562 out << " - M-level sc. dist. (min/max/mean/dev) : " << minval
1563 << " " << maxval << " " << (mean / n) << " "
1564 << sqrt(mean2 / n - mean * mean / n / n) << endl;
1565 }
1566 else
1567 {
1568 Array<OneD, NekDouble> tmpSend(1);
1569 tmpSend[0] = level;
1570 vRowComm->Send(0, tmpSend);
1571 }
1572 }
1573 else
1574 {
1575 if (isRoot)
1576 {
1577 out << " - Number of static cond. levels : " << level
1578 << endl;
1579 }
1580 }
1581
1582 if (isRoot)
1583 {
1584 out << "Stats at lowest static cond. level:" << endl;
1585 }
1586 tmp->PrintStats(out, variable, false);
1587}
1588} // namespace Nektar::MultiRegions
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:47
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Base class for constructing local to global mapping of degrees of freedom.
Definition: AssemblyMap.h:56
void PatchLocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
const Array< OneD, const int > & GetExtraDirEdges()
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:399
int GetNumPatches() const
Returns the number of patches in this static condensation level.
const Array< OneD, const int > & GetGlobalToUniversalMapUnique()
int GetNumGlobalCoeffs() const
Returns the total number of global coefficients.
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:361
Array< OneD, int > m_bndCondCoeffsToLocalTraceMap
Integer map of bnd cond coeff to local trace coeff.
Definition: AssemblyMap.h:390
const Array< OneD, const unsigned int > & GetNumLocalBndCoeffsPerPatch()
Returns the number of local boundary coefficients in each patch.
Array< OneD, int > m_bndCondCoeffsToLocalCoeffsMap
Integer map of bnd cond coeffs to local coefficients.
Definition: AssemblyMap.h:386
virtual const Array< OneD, const int > & v_GetExtraDirEdges()
int GetNumLocalCoeffs() const
Returns the total number of local coefficients.
std::string GetVariable()
Retrieves the variable string.
const Array< OneD, const unsigned int > & GetNumLocalIntCoeffsPerPatch()
Returns the number of local interior coefficients in each patch.
virtual void v_Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:375
void Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
NekDouble GetLocalToGlobalBndSign(const int i) const
Retrieve the sign change of a given local boundary mode.
virtual std::shared_ptr< AssemblyMap > v_LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
Generate a linear space mapping from existing mapping.
void AssembleBnd(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const
PatchMapSharedPtr m_patchMapFromPrevLevel
Mapping information for previous level in MultiLevel Solver.
Definition: AssemblyMap.h:505
Array< OneD, int > m_localToLocalIntMap
Integer map of local boundary coeffs to local interior system numbering.
Definition: AssemblyMap.h:384
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique()
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:372
std::string m_linSysIterSolver
Iterative solver: Conjugate Gradient, GMRES.
Definition: AssemblyMap.h:414
bool AtLastLevel() const
Returns true if this is the last level in the multi-level static condensation.
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
const Array< OneD, const int > & GetLocalToGlobalMap()
const Array< OneD, const int > & GetBndCondCoeffsToLocalTraceMap()
Retrieves the local indices corresponding to the boundary expansion modes to global trace.
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap()
Array< OneD, int > m_globalToUniversalBndMap
Integer map of process coeffs to universal space.
Definition: AssemblyMap.h:394
void PatchAssemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:411
std::shared_ptr< AssemblyMap > LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
virtual int v_GetNumNonDirVertexModes() const
void CalculateBndSystemBandWidth()
Calculates the bandwidth of the boundary system.
void LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm=true) const
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:380
void UniversalAbsMaxBnd(Array< OneD, NekDouble > &bndvals)
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:427
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap()
const Array< OneD, const int > & GetGlobalToUniversalMap()
std::string GetLinSysIterSolver() const
virtual int v_GetNumNonDirEdgeModes() const
virtual int v_GetNumDirFaces() const
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:330
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign() const
int GetNumGlobalDirBndCoeffs() const
Returns the number of global Dirichlet boundary coefficients.
const Array< OneD, NekDouble > & GetBndCondCoeffsToLocalCoeffsSign()
Returns the modal sign associated with a given boundary expansion mode.
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:342
void LocalToLocalBnd(const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locbnd) const
const Array< OneD, const int > & GetBndCondCoeffsToLocalCoeffsMap()
Retrieves the local indices corresponding to the boundary expansion modes.
void PrintStats(std::ostream &out, std::string variable, bool printHeader=true) const
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:434
const Array< OneD, const int > & GetLocalToGlobalBndMap()
Retrieve the global indices of the local boundary modes.
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:423
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:401
const Array< OneD, const int > & GetGlobalToUniversalBndMapUnique()
void GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
GlobalSysSolnType GetGlobalSysSolnType() const
Returns the method of solving global systems.
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:346
bool GetSingularSystem() const
Retrieves if the system is singular (true) or not (false)
std::string m_variable
Variable string identifier.
Definition: AssemblyMap.h:336
virtual int v_GetNumNonDirFaceModes() const
void PatchGlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
virtual int v_GetFullSystemBandWidth() const
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:348
LibUtilities::CommSharedPtr GetComm()
Retrieves the communicator.
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:429
const AssemblyMapSharedPtr GetNextLevelLocalToGlobalMap() const
Returns the local to global mapping for the next level in the multi-level static condensation.
int GetBndSystemBandWidth() const
Returns the bandwidth of the boundary system.
virtual void v_UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
int GetStaticCondLevel() const
Returns the level of static condensation for this map.
void GlobalToLocalBndWithoutSign(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc)
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:350
Array< OneD, int > m_localToGlobalBndMap
Integer map of local coeffs to global Boundary Dofs.
Definition: AssemblyMap.h:378
size_t GetHash() const
Retrieves the hash of this map.
int GetLocalToGlobalBndMap(const int i) const
Retrieve the global index of a given local boundary mode.
virtual ~AssemblyMap()
Destructor.
Gs::gs_data * m_dirBndGsh
gs gather communication to impose Dirhichlet BCs.
Definition: AssemblyMap.h:419
virtual void v_GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
const Array< OneD, NekDouble > & GetLocalToGlobalSign() const
int GetNumLocalDirBndCoeffs() const
Returns the number of local Dirichlet boundary coefficients.
Array< OneD, const NekDouble > GetLocalToGlobalBndSign() const
Retrieve the sign change for all local boundary modes.
void LocalBndToLocal(const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
AssemblyMap()
Default constructor.
Definition: AssemblyMap.cpp:79
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:396
virtual int v_GetNumNonDirEdges() const
virtual int v_GetNumDirEdges() const
const PatchMapSharedPtr & GetPatchMapFromPrevLevel(void) const
Returns the patch map from the previous level of the multi-level static condensation.
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:333
int GetNumGlobalBndCoeffs() const
Returns the total number of global boundary coefficients.
const Array< OneD, const int > & GetGlobalToUniversalBndMap()
Array< OneD, int > m_localToLocalBndMap
Integer map of local boundary coeffs to local boundary system numbering.
Definition: AssemblyMap.h:382
Array< OneD, NekDouble > m_bndCondCoeffsToLocalCoeffsSign
Integer map of sign of bnd cond coeffs to local coefficients.
Definition: AssemblyMap.h:388
const Array< OneD, const int > & GetBndCondIDToGlobalTraceID()
void LocalBndToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset, bool UseComm=true) const
void UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:425
void LocalToLocalInt(const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locint) const
virtual int v_GetNumNonDirFaces() const
bool GetSignChange()
Returns true if using a modal expansion requiring a change of sign of some modes.
void LocalIntToLocal(const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
void GlobalToLocalBnd(const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
std::string m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:404
int GetNumLocalBndCoeffs() const
Returns the total number of local boundary coefficients.
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:344
Array< OneD, int > m_bndCondIDToGlobalTraceID
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:392
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:99
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:217
static void Gather(Nektar::Array< OneD, NekDouble > pU, gs_op pOp, gs_data *pGsh, Nektar::Array< OneD, NekDouble > pBuffer=NullNekDouble1DArray)
Performs a gather-scatter operation of the provided values.
Definition: GsLib.hpp:278
@ gs_amax
Definition: GsLib.hpp:64
@ gs_add
Definition: GsLib.hpp:60
@ gs_max
Definition: GsLib.hpp:63
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:55
@ eNoSolnType
No Solution type specified.
std::shared_ptr< BottomUpSubStructuredGraph > BottomUpSubStructuredGraphSharedPtr
std::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:50
std::shared_ptr< PatchMap > PatchMapSharedPtr
int RoundNekDoubleToInt(NekDouble x)
Rounds a double precision number to an integer.
Definition: AssemblyMap.cpp:56
double NekDouble
void Gathr(I n, const T *x, const I *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.hpp:507
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.hpp:539
void Assmb(int n, const T *x, const int *y, T *z)
Assemble z[y[i]] += x[i]; z should be zero'd first.
Definition: Vmath.hpp:577
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
STL namespace.
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285