/[escript]/branches/trilinos_from_5897/finley/src/Assemble.h
ViewVC logotype

Contents of /branches/trilinos_from_5897/finley/src/Assemble.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6058 - (show annotations)
Thu Mar 10 06:51:55 2016 UTC (3 years, 1 month ago) by caltinay
File MIME type: text/plain
File size: 10725 byte(s)
added getPtr() to AbstractSystemMatrix so we can now use shared systemmatrix
pointers rather than circumventing them.

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 Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
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
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
34 struct AssembleParameters {
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 void Assemble_PDE_Points(const AssembleParameters& p, const escript::Data& d_dirac,
87 const escript::Data& y_dirac);
88
89 void Assemble_PDE_Single_1D(const AssembleParameters& p,
90 const escript::Data& A, const escript::Data& B,
91 const escript::Data& C, const escript::Data& D,
92 const escript::Data& X, const escript::Data& Y);
93
94 void Assemble_PDE_Single_2D(const AssembleParameters& p,
95 const escript::Data& A, const escript::Data& B,
96 const escript::Data& C, const escript::Data& D,
97 const escript::Data& X, const escript::Data& Y);
98
99 void Assemble_PDE_Single_3D(const AssembleParameters& p,
100 const escript::Data& A, const escript::Data& B,
101 const escript::Data& C, const escript::Data& D,
102 const escript::Data& X, const escript::Data& Y);
103
104 void Assemble_PDE_Single_C(const AssembleParameters& p, const escript::Data& D,
105 const escript::Data& Y);
106
107 void Assemble_PDE_System_1D(const AssembleParameters& p,
108 const escript::Data& A, const escript::Data& B,
109 const escript::Data& C, const escript::Data& D,
110 const escript::Data& X, const escript::Data& Y);
111
112 void Assemble_PDE_System_2D(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 void Assemble_PDE_System_3D(const AssembleParameters& p,
118 const escript::Data& A, const escript::Data& B,
119 const escript::Data& C, const escript::Data& D,
120 const escript::Data& X, const escript::Data& Y);
121
122 void Assemble_PDE_System_C(const AssembleParameters& p, const escript::Data& D,
123 const escript::Data& Y);
124
125 void Assemble_addToSystemMatrix(escript::ASM_ptr S, int NN_Equa,
126 const index_t* Nodes_Equa, int num_Equa, int NN_Sol,
127 const index_t* Nodes_Sol, int num_Sol, const double* array);
128
129 void Assemble_LumpedSystem(const NodeFile* nodes, const ElementFile* elements,
130 escript::Data& lumpedMat, const escript::Data& D,
131 bool useHRZ);
132
133 void Assemble_AverageElementData(const ElementFile* elements,
134 escript::Data& out, const escript::Data& in);
135
136 void Assemble_CopyElementData(const ElementFile* elements, escript::Data& out,
137 const escript::Data& in);
138
139 void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
140 const escript::Data& in);
141
142 void Assemble_NodeCoordinates(const NodeFile* nodes, escript::Data& out);
143
144 void Assemble_getNormal(const NodeFile* nodes, const ElementFile* elements,
145 escript::Data& normals);
146
147 void Assemble_getSize(const NodeFile* nodes, const ElementFile* elements,
148 escript::Data& size);
149
150 void Assemble_gradient(const NodeFile* nodes, const ElementFile* elements,
151 escript::Data& gradient, const escript::Data& data);
152
153 void Assemble_integrate(const NodeFile* nodes, const ElementFile* elements,
154 const escript::Data& data, double* integrals);
155
156 void Assemble_interpolate(const NodeFile* nodes, const ElementFile* elements,
157 const escript::Data& data, escript::Data& output);
158
159 void Assemble_jacobians_1D(const double* coordinates, int numQuad,
160 const double* QuadWeights, int numShape,
161 dim_t numElements, dim_t numNodes, const index_t* nodes,
162 const double* DSDv, int numTest, const double* DTDv,
163 double* dTdX, double* volume, const index_t* elementId);
164 void Assemble_jacobians_2D(const double* coordinates, int numQuad,
165 const double* QuadWeights, int numShape,
166 dim_t numElements, dim_t numNodes, const index_t* nodes,
167 const double* DSDv, int numTest, const double* DTDv,
168 double* dTdX, double* volume, const index_t* elementId);
169 void Assemble_jacobians_2D_M1D_E1D(const double* coordinates, int numQuad,
170 const double* QuadWeights, int numShape,
171 dim_t numElements, dim_t numNodes, const index_t* nodes,
172 const double* DSDv, int numTest, const double* DTDv,
173 double* dTdX, double* volume, const index_t* elementId);
174 void Assemble_jacobians_2D_M1D_E1D_C(const double* coordinates, int numQuad,
175 const double* QuadWeights, int numShape,
176 dim_t numElements, dim_t numNodes, const index_t* nodes,
177 const double* DSDv, int numTest, const double* DTDv,
178 double* dTdX, double* volume, const index_t* elementId);
179 void Assemble_jacobians_2D_M1D_E2D(const double* coordinates, int numQuad,
180 const double* QuadWeights, int numShape,
181 dim_t numElements, dim_t numNodes, const index_t* nodes,
182 const double* DSDv, int numTest, const double* DTDv,
183 double* dTdX, double* volume, const index_t* elementId);
184 void Assemble_jacobians_2D_M1D_E2D_C(const double* coordinates, int numQuad,
185 const double* QuadWeights, int numShape,
186 dim_t numElements, dim_t numNodes, const index_t* nodes,
187 const double* DSDv, int numTest, const double* DTDv,
188 double* dTdX, double* volume, const index_t* elementId);
189 void Assemble_jacobians_3D(const double* coordinates, int numQuad,
190 const double* QuadWeights, int numShape,
191 dim_t numElements, dim_t numNodes, const index_t* nodes,
192 const double* DSDv, int numTest, const double* DTDv,
193 double* dTdX, double* volume, const index_t* elementId);
194 void Assemble_jacobians_3D_M2D_E2D(const double* coordinates, int numQuad,
195 const double* QuadWeights, int numShape,
196 dim_t numElements, dim_t numNodes, const index_t* nodes,
197 const double* DSDv, int numTest, const double* DTDv,
198 double* dTdX, double* volume, const index_t* elementId);
199 void Assemble_jacobians_3D_M2D_E2D_C(const double* coordinates, int numQuad,
200 const double* QuadWeights, int numShape,
201 dim_t numElements, dim_t numNodes, const index_t* nodes,
202 const double* DSDv, int numTest, const double* DTDv,
203 double* dTdX, double* volume, const index_t* elementId);
204 void Assemble_jacobians_3D_M2D_E3D(const double* coordinates, int numQuad,
205 const double* QuadWeights, int numShape,
206 dim_t numElements, dim_t numNodes, const index_t* nodes,
207 const double* DSDv, int numTest, const double* DTDv,
208 double* dTdX, double* volume, const index_t* elementId);
209 void Assemble_jacobians_3D_M2D_E3D_C(const double* coordinates, int numQuad,
210 const double* QuadWeights, int numShape,
211 dim_t numElements, dim_t 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 } // namespace finley
216
217 #endif // __FINLEY_ASSEMBLE_H__
218

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26