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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2748 - (show annotations)
Tue Nov 17 07:32:59 2009 UTC (10 years ago) by gross
File MIME type: text/plain
File size: 7097 byte(s)
Macro elements are implemented now. VTK writer for macro elements still needs testing.
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2009 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 #ifndef INC_FINLEY_ELEMENTFILE
16 #define INC_FINLEY_ELEMENTFILE
17
18 #include "Finley.h"
19 #include "NodeFile.h"
20 #include "ReferenceElementSets.h"
21 #include "escript/DataC.h"
22
23 #ifdef PASO_MPI
24 #include "paso/Paso_MPI.h"
25 #endif
26
27
28 struct Finley_ElementFile_Jacobeans {
29 Finley_Status_t status; /* status of mesh when jacobeans where updated last time */
30 dim_t numDim; /* spatial dimension */
31 Finley_ShapeFunction* BasisFunctions; /* basis function used */
32 dim_t numQuadTotal; /* total number of quadrature nodes used to calculate jacobeans = numSub * BasisFunctions->numQuadNodes*/
33 dim_t numSides; /* number of sides (=1 normal, =2 contact) */
34 index_t* offsets; /* offset to sides (borrowed reference) */
35 dim_t numSub; /* number of subelements */
36 dim_t numShapesTotal; /* total number of shape functions = BasisFunctions->numShapes * numSides */
37 index_t* node_selection; /* local node selection list of length numSub * numShapesTotal (borrowed reference) */
38 dim_t numElements; /* number of elements */
39 double* volume; /* local volume */
40 double* DSDX; /* derivatives of shape functions in global coordinates at quadrature points*/
41 };
42
43 typedef struct Finley_ElementFile_Jacobeans Finley_ElementFile_Jacobeans;
44
45 struct Finley_ElementFile {
46 Paso_MPIInfo *MPIInfo;
47 Paso_MPI_rank *Owner;
48
49 Finley_ReferenceElementSet *referenceElementSet; /* the reference element to be used */
50
51 dim_t numElements; /* number of elements. */
52
53 index_t *Id; /* Id[i] is the id nmber of
54 node i. this number is not
55 used but useful when
56 elements are resorted. in
57 the entire code the term
58 'element id' refers to i
59 but nor to Id[i] if not
60 explicitly stated
61 otherwise. */
62
63 index_t *Tag; /* Tag[i] is the tag of element i. */
64
65 index_t *tagsInUse; /* array of tags which are actually used */
66 dim_t numTagsInUse; /* number of tags used */
67
68
69 dim_t numNodes; /* number of nodes per element */
70 index_t *Nodes; /* Nodes[INDEX(k, i, numNodes)]
71 is the k-the node in the
72 i-the element. note that
73 in the way the nodes are
74 ordered Nodes[INDEX(k, i, numNodes)
75 is k-the node of element i
76 when refering to the
77 linear version of the
78 mesh. */
79 index_t minColor; /* minimum color */
80 index_t maxColor; /* maximum color */
81 index_t *Color; /* assigns each element a color. elements with the same color
82 are don't share a node so they can be processed simultaneously
83 at anytime Color must provide a valid value. In any case one can set
84 Color[e]=e for all e */
85 index_t order; /* order of the element integration scheme*/
86 index_t reduced_order; /* order of the reduced element integration scheme*/
87
88 Finley_ElementFile_Jacobeans* jacobeans; /* jacobeans of the shape function used for solution approximation */
89 Finley_ElementFile_Jacobeans* jacobeans_reducedS; /* jacobeans of the shape function used for solution approximation for reduced order of shape function*/
90 Finley_ElementFile_Jacobeans* jacobeans_reducedQ; /* jacobeans of the shape function used for solution approximation for reduced integration order*/
91 Finley_ElementFile_Jacobeans* jacobeans_reducedS_reducedQ; /* jacobeans of the shape function used for solution approximation for reduced integration order and reduced order of shape function*/
92
93 };
94
95 typedef struct Finley_ElementFile Finley_ElementFile;
96 Finley_ElementFile* Finley_ElementFile_alloc(Finley_ReferenceElementSet* referenceElementSet, Paso_MPIInfo *MPIInfo);
97 void Finley_ElementFile_free(Finley_ElementFile*);
98 void Finley_ElementFile_allocTable(Finley_ElementFile*,dim_t);
99 void Finley_ElementFile_freeTable(Finley_ElementFile*);
100 void Finley_ElementFile_setElementDistribution(Finley_ElementFile* in, dim_t* distribution);
101 dim_t Finley_ElementFile_getGlobalNumElements(Finley_ElementFile* in);
102 dim_t Finley_ElementFile_getMyNumElements(Finley_ElementFile* in);
103 index_t Finley_ElementFile_getFirstElement(Finley_ElementFile* in);
104 void Finley_ElementFile_distributeByRankOfDOF(Finley_ElementFile* self, Paso_MPI_rank* mpiRankOfDOF, index_t *Id);
105
106 void Finley_ElementFile_createColoring(Finley_ElementFile* in,dim_t numNodes,dim_t* degreeOfFreedom);
107 void Finley_ElementFile_optimizeOrdering(Finley_ElementFile** in);
108 void Finley_ElementFile_setNodeRange(dim_t*,dim_t*,Finley_ElementFile*);
109 void Finley_ElementFile_relableNodes(dim_t*,dim_t,Finley_ElementFile*);
110 void Finley_ElementFile_markNodes(dim_t*,dim_t,dim_t,Finley_ElementFile*,dim_t);
111 void Finley_ElementFile_scatter(dim_t*,Finley_ElementFile*,Finley_ElementFile*);
112 void Finley_ElementFile_gather(dim_t*,Finley_ElementFile*,Finley_ElementFile*);
113 void Finley_ElementFile_copyTable(dim_t,Finley_ElementFile*,dim_t,dim_t,Finley_ElementFile*);
114 void Finley_ElementFile_markDOFsConnectedToRange(index_t* mask,index_t offset,index_t marker,index_t firstDOF,index_t lastDOF,index_t *dofIndex,Finley_ElementFile*in ,bool_t useLinear);
115
116 void Finley_ElementFile_setTags(Finley_ElementFile* ,const int, escriptDataC*);
117 Finley_ElementFile_Jacobeans* Finley_ElementFile_Jacobeans_alloc(Finley_ShapeFunction* );
118 void Finley_ElementFile_Jacobeans_dealloc(Finley_ElementFile_Jacobeans*);
119 Finley_ElementFile_Jacobeans* Finley_ElementFile_borrowJacobeans(Finley_ElementFile*, Finley_NodeFile*, bool_t, bool_t);
120 void Finley_ElementFile_setTagsInUse(Finley_ElementFile* in);
121
122
123 #endif /* #ifndef INC_FINLEY_ELEMENTFILE */
124

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26