escript  Revision_
dudley/src/Assemble.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16 
17 /****************************************************************************
18 
19  Assemblage routines: header file
20 
21 *****************************************************************************/
22 
23 #ifndef __DUDLEY_ASSEMBLE_H__
24 #define __DUDLEY_ASSEMBLE_H__
25 
26 #include "Dudley.h"
27 #include "ElementFile.h"
28 #include "NodeFile.h"
29 #include <escript/AbstractSystemMatrix.h>
30 
31 namespace dudley {
32 
34 {
35  AssembleParameters(const NodeFile* nodes, const ElementFile* ef,
37  bool reducedOrder);
38 
46  int numQuad;
48  int numDim;
50  int NN;
52  int numEqu;
54  const index_t* DOF;
59  int numShapes;
60  const double* shapeFns;
61 };
62 
63 void Assemble_PDE(const NodeFile* nodes, const ElementFile* elements,
65  const escript::Data& A, const escript::Data& B,
66  const escript::Data& C, const escript::Data& D,
67  const escript::Data& X, const escript::Data& Y);
68 
69 template<typename Scalar = double>
71  const escript::Data& d_dirac,
72  const escript::Data& y_dirac);
73 
74 template<typename Scalar = double>
76  const escript::Data& A, const escript::Data& B,
77  const escript::Data& C, const escript::Data& D,
78  const escript::Data& X, const escript::Data& Y);
79 
80 template<typename Scalar = double>
82  const escript::Data& A, const escript::Data& B,
83  const escript::Data& C, const escript::Data& D,
84  const escript::Data& X, const escript::Data& Y);
85 
86 template<typename Scalar = double>
88  const escript::Data& A, const escript::Data& B,
89  const escript::Data& C, const escript::Data& D,
90  const escript::Data& X, const escript::Data& Y);
91 
92 template<typename Scalar = double>
94  const escript::Data& A, const escript::Data& B,
95  const escript::Data& C, const escript::Data& D,
96  const escript::Data& X, const escript::Data& Y);
97 
98 
104 template<typename Scalar>
106  const std::vector<index_t>& Nodes, int numEq,
107  const std::vector<Scalar>& array);
108 
112 void Assemble_LumpedSystem(const NodeFile* nodes, const ElementFile* elements,
113  escript::Data& lumpedMat, const escript::Data& D,
114  bool useHRZ);
115 
118  escript::Data& out, const escript::Data& in);
119 
122  const escript::Data& in);
123 
125 void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
126  const escript::Data& in);
127 
129 void Assemble_NodeCoordinates(const NodeFile* nodes, escript::Data& x);
130 
132 void Assemble_getNormal(const NodeFile* nodes, const ElementFile* elements,
133  escript::Data& normals);
134 
137 void Assemble_getSize(const NodeFile* nodes, const ElementFile* elements,
138  escript::Data& size);
139 
142 void Assemble_gradient(const NodeFile* nodes, const ElementFile* elements,
143  escript::Data& gradient, const escript::Data& data);
144 
146 void Assemble_integrate(const NodeFile* nodes, const ElementFile* elements,
147  const escript::Data& data, std::vector<double>& integrals);
148 
150 void Assemble_interpolate(const NodeFile* nodes, const ElementFile* elements,
151  const escript::Data& data, escript::Data& output);
152 
153 void Assemble_jacobians_2D(const double* coordinates, int numQuad,
154  dim_t numElements, int numNodes,
155  const index_t* nodes, double* dTdX, double* absD,
156  double* quadWeight, const index_t* elementId);
157 
158 void Assemble_jacobians_2D_M1D_E1D(const double* coordinates, int numQuad,
159  dim_t numElements, int numNodes,
160  const index_t* nodes, double* dTdX, double* absD,
161  double* quadWeight, const index_t* elementId);
162 
163 void Assemble_jacobians_3D(const double* coordinates, int numQuad,
164  dim_t numElements, int numNodes,
165  const index_t* nodes, double* dTdX, double* abs_D,
166  double* quadWeight, const index_t* elementId);
167 
168 void Assemble_jacobians_3D_M2D_E2D(const double* coordinates, int numQuad,
169  dim_t numElements, int numNodes,
170  const index_t* nodes, double* dTdX, double* absD,
171  double* quadWeight, const index_t* elementId);
172 
173 
174 } // namespace dudley
175 
176 #endif // __DUDLEY_ASSEMBLE_H__
177 
int numShapes
Definition: dudley/src/Assemble.h:59
int numEqu
number of equations (= matrix row/column block size)
Definition: dudley/src/Assemble.h:52
void Assemble_integrate(const NodeFile *nodes, const ElementFile *elements, const escript::Data &data, std::vector< double > &integrals)
integrates data on quadrature points
Definition: dudley/src/Assemble_integrate.cpp:24
Definition: Dudley.h:52
void Assemble_gradient(const NodeFile *nodes, const ElementFile *elements, escript::Data &gradient, const escript::Data &data)
Definition: dudley/src/Assemble_gradient.cpp:27
void Assemble_jacobians_2D(const double *coordinates, int numQuad, dim_t numElements, int numNodes, const index_t *nodes, double *dTdX, double *absD, double *quadWeight, const index_t *elementId)
Definition: dudley/src/Assemble_jacobians.cpp:55
void Assemble_PDE(const NodeFile *nodes, const ElementFile *elements, escript::ASM_ptr S, escript::Data &F, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y)
Definition: dudley/src/Assemble_PDE.cpp:83
void Assemble_PDE_Points(const AssembleParameters &p, const escript::Data &d_dirac, const escript::Data &y_dirac)
Definition: dudley/src/Assemble_PDE_Points.cpp:45
dim_t DOF_UpperBound
number of local degrees of freedom
Definition: dudley/src/Assemble.h:56
escript::AbstractSystemMatrix * S
system matrix to be updated
Definition: dudley/src/Assemble.h:42
Definition: dudley/src/Assemble.h:33
void Assemble_CopyElementData(const ElementFile *elements, escript::Data &out, const escript::Data &in)
copies data between different types of elements
Definition: dudley/src/Assemble_CopyElementData.cpp:23
void Assemble_getNormal(const NodeFile *nodes, const ElementFile *elements, escript::Data &normals)
calculates the normal vector at quadrature points on face elements
Definition: dudley/src/Assemble_getNormal.cpp:25
int index_t
type for array/matrix indices used both globally and on each rank
Definition: DataTypes.h:59
void Assemble_jacobians_2D_M1D_E1D(const double *coordinates, int numQuad, dim_t numElements, int numNodes, const index_t *nodes, double *dTdX, double *absD, double *quadWeight, const index_t *elementId)
Definition: dudley/src/Assemble_jacobians.cpp:133
void Assemble_PDE_Single_2D(const AssembleParameters &p, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y)
Definition: dudley/src/Assemble_PDE_Single_2D.cpp:45
const ElementFile * elements
element file these parameters apply to
Definition: dudley/src/Assemble.h:40
int numQuad
number of quadrature nodes
Definition: dudley/src/Assemble.h:46
Data represents a collection of datapoints.
Definition: Data.h:63
const index_t * DOF
row and column degrees of freedom
Definition: dudley/src/Assemble.h:54
Definition: dudley/src/ElementFile.h:51
int NN
leading dimension of element node table
Definition: dudley/src/Assemble.h:50
A suite of factory methods for creating 2D and 3D dudley domains.
Definition: dudley/src/Assemble.h:31
void Assemble_NodeCoordinates(const NodeFile *nodes, escript::Data &x)
copies node coordinates into expanded Data object x
Definition: dudley/src/Assemble_NodeCoordinates.cpp:26
void Assemble_CopyNodalData(const NodeFile *nodes, escript::Data &out, const escript::Data &in)
copies data between different types of nodal representations
Definition: dudley/src/Assemble_CopyNodalData.cpp:22
void Assemble_PDE_System_2D(const AssembleParameters &p, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y)
Definition: dudley/src/Assemble_PDE_System_2D.cpp:48
void Assemble_jacobians_3D_M2D_E2D(const double *coordinates, int numQuad, dim_t numElements, int numNodes, const index_t *nodes, double *dTdX, double *absD, double *quadWeight, const index_t *elementId)
Definition: dudley/src/Assemble_jacobians.cpp:254
AssembleParameters(const NodeFile *nodes, const ElementFile *ef, escript::ASM_ptr sm, escript::Data &rhs, bool reducedOrder)
Definition: dudley/src/Assemble_getAssembleParameters.cpp:34
void Assemble_LumpedSystem(const NodeFile *nodes, const ElementFile *elements, escript::Data &lumpedMat, const escript::Data &D, bool useHRZ)
Definition: dudley/src/Assemble_LumpedSystem.cpp:25
escript::Data & F
right-hand side to be updated
Definition: dudley/src/Assemble.h:44
int numDim
number of spatial dimensions
Definition: dudley/src/Assemble.h:48
void Assemble_PDE_System_3D(const AssembleParameters &p, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y)
Definition: dudley/src/Assemble_PDE_System_3D.cpp:48
const ElementFile_Jacobians * jac
reference to jacobians
Definition: dudley/src/Assemble.h:58
void Assemble_jacobians_3D(const double *coordinates, int numQuad, dim_t numElements, int numNodes, const index_t *nodes, double *dTdX, double *abs_D, double *quadWeight, const index_t *elementId)
Definition: dudley/src/Assemble_jacobians.cpp:179
Definition: dudley/src/ElementFile.h:27
Base class for escript system matrices.
Definition: AbstractSystemMatrix.h:42
void Assemble_AverageElementData(const ElementFile *elements, escript::Data &out, const escript::Data &in)
averages data
Definition: dudley/src/Assemble_AverageElementData.cpp:25
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:32
void Assemble_interpolate(const NodeFile *nodes, const ElementFile *elements, const escript::Data &data, escript::Data &output)
interpolates nodal data in a data array onto elements (=integration points)
Definition: dudley/src/Assemble_interpolate.cpp:25
void Assemble_addToSystemMatrix(escript::AbstractSystemMatrix *S, const std::vector< index_t > &Nodes, int numEq, const std::vector< Scalar > &array)
Definition: dudley/src/NodeFile.h:38
void Assemble_getSize(const NodeFile *nodes, const ElementFile *elements, escript::Data &size)
Definition: dudley/src/Assemble_getSize.cpp:24
void Assemble_PDE_Single_3D(const AssembleParameters &p, const escript::Data &A, const escript::Data &B, const escript::Data &C, const escript::Data &D, const escript::Data &X, const escript::Data &Y)
Definition: dudley/src/Assemble_PDE_Single_3D.cpp:45
const double * shapeFns
Definition: dudley/src/Assemble.h:60
index_t dim_t
Definition: DataTypes.h:64