/[escript]/trunk/finley/src/Assemble.h
ViewVC logotype

Contents of /trunk/finley/src/Assemble.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6651 - (show annotations)
Wed Feb 7 02:12:08 2018 UTC (20 months ago) by jfenwick
File MIME type: text/plain
File size: 11705 byte(s)
Make everyone sad by touching all the files

Copyright dates update

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,
36 escript::ASM_ptr sm, escript::Data& rhs,
37 bool reducedOrder);
38
39 /// element file these parameters apply to
40 const ElementFile* elements;
41 /// system matrix to be updated
42 escript::ASM_ptr S;
43 /// right-hand side to be updated
44 escript::Data& F;
45 /// total number of quadrature nodes = numQuadSub * numQuadSub
46 int numQuadTotal;
47 /// number of quadrature nodes per subelements
48 int numQuadSub;
49 /// number of sides
50 int numSides;
51 /// number of sub-elements
52 int numSub;
53 /// number of spatial dimensions
54 int numDim;
55 /// leading dimension of element node table
56 int NN;
57 /// number of elements
58 dim_t numElements;
59
60 int numEqu;
61 const index_t* row_DOF;
62 index_t row_DOF_UpperBound;
63 ElementFile_Jacobians* row_jac;
64 const int* row_node;
65 int row_numShapesTotal;
66 int row_numShapes;
67 int numComp;
68 const index_t* col_DOF;
69 index_t col_DOF_UpperBound;
70 ElementFile_Jacobians* col_jac;
71 const int* col_node;
72 int col_numShapesTotal;
73 int col_numShapes;
74 };
75
76
77 /// Entry point for PDE assembly. Checks arguments, populates an
78 /// AssembleParameters structure and calls appropriate method for the actual
79 /// work.
80 void Assemble_PDE(const NodeFile* nodes, const ElementFile* elements,
81 escript::ASM_ptr S, escript::Data& F,
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>
87 void Assemble_PDE_Points(const AssembleParameters& p,
88 const escript::Data& d_dirac,
89 const escript::Data& y_dirac);
90
91 void Assemble_PDE_Single_1D(const AssembleParameters& p,
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>
97 void Assemble_PDE_Single_2D(const AssembleParameters& p,
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>
103 void Assemble_PDE_Single_3D(const AssembleParameters& p,
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>
109 void Assemble_PDE_Single_C(const AssembleParameters& p, const escript::Data& D,
110 const escript::Data& Y);
111
112 void Assemble_PDE_System_1D(const AssembleParameters& p,
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>
118 void Assemble_PDE_System_2D(const AssembleParameters& p,
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>
124 void Assemble_PDE_System_3D(const AssembleParameters& p,
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>
130 void Assemble_PDE_System_C(const AssembleParameters& p, const escript::Data& D,
131 const escript::Data& Y);
132
133 template<typename Scalar = double>
134 void Assemble_addToSystemMatrix(escript::ASM_ptr S, int NN_Equa,
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
142 /// averages data
143 template<typename Scalar>
144 void Assemble_AverageElementData(const ElementFile* elements,
145 escript::Data& out, const escript::Data& in);
146
147 /// copies data between different types of elements
148 template<typename Scalar>
149 void Assemble_CopyElementData(const ElementFile* elements, escript::Data& out,
150 const escript::Data& in);
151
152 /// copies data between different types of nodal representations
153 template<typename Scalar>
154 void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
155 const escript::Data& in);
156
157 /// copies node coordinates into expanded Data object `x`
158 void Assemble_NodeCoordinates(const NodeFile* nodes, escript::Data& x);
159
160 /// calculates the normal vector at quadrature points on face elements
161 void Assemble_getNormal(const NodeFile* nodes, const ElementFile* elements,
162 escript::Data& normals);
163
164 /// calculates the minimum distance between two vertices of elements and
165 /// assigns the value to each quadrature point in `size`
166 void Assemble_getSize(const NodeFile* nodes, const ElementFile* elements,
167 escript::Data& size);
168
169 /// Assemblage of Jacobians: calculates the gradient of nodal data at
170 /// quadrature points
171 template<typename Scalar>
172 void Assemble_gradient(const NodeFile* nodes, const ElementFile* elements,
173 escript::Data& gradient, const escript::Data& data);
174
175 /// integrates data on quadrature points
176 template<typename Scalar>
177 void Assemble_integrate(const NodeFile* nodes, const ElementFile* elements,
178 const escript::Data& data, Scalar* integrals);
179
180 /// interpolates nodal data in a data array onto elements (=integration points)
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

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26