MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_CoupledRBMFactory_def.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // MueLu: A package for multigrid based preconditioning
4 //
5 // Copyright 2012 NTESS and the MueLu contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef MUELU_COUPLEDRBMFACTORY_DEF_HPP
11 #define MUELU_COUPLEDRBMFACTORY_DEF_HPP
12 
13 #include <Xpetra_Matrix.hpp>
14 #include <Xpetra_MultiVectorFactory.hpp>
15 #include <Xpetra_VectorFactory.hpp>
16 
17 #include "MueLu_CoupledRBMFactory_decl.hpp"
18 #include "MueLu_Level.hpp"
19 #include "MueLu_Monitor.hpp"
20 
21 namespace MueLu {
22 
23 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
24 CoupledRBMFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::~CoupledRBMFactory() {}
25 
26 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
27 void CoupledRBMFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::DeclareInput(Level& currentLevel) const {
28  if (currentLevel.IsAvailable(nspName_, NoFactory::get()) == false && currentLevel.GetLevelID() == 0) {
29  Input(currentLevel, "A");
30  // Input(currentLevel,"Coordinates");
31  }
32  if (currentLevel.GetLevelID() != 0) {
33  currentLevel.DeclareInput("Nullspace", GetFactory(nspName_).get(), this); /* ! "Nullspace" and nspName_ mismatch possible here */
34  }
35 }
36 
37 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
38 void CoupledRBMFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(Level& currentLevel) const {
39  FactoryMonitor m(*this, "Structural acoustics nullspace factory", currentLevel);
40  RCP<MultiVector> nullspace;
41  if (currentLevel.GetLevelID() == 0) {
42  if (currentLevel.IsAvailable(nspName_, NoFactory::get())) {
43  nullspace = currentLevel.Get<RCP<MultiVector> >(nspName_, NoFactory::get());
44  GetOStream(Runtime1) << "Use user-given rigid body modes " << nspName_ << ": nullspace dimension=" << nullspace->getNumVectors() << " nullspace length=" << nullspace->getGlobalLength() << std::endl;
45  } else {
46  RCP<Matrix> A = Get<RCP<Matrix> >(currentLevel, "A");
47  RCP<MultiVector> Coords = Get<RCP<MultiVector> >(currentLevel, "Coordinates");
48  GetOStream(Runtime1) << "Generating nullspace for structural acoustics: dimension = " << numPDEs_ << std::endl;
49  RCP<const Map> xmap = A->getDomainMap();
50  nullspace = MultiVectorFactory::Build(xmap, 6);
51  Scalar zero = (Scalar)0.0;
52  nullspace->putScalar(zero);
53  ArrayRCP<Scalar> xnodes, ynodes, znodes;
54  Scalar cx, cy, cz;
55  ArrayRCP<Scalar> nsValues0, nsValues1, nsValues2, nsValues3, nsValues4, nsValues5;
56  int nDOFs = xmap->getLocalNumElements();
57  xnodes = Coords->getDataNonConst(0);
58  ynodes = Coords->getDataNonConst(1);
59  znodes = Coords->getDataNonConst(2);
60  cx = Coords->getVector(0)->meanValue();
61  cy = Coords->getVector(1)->meanValue();
62  cz = Coords->getVector(2)->meanValue();
63  nsValues0 = nullspace->getDataNonConst(0);
64  nsValues1 = nullspace->getDataNonConst(1);
65  nsValues2 = nullspace->getDataNonConst(2);
66  nsValues3 = nullspace->getDataNonConst(3);
67  nsValues4 = nullspace->getDataNonConst(4);
68  nsValues5 = nullspace->getDataNonConst(5);
69  for (int j = 0; j < nDOFs; j += numPDEs_) {
70  Scalar one = (Scalar)1.0;
71  if (xmap->getGlobalElement(j) >= lastAcousticDOF_) {
72  Scalar xdiff = xnodes[j] - cx;
73  Scalar ydiff = ynodes[j] - cy;
74  Scalar zdiff = znodes[j] - cz;
75  // translation
76  nsValues0[j + 0] = one;
77  nsValues1[j + 1] = one;
78  nsValues2[j + 2] = one;
79  // rotate around z-axis (x-y plane)
80  nsValues3[j + 0] = -ydiff;
81  nsValues3[j + 1] = xdiff;
82  // rotate around x-axis (y-z plane)
83  nsValues4[j + 1] = -zdiff;
84  nsValues4[j + 2] = ydiff;
85  // rotate around y-axis (x-z plane)
86  nsValues5[j + 0] = zdiff;
87  nsValues5[j + 2] = -xdiff;
88  } else {
89  // translation
90  nsValues0[j + 0] = one;
91  // insert random values and keep the top row for this node empty
92  nsValues1[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
93  nsValues1[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
94  nsValues2[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
95  nsValues2[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
96  nsValues3[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
97  nsValues3[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
98  nsValues4[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
99  nsValues4[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
100  nsValues5[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
101  nsValues5[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
102  }
103  }
104  } // end if "Nullspace" not available
105  } else {
106  nullspace = currentLevel.Get<RCP<MultiVector> >("Nullspace", GetFactory(nspName_).get());
107  }
108  Set(currentLevel, "Nullspace", nullspace);
109 }
110 
111 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
112 void CoupledRBMFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::BuildRBM(RCP<Matrix>& A, RCP<MultiVector>& Coords, RCP<MultiVector>& nullspace) const {
113  GetOStream(Runtime1) << "Generating nullspace for structural acoustics: dimension = " << numPDEs_ << std::endl;
114  RCP<const Map> xmap = A->getDomainMap();
115  nullspace = MultiVectorFactory::Build(xmap, 6);
116  Scalar zero = (Scalar)0.0;
117  nullspace->putScalar(zero);
118  ArrayRCP<Scalar> xnodes, ynodes, znodes;
119  Scalar cx, cy, cz;
120  ArrayRCP<Scalar> nsValues0, nsValues1, nsValues2, nsValues3, nsValues4, nsValues5;
121  int nDOFs = xmap->getLocalNumElements();
122  xnodes = Coords->getDataNonConst(0);
123  ynodes = Coords->getDataNonConst(1);
124  znodes = Coords->getDataNonConst(2);
125  cx = Coords->getVector(0)->meanValue();
126  cy = Coords->getVector(1)->meanValue();
127  cz = Coords->getVector(2)->meanValue();
128  nsValues0 = nullspace->getDataNonConst(0);
129  nsValues1 = nullspace->getDataNonConst(1);
130  nsValues2 = nullspace->getDataNonConst(2);
131  nsValues3 = nullspace->getDataNonConst(3);
132  nsValues4 = nullspace->getDataNonConst(4);
133  nsValues5 = nullspace->getDataNonConst(5);
134  for (int j = 0; j < nDOFs; j += numPDEs_) {
135  Scalar one = (Scalar)1.0;
136  if (xmap->getGlobalElement(j) >= lastAcousticDOF_) {
137  Scalar xdiff = xnodes[j] - cx;
138  Scalar ydiff = ynodes[j] - cy;
139  Scalar zdiff = znodes[j] - cz;
140  // translation
141  nsValues0[j + 0] = one;
142  nsValues1[j + 1] = one;
143  nsValues2[j + 2] = one;
144  // rotate around z-axis (x-y plane)
145  nsValues3[j + 0] = -ydiff;
146  nsValues3[j + 1] = xdiff;
147  // rotate around x-axis (y-z plane)
148  nsValues4[j + 1] = -zdiff;
149  nsValues4[j + 2] = ydiff;
150  // rotate around y-axis (x-z plane)
151  nsValues5[j + 0] = zdiff;
152  nsValues5[j + 2] = -xdiff;
153  } else {
154  // translation
155  nsValues0[j + 0] = one;
156  // insert random values and keep the top row for this node empty
157  nsValues1[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
158  nsValues1[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
159  nsValues2[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
160  nsValues2[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
161  nsValues3[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
162  nsValues3[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
163  nsValues4[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
164  nsValues4[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
165  nsValues5[j + 1] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
166  nsValues5[j + 2] = (Scalar)(((double)rand()) / ((double)RAND_MAX));
167  }
168  }
169 }
170 
171 } // namespace MueLu
172 
173 #define MUELU_COUPLEDRBMFACTORY_SHORT
174 #endif // MUELU_COUPLEDRBMFACTORY_DEF_HPP
MueLu::Level::Get
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level-&gt;Get&lt; RCP&lt;Matrix&gt; &gt;(&quot;A&quot;, factory) if factory == NULL =&gt; use default factory.
Definition: MueLu_Level.hpp:154
MueLu::CoupledRBMFactory::DeclareInput
void DeclareInput(Level &currentLevel) const
Specifies the data that this class needs, and the factories that generate that data.
Definition: MueLu_CoupledRBMFactory_def.hpp:27
MueLu::CoupledRBMFactory::Build
void Build(Level &currentLevel) const
Build an object with this factory.
Definition: MueLu_CoupledRBMFactory_def.hpp:38
MueLu::FactoryMonitor
Timer to be used in factories. Similar to Monitor but with additional timers.
Definition: MueLu_Monitor.hpp:165
MueLu_Level.hpp
Xpetra_Matrix.hpp
MueLu::NoFactory::get
static const NoFactory * get()
Definition: MueLu_NoFactory.cpp:34
MueLu::CoupledRBMFactory::BuildRBM
void BuildRBM(RCP< Matrix > &A, RCP< MultiVector > &Coords, RCP< MultiVector > &nullspace) const
Definition: MueLu_CoupledRBMFactory_def.hpp:112
Scalar
MueLu::DefaultScalar Scalar
Definition: MueLu_UseDefaultTypes.hpp:12
MueLu::Level
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:63
MueLu::CoupledRBMFactory::~CoupledRBMFactory
virtual ~CoupledRBMFactory()
Destructor.
Definition: MueLu_CoupledRBMFactory_def.hpp:24
MueLu_Monitor.hpp
Xpetra_VectorFactory.hpp
MueLu::Level::GetLevelID
int GetLevelID() const
Return level number.
Definition: MueLu_Level.cpp:51
MueLu::Runtime1
Description of what is happening (more verbose)
Definition: MueLu_VerbosityLevel.hpp:26
MueLu::Level::DeclareInput
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
Definition: MueLu_Level.cpp:156
MueLu_CoupledRBMFactory_decl.hpp
MueLu::Level::IsAvailable
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need&#39;s value has been saved.
Definition: MueLu_Level.cpp:281