My Project
asmhandler.hpp
Go to the documentation of this file.
1 //==============================================================================
11 //==============================================================================
12 #ifndef ASMHANDLER_HPP_
13 #define ASMHANDLER_HPP_
14 
15 #include <opm/common/utility/platform_dependent/disable_warnings.h>
16 
17 #include <dune/geometry/referenceelements.hh>
18 #include <dune/common/fmatrix.hh>
19 #include <dune/istl/bvector.hh>
20 #include <dune/common/fvector.hh>
21 
22 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
23 
24 
26 #include <opm/elasticity/mpc.hh>
28 
29 namespace Opm {
30 namespace Elasticity {
31 
32 
34  template<class GridType>
35 class ASMHandler {
36  public:
38  static const int dim = GridType::dimension;
39 
41  typedef typename GridType::LeafGridView::IndexSet LeafIndexSet;
42 
44  typedef typename GridType::LeafGridView::template Codim<0>::Iterator LeafIterator;
45 
48  ASMHandler(const GridType& gv_) : gv(gv_), maxeqn(0)
49  {
50  }
51 
54  {
55  for (MPCMap::iterator it=mpcs.begin(); it != mpcs.end();++it)
56  delete it->second;
57  }
58 
60  typedef Dune::FieldVector<double,dim> NodeValue;
61 
64  size_t getEqns() const
65  {
66  return maxeqn;
67  }
68 
73  int getEquationForDof(int node, int dof)
74  {
75  return meqn[node*dim+dof];
76  }
77 
81  {
82  return A;
83  }
84 
88  {
89  return b;
90  }
91 
94  void initForAssembly();
95 
102  template<int esize>
103  void addElement(const Dune::FieldMatrix<double,esize,esize>* K,
104  const Dune::FieldVector<double,esize>* S,
105  const LeafIterator& cell,
106  Vector* b=NULL);
107 
112  template<int comp>
113  void extractValues(Dune::FieldVector<double,comp>& v,
114  const Vector& u, const LeafIterator& it);
115 
117  void expandSolution(Vector& result, const Vector& u);
118 
122  void addMPC(MPC* mpc);
123 
128  MPC* getMPC(int node, int dof);
129 
133  void updateFixedNode(int node,
134  const std::pair<Direction,NodeValue>& entry);
135 
138  bool isFixed(int node)
139  {
140  return (fixedNodes.find(node) != fixedNodes.end());
141  }
142 
144  void printOperator() const;
145 
147  void printLoadVector() const;
148 
152  {
153  return adjacencyPattern;
154  }
155  protected:
158  {
159  for (MPCMap::iterator it = mpcs.begin();
160  it != mpcs.end();++it)
161  resolveMPCChain(it->second);
162  }
163 
166  void resolveMPCChain(MPC* mpc);
167 
169  void preprocess();
170 
175  void nodeAdjacency(const LeafIterator& it, int vertexsize, int row);
176 
179 
189  template<int esize>
190  void addDOF(int row, int erow,
191  const Dune::FieldMatrix<double,esize,esize>* K,
192  const Dune::FieldVector<double,esize>* S,
193  const LeafIndexSet& set,
194  const LeafIterator& cell,
195  Vector* b,
196  double scale=1.f);
197 
200 
202  std::vector<int> meqn;
203 
205  typedef std::pair<Direction,NodeValue> fixEntry;
206 
208  typedef std::map<int,fixEntry> fixMap;
209 
211  typedef typename fixMap::iterator fixIt;
212 
215 
218 
221 
224 
226  const GridType& gv;
227 
229  size_t maxeqn;
230  private:
232  ASMHandler(const ASMHandler&) {}
234  ASMHandler& operator=(const ASMHandler&) {}
235 };
236 
237 }
238 }
239 
240 #include "asmhandler_impl.hpp"
241 
242 #endif
Class handling finite element assembly - template implementations.
Class handling finite element assembly.
Definition: asmhandler.hpp:35
static const int dim
The dimension of the grid.
Definition: asmhandler.hpp:38
std::pair< Direction, NodeValue > fixEntry
Fixed nodes.
Definition: asmhandler.hpp:205
bool isFixed(int node)
Check if a node is marked as fixed (in any direction)
Definition: asmhandler.hpp:138
void initForAssembly()
This function needs to be called before starting the element assembly.
Definition: asmhandler_impl.hpp:22
MPCMap mpcs
The set of MPC.
Definition: asmhandler.hpp:199
void resolveMPCChain(MPC *mpc)
Internal function.
Definition: asmhandler_impl.hpp:250
~ASMHandler()
Destructor.
Definition: asmhandler.hpp:53
const GridType & gv
A reference to the grid in use.
Definition: asmhandler.hpp:226
AdjacencyPattern & getAdjacencyPattern()
Access current adjacency pattern.
Definition: asmhandler.hpp:151
fixMap fixedNodes
The map holding information about our fixed nodes.
Definition: asmhandler.hpp:214
Dune::FieldVector< double, dim > NodeValue
A vectorial node value.
Definition: asmhandler.hpp:60
size_t maxeqn
The number of equations in the system.
Definition: asmhandler.hpp:229
AdjacencyPattern adjacencyPattern
Holds the adjacency pattern of the sparse matrix.
Definition: asmhandler.hpp:217
std::vector< int > meqn
Vector of (interleaved) dof<->eqn mapping.
Definition: asmhandler.hpp:202
fixMap::iterator fixIt
Iterator over a fixmap.
Definition: asmhandler.hpp:211
GridType::LeafGridView::IndexSet LeafIndexSet
A set of indices.
Definition: asmhandler.hpp:41
int getEquationForDof(int node, int dof)
Get the equation number for a given dof on a given node.
Definition: asmhandler.hpp:73
void addElement(const Dune::FieldMatrix< double, esize, esize > *K, const Dune::FieldVector< double, esize > *S, const LeafIterator &cell, Vector *b=NULL)
Add an element matrix/vector to the system.
Definition: asmhandler_impl.hpp:85
void resolveMPCChains()
Resolve chained MPCs.
Definition: asmhandler.hpp:157
Vector b
The load vector.
Definition: asmhandler.hpp:223
ASMHandler(const GridType &gv_)
The default constructor.
Definition: asmhandler.hpp:48
std::map< int, fixEntry > fixMap
A mapping from dof to fix value info.
Definition: asmhandler.hpp:208
Matrix A
The linear operator.
Definition: asmhandler.hpp:220
void updateFixedNode(int node, const std::pair< Direction, NodeValue > &entry)
Update/add a fixed node.
Definition: asmhandler_impl.hpp:214
void extractValues(Dune::FieldVector< double, comp > &v, const Vector &u, const LeafIterator &it)
Extract values corresponding to cell.
Definition: asmhandler_impl.hpp:114
void expandSolution(Vector &result, const Vector &u)
Expand a system vector to a solution vector.
Definition: asmhandler_impl.hpp:145
void addMPC(MPC *mpc)
Add a MPC.
Definition: asmhandler_impl.hpp:187
Matrix & getOperator()
Obtain a reference to the linear operator.
Definition: asmhandler.hpp:80
void preprocess()
Internal function. Generate meqn for registered MPC/fixed nodes.
Definition: asmhandler_impl.hpp:315
void nodeAdjacency(const LeafIterator &it, int vertexsize, int row)
Internal function.
Definition: asmhandler_impl.hpp:351
GridType::LeafGridView::template Codim< 0 >::Iterator LeafIterator
An iterator over grid cells.
Definition: asmhandler.hpp:44
void printLoadVector() const
Print the current load vector.
Definition: asmhandler_impl.hpp:241
Vector & getLoadVector()
Obtain a reference to the load vector.
Definition: asmhandler.hpp:87
MPC * getMPC(int node, int dof)
Look for and return a MPC for a specified node+dof.
Definition: asmhandler_impl.hpp:205
size_t getEqns() const
Get the number of equations in the system.
Definition: asmhandler.hpp:64
void determineAdjacencyPattern()
Internal function. Calculate adjacency pattern.
Definition: asmhandler_impl.hpp:375
void addDOF(int row, int erow, const Dune::FieldMatrix< double, esize, esize > *K, const Dune::FieldVector< double, esize > *S, const LeafIndexSet &set, const LeafIterator &cell, Vector *b, double scale=1.f)
Internal function.
Definition: asmhandler_impl.hpp:52
void printOperator() const
Print the current operator.
Definition: asmhandler_impl.hpp:235
A class for representing a general multi-point constraint equation.
Definition: mpc.hh:59
Logging helper utilities.
Helper class with some matrix operations.
Dune::BCRSMatrix< Dune::FieldMatrix< double, 1, 1 > > Matrix
A sparse matrix holding our operator.
Definition: matrixops.hpp:27
Dune::BlockVector< Dune::FieldVector< double, 1 > > Vector
A vector holding our RHS.
Definition: matrixops.hpp:33
std::vector< std::set< int > > AdjacencyPattern
For storing matrix adjacency/sparsity patterns.
Definition: matrixops.hpp:30
Representation of multi-point constraint (MPC) equations.
std::map< int, MPC * > MPCMap
A mapping from dof to MPCs.
Definition: mpc.hh:159
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43