escript  Revision_
finley/src/Assemble.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2020 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-2017 by Centre for Geoscience Computing (GeoComp)
14 * Development from 2019 by School of Earth and Environmental Sciences
15 **
16 *****************************************************************************/
17 
18 /****************************************************************************
19 
20  Assemblage routines: header file
21 
22 *****************************************************************************/
23 
24 #ifndef __FINLEY_ASSEMBLE_H__
25 #define __FINLEY_ASSEMBLE_H__
26 
27 #include "Finley.h"
28 #include "ElementFile.h"
29 #include "NodeFile.h"
30 #include <escript/AbstractSystemMatrix.h>
31 
32 namespace finley {
33 
35 {
36  AssembleParameters(const NodeFile* nodes, const ElementFile* ef,
38  bool reducedOrder);
39 
51  int numSides;
53  int numSub;
55  int numDim;
57  int NN;
60 
61  int numEqu;
62  const index_t* row_DOF;
65  const int* row_node;
68  int numComp;
69  const index_t* col_DOF;
72  const int* col_node;
75 };
76 
77 
81 void Assemble_PDE(const NodeFile* nodes, const ElementFile* elements,
83  const escript::Data& A, const escript::Data& B,
84  const escript::Data& C, const escript::Data& D,
85  const escript::Data& X, const escript::Data& Y);
86 
87 template<typename Scalar>
89  const escript::Data& d_dirac,
90  const escript::Data& y_dirac);
91 
93  const escript::Data& A, const escript::Data& B,
94  const escript::Data& C, const escript::Data& D,
95  const escript::Data& X, const escript::Data& Y);
96 
97 template<typename Scalar>
99  const escript::Data& A, const escript::Data& B,
100  const escript::Data& C, const escript::Data& D,
101  const escript::Data& X, const escript::Data& Y);
102 
103 template<typename Scalar>
105  const escript::Data& A, const escript::Data& B,
106  const escript::Data& C, const escript::Data& D,
107  const escript::Data& X, const escript::Data& Y);
108 
109 template<typename Scalar>
111  const escript::Data& Y);
112 
114  const escript::Data& A, const escript::Data& B,
115  const escript::Data& C, const escript::Data& D,
116  const escript::Data& X, const escript::Data& Y);
117 
118 template<typename Scalar>
120  const escript::Data& A, const escript::Data& B,
121  const escript::Data& C, const escript::Data& D,
122  const escript::Data& X, const escript::Data& Y);
123 
124 template<typename Scalar>
126  const escript::Data& A, const escript::Data& B,
127  const escript::Data& C, const escript::Data& D,
128  const escript::Data& X, const escript::Data& Y);
129 
130 template<typename Scalar>
132  const escript::Data& Y);
133 
134 template<typename Scalar = double>
136  const index_t* Nodes_Equa, int num_Equa, int NN_Sol,
137  const index_t* Nodes_Sol, int num_Sol, const Scalar* array);
138 
139 void Assemble_LumpedSystem(const NodeFile* nodes, const ElementFile* elements,
140  escript::Data& lumpedMat, const escript::Data& D,
141  bool useHRZ);
142 
144 template<typename Scalar>
145 void Assemble_AverageElementData(const ElementFile* elements,
146  escript::Data& out, const escript::Data& in);
147 
149 template<typename Scalar>
150 void Assemble_CopyElementData(const ElementFile* elements, escript::Data& out,
151  const escript::Data& in);
152 
154 template<typename Scalar>
155 void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
156  const escript::Data& in);
157 
159 void Assemble_NodeCoordinates(const NodeFile* nodes, escript::Data& x);
160 
162 void Assemble_getNormal(const NodeFile* nodes, const ElementFile* elements,
163  escript::Data& normals);
164 
167 void Assemble_getSize(const NodeFile* nodes, const ElementFile* elements,
168  escript::Data& size);
169 
172 template<typename Scalar>
173 void Assemble_gradient(const NodeFile* nodes, const ElementFile* elements,
174  escript::Data& gradient, const escript::Data& data);
175 
177 template<typename Scalar>
178 void Assemble_integrate(const NodeFile* nodes, const ElementFile* elements,
179  const escript::Data& data, Scalar* integrals);
180 
182 template<typename Scalar>
183 void Assemble_interpolate(const NodeFile* nodes, const ElementFile* elements,
184  const escript::Data& data, escript::Data& output);
185 
186 void Assemble_jacobians_1D(const double* coordinates, int numQuad,
187  const double* QuadWeights, int numShape,
188  dim_t numElements, int numNodes, const index_t* nodes,
189  const double* DSDv, int numTest, const double* DTDv,
190  double* dTdX, double* volume, const index_t* elementId);
191 
192 void Assemble_jacobians_2D(const double* coordinates, int numQuad,
193  const double* QuadWeights, int numShape,
194  dim_t numElements, int numNodes, const index_t* nodes,
195  const double* DSDv, int numTest, const double* DTDv,
196  double* dTdX, double* volume, const index_t* elementId);
197 
198 void Assemble_jacobians_2D_M1D_E1D(const double* coordinates, int numQuad,
199  const double* QuadWeights, int numShape,
200  dim_t numElements, int numNodes, const index_t* nodes,
201  const double* DSDv, int numTest, const double* DTDv,
202  double* dTdX, double* volume, const index_t* elementId);
203 
204 void Assemble_jacobians_2D_M1D_E1D_C(const double* coordinates, int numQuad,
205  const double* QuadWeights, int numShape,
206  dim_t numElements, int numNodes, const index_t* nodes,
207  const double* DSDv, int numTest, const double* DTDv,
208  double* dTdX, double* volume, const index_t* elementId);
209 
210 void Assemble_jacobians_2D_M1D_E2D(const double* coordinates, int numQuad,
211  const double* QuadWeights, int numShape,
212  dim_t numElements, int numNodes, const index_t* nodes,
213  const double* DSDv, int numTest, const double* DTDv,
214  double* dTdX, double* volume, const index_t* elementId);
215 
216 void Assemble_jacobians_2D_M1D_E2D_C(const double* coordinates, int numQuad,
217  const double* QuadWeights, int numShape,
218  dim_t numElements, int numNodes, const index_t* nodes,
219  const double* DSDv, int numTest, const double* DTDv,
220  double* dTdX, double* volume, const index_t* elementId);
221 
222 void Assemble_jacobians_3D(const double* coordinates, int numQuad,
223  const double* QuadWeights, int numShape,
224  dim_t numElements, int numNodes, const index_t* nodes,
225  const double* DSDv, int numTest, const double* DTDv,
226  double* dTdX, double* volume, const index_t* elementId);
227 
228 void Assemble_jacobians_3D_M2D_E2D(const double* coordinates, int numQuad,
229  const double* QuadWeights, int numShape,
230  dim_t numElements, int numNodes, const index_t* nodes,
231  const double* DSDv, int numTest, const double* DTDv,
232  double* dTdX, double* volume, const index_t* elementId);
233 
234 void Assemble_jacobians_3D_M2D_E2D_C(const double* coordinates, int numQuad,
235  const double* QuadWeights, int numShape,
236  dim_t numElements, int numNodes, const index_t* nodes,
237  const double* DSDv, int numTest, const double* DTDv,
238  double* dTdX, double* volume, const index_t* elementId);
239 
240 void Assemble_jacobians_3D_M2D_E3D(const double* coordinates, int numQuad,
241  const double* QuadWeights, int numShape,
242  dim_t numElements, int numNodes, const index_t* nodes,
243  const double* DSDv, int numTest, const double* DTDv,
244  double* dTdX, double* volume, const index_t* elementId);
245 
246 void Assemble_jacobians_3D_M2D_E3D_C(const double* coordinates, int numQuad,
247  const double* QuadWeights, int numShape,
248  dim_t numElements, int numNodes, const index_t* nodes,
249  const double* DSDv, int numTest, const double* DTDv,
250  double* dTdX, double* volume, const index_t* elementId);
251 
252 } // namespace finley
253 
254 #endif // __FINLEY_ASSEMBLE_H__
255 
#define S(_J_, _I_)
Definition: ShapeFunctions.cpp:122
Data represents a collection of datapoints.
Definition: Data.h:64
Definition: finley/src/ElementFile.h:63
Definition: finley/src/NodeFile.h:42
index_t dim_t
Definition: DataTypes.h:66
int index_t
type for array/matrix indices used both globally and on each rank
Definition: DataTypes.h:61
boost::shared_ptr< AbstractSystemMatrix > ASM_ptr
Definition: AbstractSystemMatrix.h:33
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:49
A suite of factory methods for creating various finley domains.
Definition: finley/src/Assemble.h:32
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)
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:50
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:376
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:46
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:46
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:27
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:33
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:140
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:733
void Assemble_NodeCoordinates(const NodeFile *nodes, escript::Data &x)
copies node coordinates into expanded Data object x
Definition: finley/src/Assemble_NodeCoordinates.cpp:34
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:46
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:34
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:31
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:186
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:34
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:452
void Assemble_PDE_System_C(const AssembleParameters &p, const escript::Data &D, const escript::Data &Y)
Definition: Assemble_PDE_System_C.cpp:43
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:49
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:300
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:49
void Assemble_LumpedSystem(const NodeFile *nodes, const ElementFile *elements, escript::Data &lumpedMat, const escript::Data &D, bool useHRZ)
Definition: finley/src/Assemble_LumpedSystem.cpp:37
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:86
void Assemble_getSize(const NodeFile *nodes, const ElementFile *elements, escript::Data &size)
Definition: finley/src/Assemble_getSize.cpp:33
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:87
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:246
void Assemble_PDE_Points(const AssembleParameters &p, const escript::Data &d_dirac, const escript::Data &y_dirac)
Definition: finley/src/Assemble_PDE_Points.cpp:44
void Assemble_PDE_Single_C(const AssembleParameters &p, const escript::Data &D, const escript::Data &Y)
Definition: Assemble_PDE_Single_C.cpp:42
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:49
void Assemble_gradient(const NodeFile *nodes, const ElementFile *elements, escript::Data &gradient, const escript::Data &data)
Definition: finley/src/Assemble_gradient.cpp:34
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:533
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:666
void Assemble_AverageElementData(const ElementFile *elements, escript::Data &out, const escript::Data &in)
averages data
Definition: finley/src/Assemble_AverageElementData.cpp:33
Definition: finley/src/Assemble.h:35
int NN
leading dimension of element node table
Definition: finley/src/Assemble.h:57
ElementFile_Jacobians * col_jac
Definition: finley/src/Assemble.h:71
int row_numShapesTotal
Definition: finley/src/Assemble.h:66
int numComp
Definition: finley/src/Assemble.h:68
const ElementFile * elements
element file these parameters apply to
Definition: finley/src/Assemble.h:41
index_t col_DOF_UpperBound
Definition: finley/src/Assemble.h:70
escript::ASM_ptr S
system matrix to be updated
Definition: finley/src/Assemble.h:43
const index_t * row_DOF
Definition: finley/src/Assemble.h:62
int numQuadTotal
total number of quadrature nodes = numQuadSub * numQuadSub
Definition: finley/src/Assemble.h:47
int numSub
number of sub-elements
Definition: finley/src/Assemble.h:53
int row_numShapes
Definition: finley/src/Assemble.h:67
int col_numShapes
Definition: finley/src/Assemble.h:74
int numQuadSub
number of quadrature nodes per subelements
Definition: finley/src/Assemble.h:49
int col_numShapesTotal
Definition: finley/src/Assemble.h:73
ElementFile_Jacobians * row_jac
Definition: finley/src/Assemble.h:64
int numEqu
Definition: finley/src/Assemble.h:61
escript::Data & F
right-hand side to be updated
Definition: finley/src/Assemble.h:45
const int * col_node
Definition: finley/src/Assemble.h:72
index_t row_DOF_UpperBound
Definition: finley/src/Assemble.h:63
int numSides
number of sides
Definition: finley/src/Assemble.h:51
dim_t numElements
number of elements
Definition: finley/src/Assemble.h:59
int numDim
number of spatial dimensions
Definition: finley/src/Assemble.h:55
const index_t * col_DOF
Definition: finley/src/Assemble.h:69
AssembleParameters(const NodeFile *nodes, const ElementFile *ef, escript::ASM_ptr sm, escript::Data &rhs, bool reducedOrder)
Definition: finley/src/Assemble_getAssembleParameters.cpp:34
const int * row_node
Definition: finley/src/Assemble.h:65
Definition: finley/src/ElementFile.h:29