escript  Revision_
finley/src/Assemble.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2018 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 __FINLEY_ASSEMBLE_H__
24 #define __FINLEY_ASSEMBLE_H__
25 
26 #include "Finley.h"
27 #include "ElementFile.h"
28 #include "NodeFile.h"
29 #include <escript/AbstractSystemMatrix.h>
30 
31 namespace finley {
32 
33 struct AssembleParameters
34 {
35  AssembleParameters(const NodeFile* nodes, const ElementFile* ef,
37  bool reducedOrder);
38 
40  const ElementFile* elements;
46  int numQuadTotal;
50  int numSides;
52  int numSub;
54  int numDim;
56  int NN;
59 
60  int numEqu;
61  const index_t* row_DOF;
64  const int* row_node;
66  int row_numShapes;
67  int numComp;
68  const index_t* col_DOF;
71  const int* col_node;
73  int col_numShapes;
74 };
75 
76 
80 void Assemble_PDE(const NodeFile* nodes, const ElementFile* elements,
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>
88  const escript::Data& d_dirac,
89  const escript::Data& y_dirac);
90 
92  const escript::Data& A, const escript::Data& B,
93  const escript::Data& C, const escript::Data& D,
94  const escript::Data& X, const escript::Data& Y);
95 
96 template<typename Scalar>
98  const escript::Data& A, const escript::Data& B,
99  const escript::Data& C, const escript::Data& D,
100  const escript::Data& X, const escript::Data& Y);
101 
102 template<typename Scalar>
104  const escript::Data& A, const escript::Data& B,
105  const escript::Data& C, const escript::Data& D,
106  const escript::Data& X, const escript::Data& Y);
107 
108 template<typename Scalar>
110  const escript::Data& Y);
111 
113  const escript::Data& A, const escript::Data& B,
114  const escript::Data& C, const escript::Data& D,
115  const escript::Data& X, const escript::Data& Y);
116 
117 template<typename Scalar>
119  const escript::Data& A, const escript::Data& B,
120  const escript::Data& C, const escript::Data& D,
121  const escript::Data& X, const escript::Data& Y);
122 
123 template<typename Scalar>
125  const escript::Data& A, const escript::Data& B,
126  const escript::Data& C, const escript::Data& D,
127  const escript::Data& X, const escript::Data& Y);
128 
129 template<typename Scalar>
131  const escript::Data& Y);
132 
133 template<typename Scalar = double>
135  const index_t* Nodes_Equa, int num_Equa, int NN_Sol,
136  const index_t* Nodes_Sol, int num_Sol, const Scalar* array);
137 
138 void Assemble_LumpedSystem(const NodeFile* nodes, const ElementFile* elements,
139  escript::Data& lumpedMat, const escript::Data& D,
140  bool useHRZ);
141 
143 template<typename Scalar>
144 void Assemble_AverageElementData(const ElementFile* elements,
145  escript::Data& out, const escript::Data& in);
146 
148 template<typename Scalar>
149 void Assemble_CopyElementData(const ElementFile* elements, escript::Data& out,
150  const escript::Data& in);
151 
153 template<typename Scalar>
154 void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
155  const escript::Data& in);
156 
158 void Assemble_NodeCoordinates(const NodeFile* nodes, escript::Data& x);
159 
161 void Assemble_getNormal(const NodeFile* nodes, const ElementFile* elements,
162  escript::Data& normals);
163 
166 void Assemble_getSize(const NodeFile* nodes, const ElementFile* elements,
167  escript::Data& size);
168 
171 template<typename Scalar>
172 void Assemble_gradient(const NodeFile* nodes, const ElementFile* elements,
173  escript::Data& gradient, const escript::Data& data);
174 
176 template<typename Scalar>
177 void Assemble_integrate(const NodeFile* nodes, const ElementFile* elements,
178  const escript::Data& data, Scalar* integrals);
179 
181 template<typename Scalar>
182 void Assemble_interpolate(const NodeFile* nodes, const ElementFile* elements,
183  const escript::Data& data, escript::Data& output);
184 
185 void Assemble_jacobians_1D(const double* coordinates, int numQuad,
186  const double* QuadWeights, int numShape,
187  dim_t numElements, int numNodes, const index_t* nodes,
188  const double* DSDv, int numTest, const double* DTDv,
189  double* dTdX, double* volume, const index_t* elementId);
190 
191 void Assemble_jacobians_2D(const double* coordinates, int numQuad,
192  const double* QuadWeights, int numShape,
193  dim_t numElements, int numNodes, const index_t* nodes,
194  const double* DSDv, int numTest, const double* DTDv,
195  double* dTdX, double* volume, const index_t* elementId);
196 
197 void Assemble_jacobians_2D_M1D_E1D(const double* coordinates, int numQuad,
198  const double* QuadWeights, int numShape,
199  dim_t numElements, int numNodes, const index_t* nodes,
200  const double* DSDv, int numTest, const double* DTDv,
201  double* dTdX, double* volume, const index_t* elementId);
202 
203 void Assemble_jacobians_2D_M1D_E1D_C(const double* coordinates, int numQuad,
204  const double* QuadWeights, int numShape,
205  dim_t numElements, int numNodes, const index_t* nodes,
206  const double* DSDv, int numTest, const double* DTDv,
207  double* dTdX, double* volume, const index_t* elementId);
208 
209 void Assemble_jacobians_2D_M1D_E2D(const double* coordinates, int numQuad,
210  const double* QuadWeights, int numShape,
211  dim_t numElements, int numNodes, const index_t* nodes,
212  const double* DSDv, int numTest, const double* DTDv,
213  double* dTdX, double* volume, const index_t* elementId);
214 
215 void Assemble_jacobians_2D_M1D_E2D_C(const double* coordinates, int numQuad,
216  const double* QuadWeights, int numShape,
217  dim_t numElements, int numNodes, const index_t* nodes,
218  const double* DSDv, int numTest, const double* DTDv,
219  double* dTdX, double* volume, const index_t* elementId);
220 
221 void Assemble_jacobians_3D(const double* coordinates, int numQuad,
222  const double* QuadWeights, int numShape,
223  dim_t numElements, int numNodes, const index_t* nodes,
224  const double* DSDv, int numTest, const double* DTDv,
225  double* dTdX, double* volume, const index_t* elementId);
226 
227 void Assemble_jacobians_3D_M2D_E2D(const double* coordinates, int numQuad,
228  const double* QuadWeights, int numShape,
229  dim_t numElements, int numNodes, const index_t* nodes,
230  const double* DSDv, int numTest, const double* DTDv,
231  double* dTdX, double* volume, const index_t* elementId);
232 
233 void Assemble_jacobians_3D_M2D_E2D_C(const double* coordinates, int numQuad,
234  const double* QuadWeights, int numShape,
235  dim_t numElements, int numNodes, const index_t* nodes,
236  const double* DSDv, int numTest, const double* DTDv,
237  double* dTdX, double* volume, const index_t* elementId);
238 
239 void Assemble_jacobians_3D_M2D_E3D(const double* coordinates, int numQuad,
240  const double* QuadWeights, int numShape,
241  dim_t numElements, int numNodes, const index_t* nodes,
242  const double* DSDv, int numTest, const double* DTDv,
243  double* dTdX, double* volume, const index_t* elementId);
244 
245 void Assemble_jacobians_3D_M2D_E3D_C(const double* coordinates, int numQuad,
246  const double* QuadWeights, int numShape,
247  dim_t numElements, int numNodes, const index_t* nodes,
248  const double* DSDv, int numTest, const double* DTDv,
249  double* dTdX, double* volume, const index_t* elementId);
250 
251 } // namespace finley
252 
253 #endif // __FINLEY_ASSEMBLE_H__
254 
finley::Assemble_jacobians_3D_M2D_E2D
void Assemble_jacobians_3D_M2D_E2D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:676
finley::Assemble_NodeCoordinates
void Assemble_NodeCoordinates(const NodeFile *nodes, escript::Data &x)
copies node coordinates into expanded Data object x
Definition: finley/src/Assemble_NodeCoordinates.cpp:46
finley::AssembleParameters::NN
int NN
leading dimension of element node table
Definition: finley/src/Assemble.h:86
finley::AssembleParameters::numDim
int numDim
number of spatial dimensions
Definition: finley/src/Assemble.h:84
finley::Assemble_PDE_Points
void Assemble_PDE_Points(const AssembleParameters &p, const escript::Data &d_dirac, const escript::Data &y_dirac)
Definition: finley/src/Assemble_PDE_Points.cpp:68
finley::Assemble_integrate
void Assemble_integrate(const NodeFile *nodes, const ElementFile *elements, const escript::Data &data, Scalar *integrals)
integrates data on quadrature points
Definition: finley/src/Assemble_integrate.cpp:45
finley::AssembleParameters::F
escript::Data & F
right-hand side to be updated
Definition: finley/src/Assemble.h:74
finley::Assemble_AverageElementData
void Assemble_AverageElementData(const ElementFile *elements, escript::Data &out, const escript::Data &in)
averages data
Definition: finley/src/Assemble_AverageElementData.cpp:45
finley::AssembleParameters::col_numShapes
int col_numShapes
Definition: finley/src/Assemble.h:103
NodeFile.h
finley::ElementFile_Jacobians
Definition: finley/src/ElementFile.h:38
finley::AssembleParameters::numElements
dim_t numElements
number of elements
Definition: finley/src/Assemble.h:88
finley::Assemble_jacobians_2D_M1D_E1D
void Assemble_jacobians_2D_M1D_E1D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:150
finley::AssembleParameters::col_jac
ElementFile_Jacobians * col_jac
Definition: finley/src/Assemble.h:100
finley::Assemble_jacobians_1D
void Assemble_jacobians_1D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:60
finley::AssembleParameters::row_DOF_UpperBound
index_t row_DOF_UpperBound
Definition: finley/src/Assemble.h:92
finley::AssembleParameters::row_numShapes
int row_numShapes
Definition: finley/src/Assemble.h:96
finley::Assemble_PDE_System_2D
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: finley/src/Assemble_PDE_System_2D.cpp:78
finley::Assemble_jacobians_3D_M2D_E2D_C
void Assemble_jacobians_3D_M2D_E2D_C(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:743
finley::Assemble_interpolate
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: finley/src/Assemble_interpolate.cpp:47
finley::Assemble_PDE_Single_2D
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: finley/src/Assemble_PDE_Single_2D.cpp:72
finley::AssembleParameters::row_numShapesTotal
int row_numShapesTotal
Definition: finley/src/Assemble.h:95
finley::AssembleParameters::S
escript::ASM_ptr S
system matrix to be updated
Definition: finley/src/Assemble.h:72
finley::AssembleParameters::col_DOF
const index_t * col_DOF
Definition: finley/src/Assemble.h:98
finley::Assemble_jacobians_3D_M2D_E3D_C
void Assemble_jacobians_3D_M2D_E3D_C(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:543
finley::AssembleParameters::col_DOF_UpperBound
index_t col_DOF_UpperBound
Definition: finley/src/Assemble.h:99
finley::Assemble_jacobians_3D
void Assemble_jacobians_3D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:386
escript::Data
Data represents a collection of datapoints.
Definition: Data.h:62
finley::AssembleParameters::numSides
int numSides
number of sides
Definition: finley/src/Assemble.h:80
finley::Assemble_CopyElementData
void Assemble_CopyElementData(const ElementFile *elements, escript::Data &out, const escript::Data &in)
copies data between different types of elements
Definition: finley/src/Assemble_CopyElementData.cpp:43
escript::DataTypes::dim_t
index_t dim_t
Definition: DataTypes.h:87
finley::AssembleParameters::col_node
const int * col_node
Definition: finley/src/Assemble.h:101
finley::AssembleParameters::numComp
int numComp
Definition: finley/src/Assemble.h:97
ElementFile.h
finley::Assemble_jacobians_3D_M2D_E3D
void Assemble_jacobians_3D_M2D_E3D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:462
finley::Assemble_addToSystemMatrix
void Assemble_addToSystemMatrix(escript::ASM_ptr S, int NN_Equa, const index_t *Nodes_Equa, int num_Equa, int NN_Sol, const index_t *Nodes_Sol, int num_Sol, const Scalar *array)
finley::AssembleParameters::numQuadTotal
int numQuadTotal
total number of quadrature nodes = numQuadSub * numQuadSub
Definition: finley/src/Assemble.h:76
finley::Assemble_jacobians_2D
void Assemble_jacobians_2D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:97
finley::Assemble_getSize
void Assemble_getSize(const NodeFile *nodes, const ElementFile *elements, escript::Data &size)
Definition: finley/src/Assemble_getSize.cpp:46
finley::AssembleParameters
Definition: finley/src/Assemble.h:48
finley::Assemble_PDE_System_C
void Assemble_PDE_System_C(const AssembleParameters &p, const escript::Data &D, const escript::Data &Y)
Definition: Assemble_PDE_System_C.cpp:66
finley::AssembleParameters::numSub
int numSub
number of sub-elements
Definition: finley/src/Assemble.h:82
finley::Assemble_gradient
void Assemble_gradient(const NodeFile *nodes, const ElementFile *elements, escript::Data &gradient, const escript::Data &data)
Definition: finley/src/Assemble_gradient.cpp:47
finley::AssembleParameters::elements
const ElementFile * elements
element file these parameters apply to
Definition: finley/src/Assemble.h:70
finley::AssembleParameters::row_jac
ElementFile_Jacobians * row_jac
Definition: finley/src/Assemble.h:93
finley::Assemble_LumpedSystem
void Assemble_LumpedSystem(const NodeFile *nodes, const ElementFile *elements, escript::Data &lumpedMat, const escript::Data &D, bool useHRZ)
Definition: finley/src/Assemble_LumpedSystem.cpp:52
escript::DataTypes::index_t
int index_t
type for array/matrix indices used both globally and on each rank
Definition: DataTypes.h:82
finley::NodeFile
Definition: finley/src/NodeFile.h:39
finley::Assemble_PDE_System_3D
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: finley/src/Assemble_PDE_System_3D.cpp:78
finley::Assemble_jacobians_2D_M1D_E2D
void Assemble_jacobians_2D_M1D_E2D(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:256
finley::Assemble_PDE
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: finley/src/Assemble_PDE.cpp:122
finley::AssembleParameters::row_node
const int * row_node
Definition: finley/src/Assemble.h:94
finley::AssembleParameters::AssembleParameters
AssembleParameters(const NodeFile *nodes, const ElementFile *ef, escript::ASM_ptr sm, escript::Data &rhs, bool reducedOrder)
Definition: finley/src/Assemble_getAssembleParameters.cpp:31
Finley.h
finley::ElementFile
Definition: finley/src/ElementFile.h:72
finley::Assemble_jacobians_2D_M1D_E1D_C
void Assemble_jacobians_2D_M1D_E1D_C(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:196
finley::AssembleParameters::row_DOF
const index_t * row_DOF
Definition: finley/src/Assemble.h:91
finley::Assemble_getNormal
void Assemble_getNormal(const NodeFile *nodes, const ElementFile *elements, escript::Data &normals)
calculates the normal vector at quadrature points on face elements
Definition: finley/src/Assemble_getNormal.cpp:46
finley::Assemble_PDE_Single_C
void Assemble_PDE_Single_C(const AssembleParameters &p, const escript::Data &D, const escript::Data &Y)
Definition: Assemble_PDE_Single_C.cpp:64
S
#define S(_J_, _I_)
Definition: ShapeFunctions.cpp:134
finley::Assemble_jacobians_2D_M1D_E2D_C
void Assemble_jacobians_2D_M1D_E2D_C(const double *coordinates, int numQuad, const double *QuadWeights, int numShape, dim_t numElements, int numNodes, const index_t *nodes, const double *DSDv, int numTest, const double *DTDv, double *dTdX, double *volume, const index_t *elementId)
Definition: finley/src/Assemble_jacobians.cpp:310
escript::ASM_ptr
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:43
finley::AssembleParameters::numQuadSub
int numQuadSub
number of quadrature nodes per subelements
Definition: finley/src/Assemble.h:78
finley::Assemble_CopyNodalData
void Assemble_CopyNodalData(const NodeFile *nodes, escript::Data &out, const escript::Data &in)
copies data between different types of nodal representations
Definition: finley/src/Assemble_CopyNodalData.cpp:37
escript::Scalar
Data Scalar(double value, const FunctionSpace &what, bool expanded)
A collection of factory functions for creating Data objects which contain data points of various shap...
Definition: DataFactory.cpp:60
finley::AssembleParameters::numEqu
int numEqu
Definition: finley/src/Assemble.h:90
finley
A suite of factory methods for creating various finley domains.
Definition: finley/src/Assemble.h:31
finley::Assemble_PDE_System_1D
void Assemble_PDE_System_1D(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: Assemble_PDE_System_1D.cpp:78
finley::AssembleParameters::col_numShapesTotal
int col_numShapesTotal
Definition: finley/src/Assemble.h:102
finley::Assemble_PDE_Single_3D
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: finley/src/Assemble_PDE_Single_3D.cpp:72
finley::Assemble_PDE_Single_1D
void Assemble_PDE_Single_1D(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: Assemble_PDE_Single_1D.cpp:72