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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6009 - (show annotations)
Wed Mar 2 04:13:26 2016 UTC (3 years, 1 month ago) by caltinay
File MIME type: text/plain
File size: 6051 byte(s)
Much needed sync with trunk...

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 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
33 struct Assemble_Parameters {
34 // number of quadrature nodes
35 int numQuad;
36 // number of spatial dimensions
37 int numDim;
38 // leading dimension of element node table
39 int NN;
40 // number of elements
41 dim_t numElements;
42
43 int numEqu;
44 const index_t* row_DOF;
45 dim_t row_DOF_UpperBound;
46 Dudley_ElementFile_Jacobians* row_jac;
47 int numShapes;
48
49 int numComp;
50 const index_t* col_DOF;
51 dim_t col_DOF_UpperBound;
52
53 const double* shapeFns;
54 };
55
56 inline bool Assemble_reducedIntegrationOrder(const escript::Data* in)
57 {
58 const int fs = in->getFunctionSpace().getTypeCode();
59 return (fs == DUDLEY_REDUCED_ELEMENTS || fs == DUDLEY_REDUCED_FACE_ELEMENTS);
60 }
61
62 void Assemble_getAssembleParameters(const Dudley_NodeFile*, const Dudley_ElementFile*,
63 escript::ASM_ptr, const escript::Data&,
64 bool, Assemble_Parameters*);
65
66 void Assemble_PDE(const Dudley_NodeFile* nodes, const Dudley_ElementFile* elements,
67 escript::ASM_ptr S, escript::Data& F,
68 const escript::Data& A, const escript::Data& B,
69 const escript::Data& C, const escript::Data& D,
70 const escript::Data& X, const escript::Data& Y);
71
72 void Assemble_PDE_Points(const Assemble_Parameters& p, const Dudley_ElementFile*,
73 escript::ASM_ptr S, escript::Data& F,
74 const escript::Data& d_dirac,
75 const escript::Data& y_dirac);
76
77 void Assemble_PDE_Single_2D(const Assemble_Parameters& p, const Dudley_ElementFile*,
78 escript::ASM_ptr S, escript::Data& F,
79 const escript::Data& A, const escript::Data& B,
80 const escript::Data& C, const escript::Data& D,
81 const escript::Data& X, const escript::Data& Y);
82
83 void Assemble_PDE_Single_3D(const Assemble_Parameters& p, const Dudley_ElementFile*,
84 escript::ASM_ptr S, escript::Data& F,
85 const escript::Data& A, const escript::Data& B,
86 const escript::Data& C, const escript::Data& D,
87 const escript::Data& X, const escript::Data& Y);
88
89 void Assemble_PDE_System_2D(const Assemble_Parameters& p, const Dudley_ElementFile*,
90 escript::ASM_ptr S, escript::Data& F,
91 const escript::Data& A, const escript::Data& B,
92 const escript::Data& C, const escript::Data& D,
93 const escript::Data& X, const escript::Data& Y);
94
95 void Assemble_PDE_System_3D(const Assemble_Parameters& p, const Dudley_ElementFile*,
96 escript::ASM_ptr S, escript::Data& F,
97 const escript::Data& A, const escript::Data& B,
98 const escript::Data& C, const escript::Data& D,
99 const escript::Data& X, const escript::Data& Y);
100
101 void Assemble_NodeCoordinates(Dudley_NodeFile *, escript::Data *);
102 void Assemble_setNormal(Dudley_NodeFile *, Dudley_ElementFile *, escript::Data *);
103 void Assemble_interpolate(Dudley_NodeFile *, Dudley_ElementFile *, const escript::Data *, escript::Data *);
104 void Assemble_gradient(Dudley_NodeFile *, Dudley_ElementFile *, escript::Data *, const escript::Data *);
105 void Assemble_integrate(Dudley_NodeFile *, Dudley_ElementFile *, const escript::Data *, double *);
106 void Assemble_getSize(Dudley_NodeFile *, Dudley_ElementFile *, escript::Data *);
107 void Assemble_CopyNodalData(Dudley_NodeFile * nodes, escript::Data * out, const escript::Data * in);
108 void Assemble_CopyElementData(Dudley_ElementFile * elements, escript::Data * out, const escript::Data * in);
109 void Assemble_AverageElementData(Dudley_ElementFile * elements, escript::Data * out, const escript::Data * in);
110 void Assemble_addToSystemMatrix(escript::ASM_ptr in, const dim_t NN_Equa, const index_t * Nodes_Equa, const dim_t num_Equa,
111 const dim_t NN_Sol, const index_t * Nodes_Sol, const dim_t num_Sol, const double *array);
112
113 void Assemble_jacobians_2D(double *, dim_t, dim_t, dim_t, index_t *, double *, double *abs_D, double *quadweight,
114 index_t *);
115 void Assemble_jacobians_2D_M1D_E1D(double *, dim_t, dim_t, dim_t, index_t *, double *, double *abs_D,
116 double *quadweight, index_t *);
117 void Assemble_jacobians_3D(double *, dim_t, dim_t, dim_t, index_t *, double *, double *abs_D, double *quadweight,
118 index_t *);
119 void Assemble_jacobians_3D_M2D_E2D(double *, dim_t, dim_t, dim_t, index_t *, double *, double *abs_D,
120 double *quadweight, index_t *);
121
122 void Assemble_LumpedSystem(Dudley_NodeFile* nodes, Dudley_ElementFile* elements,
123 escript::Data& lumpedMat, const escript::Data& D,
124 bool useHRZ);
125
126 } // namespace dudley
127
128 #endif // __DUDLEY_ASSEMBLE_H__
129

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26