/[escript]/branches/trilinos_from_5897/dudley/src/ElementFile.cpp
ViewVC logotype

Contents of /branches/trilinos_from_5897/dudley/src/ElementFile.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6079 - (show annotations)
Mon Mar 21 12:22:38 2016 UTC (2 years, 10 months ago) by caltinay
File size: 6516 byte(s)
Big commit - making dudley much more like finley to make it more
managable. Fixed quite a few issues that had been fixed in finley.
Disposed of all ReducedNode/ReducedDOF entities that dudley never supported.
Compiles and passes tests.

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 #include "ElementFile.h"
18 #include "ShapeTable.h"
19
20 namespace dudley {
21
22 ElementFile::ElementFile(ElementTypeId type, escript::JMPI mpiInfo) :
23 MPIInfo(mpiInfo),
24 numElements(0),
25 Id(NULL),
26 Tag(NULL),
27 Owner(NULL),
28 Nodes(NULL),
29 Color(NULL),
30 minColor(0),
31 maxColor(-1),
32 etype(type)
33 {
34 jacobians = new ElementFile_Jacobians();
35 jacobians_reducedQ = new ElementFile_Jacobians();
36
37 numDim = Dims[type];
38 numNodes = numDim + 1;
39 numLocalDim = localDims[type];
40 numShapes = numLocalDim + 1;
41 ename = getElementName(type);
42 }
43
44 ElementFile::~ElementFile()
45 {
46 freeTable();
47 delete jacobians;
48 delete jacobians_reducedQ;
49 }
50
51 void ElementFile::allocTable(dim_t NE)
52 {
53 if (numElements > 0)
54 freeTable();
55
56 numElements = NE;
57 Owner = new int[numElements];
58 Id = new index_t[numElements];
59 Nodes = new index_t[numElements * numNodes];
60 Tag = new int[numElements];
61 Color = new index_t[numElements];
62
63 // this initialization makes sure that data are located on the right
64 // processor
65 #pragma omp parallel for
66 for (index_t e = 0; e < numElements; e++) {
67 for (int i = 0; i < numNodes; i++)
68 Nodes[INDEX2(i, e, numNodes)] = -1;
69 Owner[e] = -1;
70 Id[e] = -1;
71 Tag[e] = -1;
72 Color[e] = -1;
73 }
74 maxColor = -1;
75 minColor = 0;
76 }
77
78 void ElementFile::freeTable()
79 {
80 delete[] Owner;
81 delete[] Id;
82 delete[] Nodes;
83 delete[] Tag;
84 delete[] Color;
85 tagsInUse.clear();
86 numElements = 0;
87 maxColor = -1;
88 minColor = 0;
89 }
90
91 void ElementFile::copyTable(index_t offset, index_t nodeOffset,
92 index_t idOffset, const ElementFile* in)
93 {
94 const int NN_in = in->numNodes;
95 if (NN_in > numNodes) {
96 throw DudleyException("ElementFile::copyTable: dimensions of element files don't match.");
97 }
98
99 if (MPIInfo->comm != in->MPIInfo->comm) {
100 throw DudleyException("ElementFile::copyTable: MPI communicators of element files don't match.");
101 }
102
103 #pragma omp parallel for
104 for (index_t n = 0; n < in->numElements; n++) {
105 Owner[offset + n] = in->Owner[n];
106 Id[offset + n] = in->Id[n] + idOffset;
107 Tag[offset + n] = in->Tag[n];
108 for (int i = 0; i < numNodes; i++)
109 Nodes[INDEX2(i, offset + n, numNodes)] =
110 in->Nodes[INDEX2(i, n, NN_in)] + nodeOffset;
111 }
112 }
113
114 void ElementFile::gather(const index_t* index, const ElementFile* in)
115 {
116 const int NN_in = in->numNodes;
117 #pragma omp parallel for
118 for (index_t e = 0; e < numElements; e++) {
119 const index_t k = index[e];
120 Id[e] = in->Id[k];
121 Tag[e] = in->Tag[k];
122 Owner[e] = in->Owner[k];
123 Color[e] = in->Color[k] + maxColor + 1;
124 for (int j = 0; j < std::min(numNodes, NN_in); j++)
125 Nodes[INDEX2(j, e, numNodes)] = in->Nodes[INDEX2(j, k, NN_in)];
126 }
127 minColor = std::min(minColor, in->minColor + maxColor + 1);
128 maxColor = std::max(maxColor, in->maxColor + maxColor + 1);
129 }
130
131 void ElementFile::swapTable(ElementFile* other)
132 {
133 std::swap(numElements, other->numElements);
134 std::swap(Owner, other->Owner);
135 std::swap(Id, other->Id);
136 std::swap(Nodes, other->Nodes);
137 std::swap(Tag, other->Tag);
138 std::swap(Color, other->Color);
139 std::swap(minColor, other->minColor);
140 std::swap(maxColor, other->maxColor);
141 std::swap(tagsInUse, other->tagsInUse);
142 }
143
144 void ElementFile::optimizeOrdering()
145 {
146 if (numElements < 1)
147 return;
148
149 util::ValueAndIndexList item_list(numElements);
150 index_t* index = new index_t[numElements];
151 ElementFile* out = new ElementFile(etype, MPIInfo);
152 out->allocTable(numElements);
153
154 #pragma omp parallel for
155 for (index_t e = 0; e < numElements; e++) {
156 std::pair<index_t,index_t> entry(Nodes[INDEX2(0, e, numNodes)], e);
157 for (int i = 1; i < numNodes; i++)
158 entry.first = std::min(entry.first, Nodes[INDEX2(i, e, numNodes)]);
159 item_list[e] = entry;
160 }
161 util::sortValueAndIndex(item_list);
162
163 #pragma omp parallel for
164 for (index_t e = 0; e < numElements; e++)
165 index[e] = item_list[e].second;
166
167 out->gather(index, this);
168 swapTable(out);
169 delete out;
170 delete[] index;
171 }
172
173 void ElementFile::setTags(int newTag, const escript::Data& mask)
174 {
175 const int numQuad = hasReducedIntegrationOrder(mask) ? 1 : numNodes;
176
177 if (1 != mask.getDataPointSize()) {
178 throw DudleyException("ElementFile::setTags: number of components of mask must be 1.");
179 } else if (!mask.numSamplesEqual(numQuad, numElements)) {
180 throw DudleyException("ElementFile::setTags: illegal number of samples of mask Data object");
181 }
182
183 if (mask.actsExpanded()) {
184 #pragma omp parallel for
185 for (index_t n = 0; n < numElements; n++) {
186 if (mask.getSampleDataRO(n)[0] > 0)
187 Tag[n] = newTag;
188 }
189 } else {
190 #pragma omp parallel for
191 for (index_t n = 0; n < numElements; n++) {
192 const double* mask_array = mask.getSampleDataRO(n);
193 bool check = false;
194 for (int q = 0; q < numQuad; q++)
195 check = check || mask_array[q];
196 if (check)
197 Tag[n] = newTag;
198 }
199 }
200 updateTagList();
201 }
202
203 void ElementFile::markNodes(std::vector<short>& mask, index_t offset) const
204 {
205 #pragma omp parallel for
206 for (index_t e = 0; e < numElements; e++) {
207 for (int i = 0; i < numNodes; i++) {
208 mask[Nodes[INDEX2(i, e, numNodes)] - offset] = 1;
209 }
210 }
211 }
212
213 void ElementFile::relabelNodes(const index_t* newNode, index_t offset)
214 {
215 #pragma omp parallel for
216 for (index_t j = 0; j < numElements; j++) {
217 for (int i = 0; i < numNodes; i++) {
218 Nodes[INDEX2(i, j, numNodes)] =
219 newNode[Nodes[INDEX2(i, j, numNodes)] - offset];
220 }
221 }
222 }
223
224 } // namespace dudley
225

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision
svn:mergeinfo /branches/4.0fordebian/dudley/src/ElementFile.cpp:5567-5588 /branches/lapack2681/finley/src/ElementFile.cpp:2682-2741 /branches/pasowrap/dudley/src/ElementFile.cpp:3661-3674 /branches/py3_attempt2/dudley/src/ElementFile.cpp:3871-3891 /branches/restext/finley/src/ElementFile.cpp:2610-2624 /branches/ripleygmg_from_3668/dudley/src/ElementFile.cpp:3669-3791 /branches/stage3.0/finley/src/ElementFile.cpp:2569-2590 /branches/symbolic_from_3470/dudley/src/ElementFile.cpp:3471-3974 /branches/symbolic_from_3470/ripley/test/python/dudley/src/ElementFile.cpp:3517-3974 /release/3.0/finley/src/ElementFile.cpp:2591-2601 /release/4.0/dudley/src/ElementFile.cpp:5380-5406 /trunk/dudley/src/ElementFile.cpp:4257-4344,5898-6007 /trunk/ripley/test/python/dudley/src/ElementFile.cpp:3480-3515

  ViewVC Help
Powered by ViewVC 1.1.26