/[escript]/trunk/finley/src/CPPAdapter/MeshAdapter.cpp
ViewVC logotype

Contents of /trunk/finley/src/CPPAdapter/MeshAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1796 - (show annotations)
Wed Sep 17 01:45:46 2008 UTC (10 years, 7 months ago) by jfenwick
File size: 83416 byte(s)
Merged noarrayview branch onto trunk.


1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 #include "MeshAdapter.h"
17 #include "escript/Data.h"
18 #include "escript/DataFactory.h"
19 #ifdef USE_NETCDF
20 #include <netcdfcpp.h>
21 #endif
22 extern "C" {
23 #include "escript/blocktimer.h"
24 }
25 #include <vector>
26
27 #define IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH 256
28
29 using namespace std;
30 using namespace escript;
31
32 namespace finley {
33
34 //
35 // define the static constants
36 MeshAdapter::FunctionSpaceNamesMapType MeshAdapter::m_functionSpaceTypeNames;
37 const int MeshAdapter::DegreesOfFreedom=FINLEY_DEGREES_OF_FREEDOM;
38 const int MeshAdapter::ReducedDegreesOfFreedom=FINLEY_REDUCED_DEGREES_OF_FREEDOM;
39 const int MeshAdapter::Nodes=FINLEY_NODES;
40 const int MeshAdapter::ReducedNodes=FINLEY_REDUCED_NODES;
41 const int MeshAdapter::Elements=FINLEY_ELEMENTS;
42 const int MeshAdapter::ReducedElements=FINLEY_REDUCED_ELEMENTS;
43 const int MeshAdapter::FaceElements=FINLEY_FACE_ELEMENTS;
44 const int MeshAdapter::ReducedFaceElements=FINLEY_REDUCED_FACE_ELEMENTS;
45 const int MeshAdapter::Points=FINLEY_POINTS;
46 const int MeshAdapter::ContactElementsZero=FINLEY_CONTACT_ELEMENTS_1;
47 const int MeshAdapter::ReducedContactElementsZero=FINLEY_REDUCED_CONTACT_ELEMENTS_1;
48 const int MeshAdapter::ContactElementsOne=FINLEY_CONTACT_ELEMENTS_2;
49 const int MeshAdapter::ReducedContactElementsOne=FINLEY_REDUCED_CONTACT_ELEMENTS_2;
50
51 MeshAdapter::MeshAdapter(Finley_Mesh* finleyMesh)
52 {
53 setFunctionSpaceTypeNames();
54 //
55 // need to use a null_deleter as Finley_Mesh_free deletes the pointer
56 // for us.
57 m_finleyMesh.reset(finleyMesh,null_deleter());
58 }
59
60 //
61 // The copy constructor should just increment the use count
62 MeshAdapter::MeshAdapter(const MeshAdapter& in):
63 m_finleyMesh(in.m_finleyMesh)
64 {
65 setFunctionSpaceTypeNames();
66 }
67
68 MeshAdapter::~MeshAdapter()
69 {
70 //
71 // I hope the case for the pointer being zero has been taken care of.
72 // cout << "In MeshAdapter destructor." << endl;
73 if (m_finleyMesh.unique()) {
74 Finley_Mesh_free(m_finleyMesh.get());
75 }
76 }
77
78 int MeshAdapter::getMPISize() const
79 {
80 return m_finleyMesh.get()->MPIInfo->size;
81 }
82 int MeshAdapter::getMPIRank() const
83 {
84 return m_finleyMesh.get()->MPIInfo->rank;
85 }
86
87
88 Finley_Mesh* MeshAdapter::getFinley_Mesh() const {
89 return m_finleyMesh.get();
90 }
91
92 void MeshAdapter::write(const std::string& fileName) const
93 {
94 char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
95 strcpy(fName,fileName.c_str());
96 Finley_Mesh_write(m_finleyMesh.get(),fName);
97 checkFinleyError();
98 TMPMEMFREE(fName);
99 }
100
101 void MeshAdapter::Print_Mesh_Info(const bool full=false) const
102 {
103 Finley_PrintMesh_Info(m_finleyMesh.get(), full);
104 }
105
106 void MeshAdapter::dump(const std::string& fileName) const
107 {
108 #ifdef USE_NETCDF
109 const NcDim* ncdims[12];
110 NcVar *ids;
111 int *int_ptr;
112 Finley_Mesh *mesh = m_finleyMesh.get();
113 Finley_TagMap* tag_map;
114 int num_Tags = 0;
115 int mpi_size = mesh->MPIInfo->size;
116 int mpi_rank = mesh->MPIInfo->rank;
117 int numDim = mesh->Nodes->numDim;
118 int numNodes = mesh->Nodes->numNodes;
119 int num_Elements = mesh->Elements->numElements;
120 int num_FaceElements = mesh->FaceElements->numElements;
121 int num_ContactElements = mesh->ContactElements->numElements;
122 int num_Points = mesh->Points->numElements;
123 int num_Elements_numNodes = mesh->Elements->numNodes;
124 int num_FaceElements_numNodes = mesh->FaceElements->numNodes;
125 int num_ContactElements_numNodes = mesh->ContactElements->numNodes;
126
127 char *newFileName = Paso_MPI_appendRankToFileName(fileName.c_str(),
128 mpi_size, mpi_rank);
129
130 /* Figure out how much storage is required for tags */
131 tag_map = mesh->TagMap;
132 if (tag_map) {
133 while (tag_map) {
134 num_Tags++;
135 tag_map=tag_map->next;
136 }
137 }
138
139 // NetCDF error handler
140 NcError err(NcError::verbose_nonfatal);
141 // Create the file.
142 NcFile dataFile(newFileName, NcFile::Replace);
143 // check if writing was successful
144 if (!dataFile.is_valid())
145 throw DataException("Error - MeshAdapter::dump: opening of NetCDF file for output failed: " + *newFileName);
146
147 // Define dimensions (num_Elements and dim_Elements are identical, dim_Elements only appears if > 0)
148 if (! (ncdims[0] = dataFile.add_dim("numNodes", numNodes)) )
149 throw DataException("Error - MeshAdapter::dump: appending dimension numNodes to netCDF file failed: " + *newFileName);
150 if (! (ncdims[1] = dataFile.add_dim("numDim", numDim)) )
151 throw DataException("Error - MeshAdapter::dump: appending dimension numDim to netCDF file failed: " + *newFileName);
152 if (! (ncdims[2] = dataFile.add_dim("mpi_size_plus_1", mpi_size+1)) )
153 throw DataException("Error - MeshAdapter::dump: appending dimension mpi_size to netCDF file failed: " + *newFileName);
154 if (num_Elements>0)
155 if (! (ncdims[3] = dataFile.add_dim("dim_Elements", num_Elements)) )
156 throw DataException("Error - MeshAdapter::dump: appending dimension dim_Elements to netCDF file failed: " + *newFileName);
157 if (num_FaceElements>0)
158 if (! (ncdims[4] = dataFile.add_dim("dim_FaceElements", num_FaceElements)) )
159 throw DataException("Error - MeshAdapter::dump: appending dimension dim_FaceElements to netCDF file failed: " + *newFileName);
160 if (num_ContactElements>0)
161 if (! (ncdims[5] = dataFile.add_dim("dim_ContactElements", num_ContactElements)) )
162 throw DataException("Error - MeshAdapter::dump: appending dimension dim_ContactElements to netCDF file failed: " + *newFileName);
163 if (num_Points>0)
164 if (! (ncdims[6] = dataFile.add_dim("dim_Points", num_Points)) )
165 throw DataException("Error - MeshAdapter::dump: appending dimension dim_Points to netCDF file failed: " + *newFileName);
166 if (num_Elements>0)
167 if (! (ncdims[7] = dataFile.add_dim("dim_Elements_Nodes", num_Elements_numNodes)) )
168 throw DataException("Error - MeshAdapter::dump: appending dimension dim_Elements_Nodes to netCDF file failed: " + *newFileName);
169 if (num_FaceElements>0)
170 if (! (ncdims[8] = dataFile.add_dim("dim_FaceElements_numNodes", num_FaceElements_numNodes)) )
171 throw DataException("Error - MeshAdapter::dump: appending dimension dim_FaceElements_numNodes to netCDF file failed: " + *newFileName);
172 if (num_ContactElements>0)
173 if (! (ncdims[9] = dataFile.add_dim("dim_ContactElements_numNodes", num_ContactElements_numNodes)) )
174 throw DataException("Error - MeshAdapter::dump: appending dimension dim_ContactElements_numNodes to netCDF file failed: " + *newFileName);
175 if (num_Tags>0)
176 if (! (ncdims[10] = dataFile.add_dim("dim_Tags", num_Tags)) )
177 throw DataException("Error - MeshAdapter::dump: appending dimension dim_Tags to netCDF file failed: " + *newFileName);
178
179 // Attributes: MPI size, MPI rank, Name, order, reduced_order
180 if (!dataFile.add_att("mpi_size", mpi_size) )
181 throw DataException("Error - MeshAdapter::dump: appending mpi_size to NetCDF file failed: " + *newFileName);
182 if (!dataFile.add_att("mpi_rank", mpi_rank) )
183 throw DataException("Error - MeshAdapter::dump: appending mpi_rank to NetCDF file failed: " + *newFileName);
184 if (!dataFile.add_att("Name",mesh->Name) )
185 throw DataException("Error - MeshAdapter::dump: appending Name to NetCDF file failed: " + *newFileName);
186 if (!dataFile.add_att("numDim",numDim) )
187 throw DataException("Error - MeshAdapter::dump: appending order to NetCDF file failed: " + *newFileName);
188 if (!dataFile.add_att("order",mesh->order) )
189 throw DataException("Error - MeshAdapter::dump: appending order to NetCDF file failed: " + *newFileName);
190 if (!dataFile.add_att("reduced_order",mesh->reduced_order) )
191 throw DataException("Error - MeshAdapter::dump: appending reduced_order to NetCDF file failed: " + *newFileName);
192 if (!dataFile.add_att("numNodes",numNodes) )
193 throw DataException("Error - MeshAdapter::dump: appending numNodes to NetCDF file failed: " + *newFileName);
194 if (!dataFile.add_att("num_Elements",num_Elements) )
195 throw DataException("Error - MeshAdapter::dump: appending num_Elements to NetCDF file failed: " + *newFileName);
196 if (!dataFile.add_att("num_FaceElements",num_FaceElements) )
197 throw DataException("Error - MeshAdapter::dump: appending num_FaceElements to NetCDF file failed: " + *newFileName);
198 if (!dataFile.add_att("num_ContactElements",num_ContactElements) )
199 throw DataException("Error - MeshAdapter::dump: appending num_ContactElements to NetCDF file failed: " + *newFileName);
200 if (!dataFile.add_att("num_Points",num_Points) )
201 throw DataException("Error - MeshAdapter::dump: appending num_Points to NetCDF file failed: " + *newFileName);
202 if (!dataFile.add_att("num_Elements_numNodes",num_Elements_numNodes) )
203 throw DataException("Error - MeshAdapter::dump: appending num_Elements_numNodes to NetCDF file failed: " + *newFileName);
204 if (!dataFile.add_att("num_FaceElements_numNodes",num_FaceElements_numNodes) )
205 throw DataException("Error - MeshAdapter::dump: appending num_FaceElements_numNodes to NetCDF file failed: " + *newFileName);
206 if (!dataFile.add_att("num_ContactElements_numNodes",num_ContactElements_numNodes) )
207 throw DataException("Error - MeshAdapter::dump: appending num_ContactElements_numNodes to NetCDF file failed: " + *newFileName);
208 if (!dataFile.add_att("Elements_TypeId", mesh->Elements->ReferenceElement->Type->TypeId) )
209 throw DataException("Error - MeshAdapter::dump: appending Elements_TypeId to NetCDF file failed: " + *newFileName);
210 if (!dataFile.add_att("FaceElements_TypeId", mesh->FaceElements->ReferenceElement->Type->TypeId) )
211 throw DataException("Error - MeshAdapter::dump: appending FaceElements_TypeId to NetCDF file failed: " + *newFileName);
212 if (!dataFile.add_att("ContactElements_TypeId", mesh->ContactElements->ReferenceElement->Type->TypeId) )
213 throw DataException("Error - MeshAdapter::dump: appending ContactElements_TypeId to NetCDF file failed: " + *newFileName);
214 if (!dataFile.add_att("Points_TypeId", mesh->Points->ReferenceElement->Type->TypeId) )
215 throw DataException("Error - MeshAdapter::dump: appending Points_TypeId to NetCDF file failed: " + *newFileName);
216 if (!dataFile.add_att("num_Tags", num_Tags) )
217 throw DataException("Error - MeshAdapter::dump: appending num_Tags to NetCDF file failed: " + *newFileName);
218
219 // // // // // Nodes // // // // //
220
221 // Only write nodes if non-empty because NetCDF doesn't like empty arrays (it treats them as NC_UNLIMITED)
222 if (numNodes>0) {
223
224 // Nodes Id
225 if (! ( ids = dataFile.add_var("Nodes_Id", ncInt, ncdims[0])) )
226 throw DataException("Error - MeshAdapter::dump: appending Nodes_Id to netCDF file failed: " + *newFileName);
227 int_ptr = &mesh->Nodes->Id[0];
228 if (! (ids->put(int_ptr, numNodes)) )
229 throw DataException("Error - MeshAdapter::dump: copy Nodes_Id to netCDF buffer failed: " + *newFileName);
230
231 // Nodes Tag
232 if (! ( ids = dataFile.add_var("Nodes_Tag", ncInt, ncdims[0])) )
233 throw DataException("Error - MeshAdapter::dump: appending Nodes_Tag to netCDF file failed: " + *newFileName);
234 int_ptr = &mesh->Nodes->Tag[0];
235 if (! (ids->put(int_ptr, numNodes)) )
236 throw DataException("Error - MeshAdapter::dump: copy Nodes_Tag to netCDF buffer failed: " + *newFileName);
237
238 // Nodes gDOF
239 if (! ( ids = dataFile.add_var("Nodes_gDOF", ncInt, ncdims[0])) )
240 throw DataException("Error - MeshAdapter::dump: appending Nodes_gDOF to netCDF file failed: " + *newFileName);
241 int_ptr = &mesh->Nodes->globalDegreesOfFreedom[0];
242 if (! (ids->put(int_ptr, numNodes)) )
243 throw DataException("Error - MeshAdapter::dump: copy Nodes_gDOF to netCDF buffer failed: " + *newFileName);
244
245 // Nodes global node index
246 if (! ( ids = dataFile.add_var("Nodes_gNI", ncInt, ncdims[0])) )
247 throw DataException("Error - MeshAdapter::dump: appending Nodes_gNI to netCDF file failed: " + *newFileName);
248 int_ptr = &mesh->Nodes->globalNodesIndex[0];
249 if (! (ids->put(int_ptr, numNodes)) )
250 throw DataException("Error - MeshAdapter::dump: copy Nodes_gNI to netCDF buffer failed: " + *newFileName);
251
252 // Nodes grDof
253 if (! ( ids = dataFile.add_var("Nodes_grDfI", ncInt, ncdims[0])) )
254 throw DataException("Error - MeshAdapter::dump: appending Nodes_grDfI to netCDF file failed: " + *newFileName);
255 int_ptr = &mesh->Nodes->globalReducedDOFIndex[0];
256 if (! (ids->put(int_ptr, numNodes)) )
257 throw DataException("Error - MeshAdapter::dump: copy Nodes_grDfI to netCDF buffer failed: " + *newFileName);
258
259 // Nodes grNI
260 if (! ( ids = dataFile.add_var("Nodes_grNI", ncInt, ncdims[0])) )
261 throw DataException("Error - MeshAdapter::dump: appending Nodes_grNI to netCDF file failed: " + *newFileName);
262 int_ptr = &mesh->Nodes->globalReducedNodesIndex[0];
263 if (! (ids->put(int_ptr, numNodes)) )
264 throw DataException("Error - MeshAdapter::dump: copy Nodes_grNI to netCDF buffer failed: " + *newFileName);
265
266 // Nodes Coordinates
267 if (! ( ids = dataFile.add_var("Nodes_Coordinates", ncDouble, ncdims[0], ncdims[1]) ) )
268 throw DataException("Error - MeshAdapter::dump: appending Nodes_Coordinates to netCDF file failed: " + *newFileName);
269 if (! (ids->put(&(mesh->Nodes->Coordinates[INDEX2(0,0,numDim)]), numNodes, numDim)) )
270 throw DataException("Error - MeshAdapter::dump: copy Nodes_Coordinates to netCDF buffer failed: " + *newFileName);
271
272 // Nodes degreesOfFreedomDistribution
273 if (! ( ids = dataFile.add_var("Nodes_DofDistribution", ncInt, ncdims[2])) )
274 throw DataException("Error - MeshAdapter::dump: appending Nodes_DofDistribution to netCDF file failed: " + *newFileName);
275 int_ptr = &mesh->Nodes->degreesOfFreedomDistribution->first_component[0];
276 if (! (ids->put(int_ptr, mpi_size+1)) )
277 throw DataException("Error - MeshAdapter::dump: copy Nodes_DofDistribution to netCDF buffer failed: " + *newFileName);
278
279 // Nodes nodeDistribution
280 if (! ( ids = dataFile.add_var("Nodes_NodeDistribution", ncInt, ncdims[2])) )
281 throw DataException("Error - MeshAdapter::dump: appending Nodes_NodeDistribution to netCDF file failed: " + *newFileName);
282 int_ptr = &mesh->Nodes->nodesDistribution->first_component[0];
283 if (! (ids->put(int_ptr, mpi_size+1)) )
284 throw DataException("Error - MeshAdapter::dump: copy Nodes_NodeDistribution to netCDF buffer failed: " + *newFileName);
285
286 }
287
288 // // // // // Elements // // // // //
289
290 if (num_Elements>0) {
291
292 // Elements_Id
293 if (! ( ids = dataFile.add_var("Elements_Id", ncInt, ncdims[3])) )
294 throw DataException("Error - MeshAdapter::dump: appending Elements_Id to netCDF file failed: " + *newFileName);
295 int_ptr = &mesh->Elements->Id[0];
296 if (! (ids->put(int_ptr, num_Elements)) )
297 throw DataException("Error - MeshAdapter::dump: copy Elements_Id to netCDF buffer failed: " + *newFileName);
298
299 // Elements_Tag
300 if (! ( ids = dataFile.add_var("Elements_Tag", ncInt, ncdims[3])) )
301 throw DataException("Error - MeshAdapter::dump: appending Elements_Tag to netCDF file failed: " + *newFileName);
302 int_ptr = &mesh->Elements->Tag[0];
303 if (! (ids->put(int_ptr, num_Elements)) )
304 throw DataException("Error - MeshAdapter::dump: copy Elements_Tag to netCDF buffer failed: " + *newFileName);
305
306 // Elements_Owner
307 if (! ( ids = dataFile.add_var("Elements_Owner", ncInt, ncdims[3])) )
308 throw DataException("Error - MeshAdapter::dump: appending Elements_Owner to netCDF file failed: " + *newFileName);
309 int_ptr = &mesh->Elements->Owner[0];
310 if (! (ids->put(int_ptr, num_Elements)) )
311 throw DataException("Error - MeshAdapter::dump: copy Elements_Owner to netCDF buffer failed: " + *newFileName);
312
313 // Elements_Color
314 if (! ( ids = dataFile.add_var("Elements_Color", ncInt, ncdims[3])) )
315 throw DataException("Error - MeshAdapter::dump: appending Elements_Color to netCDF file failed: " + *newFileName);
316 int_ptr = &mesh->Elements->Color[0];
317 if (! (ids->put(int_ptr, num_Elements)) )
318 throw DataException("Error - MeshAdapter::dump: copy Elements_Color to netCDF buffer failed: " + *newFileName);
319
320 // Elements_Nodes
321 if (! ( ids = dataFile.add_var("Elements_Nodes", ncInt, ncdims[3], ncdims[7]) ) )
322 throw DataException("Error - MeshAdapter::dump: appending Elements_Nodes to netCDF file failed: " + *newFileName);
323 if (! (ids->put(&(mesh->Elements->Nodes[0]), num_Elements, num_Elements_numNodes)) )
324 throw DataException("Error - MeshAdapter::dump: copy Elements_Nodes to netCDF buffer failed: " + *newFileName);
325
326 }
327
328 // // // // // Face_Elements // // // // //
329
330 if (num_FaceElements>0) {
331
332 // FaceElements_Id
333 if (! ( ids = dataFile.add_var("FaceElements_Id", ncInt, ncdims[4])) )
334 throw DataException("Error - MeshAdapter::dump: appending FaceElements_Id to netCDF file failed: " + *newFileName);
335 int_ptr = &mesh->FaceElements->Id[0];
336 if (! (ids->put(int_ptr, num_FaceElements)) )
337 throw DataException("Error - MeshAdapter::dump: copy FaceElements_Id to netCDF buffer failed: " + *newFileName);
338
339 // FaceElements_Tag
340 if (! ( ids = dataFile.add_var("FaceElements_Tag", ncInt, ncdims[4])) )
341 throw DataException("Error - MeshAdapter::dump: appending FaceElements_Tag to netCDF file failed: " + *newFileName);
342 int_ptr = &mesh->FaceElements->Tag[0];
343 if (! (ids->put(int_ptr, num_FaceElements)) )
344 throw DataException("Error - MeshAdapter::dump: copy FaceElements_Tag to netCDF buffer failed: " + *newFileName);
345
346 // FaceElements_Owner
347 if (! ( ids = dataFile.add_var("FaceElements_Owner", ncInt, ncdims[4])) )
348 throw DataException("Error - MeshAdapter::dump: appending FaceElements_Owner to netCDF file failed: " + *newFileName);
349 int_ptr = &mesh->FaceElements->Owner[0];
350 if (! (ids->put(int_ptr, num_FaceElements)) )
351 throw DataException("Error - MeshAdapter::dump: copy FaceElements_Owner to netCDF buffer failed: " + *newFileName);
352
353 // FaceElements_Color
354 if (! ( ids = dataFile.add_var("FaceElements_Color", ncInt, ncdims[4])) )
355 throw DataException("Error - MeshAdapter::dump: appending FaceElements_Color to netCDF file failed: " + *newFileName);
356 int_ptr = &mesh->FaceElements->Color[0];
357 if (! (ids->put(int_ptr, num_FaceElements)) )
358 throw DataException("Error - MeshAdapter::dump: copy FaceElements_Color to netCDF buffer failed: " + *newFileName);
359
360 // FaceElements_Nodes
361 if (! ( ids = dataFile.add_var("FaceElements_Nodes", ncInt, ncdims[4], ncdims[8]) ) )
362 throw DataException("Error - MeshAdapter::dump: appending FaceElements_Nodes to netCDF file failed: " + *newFileName);
363 if (! (ids->put(&(mesh->FaceElements->Nodes[0]), num_FaceElements, num_FaceElements_numNodes)) )
364 throw DataException("Error - MeshAdapter::dump: copy FaceElements_Nodes to netCDF buffer failed: " + *newFileName);
365
366 }
367
368 // // // // // Contact_Elements // // // // //
369
370 if (num_ContactElements>0) {
371
372 // ContactElements_Id
373 if (! ( ids = dataFile.add_var("ContactElements_Id", ncInt, ncdims[5])) )
374 throw DataException("Error - MeshAdapter::dump: appending ContactElements_Id to netCDF file failed: " + *newFileName);
375 int_ptr = &mesh->ContactElements->Id[0];
376 if (! (ids->put(int_ptr, num_ContactElements)) )
377 throw DataException("Error - MeshAdapter::dump: copy ContactElements_Id to netCDF buffer failed: " + *newFileName);
378
379 // ContactElements_Tag
380 if (! ( ids = dataFile.add_var("ContactElements_Tag", ncInt, ncdims[5])) )
381 throw DataException("Error - MeshAdapter::dump: appending ContactElements_Tag to netCDF file failed: " + *newFileName);
382 int_ptr = &mesh->ContactElements->Tag[0];
383 if (! (ids->put(int_ptr, num_ContactElements)) )
384 throw DataException("Error - MeshAdapter::dump: copy ContactElements_Tag to netCDF buffer failed: " + *newFileName);
385
386 // ContactElements_Owner
387 if (! ( ids = dataFile.add_var("ContactElements_Owner", ncInt, ncdims[5])) )
388 throw DataException("Error - MeshAdapter::dump: appending ContactElements_Owner to netCDF file failed: " + *newFileName);
389 int_ptr = &mesh->ContactElements->Owner[0];
390 if (! (ids->put(int_ptr, num_ContactElements)) )
391 throw DataException("Error - MeshAdapter::dump: copy ContactElements_Owner to netCDF buffer failed: " + *newFileName);
392
393 // ContactElements_Color
394 if (! ( ids = dataFile.add_var("ContactElements_Color", ncInt, ncdims[5])) )
395 throw DataException("Error - MeshAdapter::dump: appending ContactElements_Color to netCDF file failed: " + *newFileName);
396 int_ptr = &mesh->ContactElements->Color[0];
397 if (! (ids->put(int_ptr, num_ContactElements)) )
398 throw DataException("Error - MeshAdapter::dump: copy ContactElements_Color to netCDF buffer failed: " + *newFileName);
399
400 // ContactElements_Nodes
401 if (! ( ids = dataFile.add_var("ContactElements_Nodes", ncInt, ncdims[5], ncdims[9]) ) )
402 throw DataException("Error - MeshAdapter::dump: appending ContactElements_Nodes to netCDF file failed: " + *newFileName);
403 if (! (ids->put(&(mesh->ContactElements->Nodes[0]), num_ContactElements, num_ContactElements_numNodes)) )
404 throw DataException("Error - MeshAdapter::dump: copy ContactElements_Nodes to netCDF buffer failed: " + *newFileName);
405
406 }
407
408 // // // // // Points // // // // //
409
410 if (num_Points>0) {
411
412 fprintf(stderr, "\n\n\nWARNING: MeshAdapter::dump has not been tested with Point elements\n\n\n");
413
414 // Points_Id
415 if (! ( ids = dataFile.add_var("Points_Id", ncInt, ncdims[6])) )
416 throw DataException("Error - MeshAdapter::dump: appending Points_Id to netCDF file failed: " + *newFileName);
417 int_ptr = &mesh->Points->Id[0];
418 if (! (ids->put(int_ptr, num_Points)) )
419 throw DataException("Error - MeshAdapter::dump: copy Points_Id to netCDF buffer failed: " + *newFileName);
420
421 // Points_Tag
422 if (! ( ids = dataFile.add_var("Points_Tag", ncInt, ncdims[6])) )
423 throw DataException("Error - MeshAdapter::dump: appending Points_Tag to netCDF file failed: " + *newFileName);
424 int_ptr = &mesh->Points->Tag[0];
425 if (! (ids->put(int_ptr, num_Points)) )
426 throw DataException("Error - MeshAdapter::dump: copy Points_Tag to netCDF buffer failed: " + *newFileName);
427
428 // Points_Owner
429 if (! ( ids = dataFile.add_var("Points_Owner", ncInt, ncdims[6])) )
430 throw DataException("Error - MeshAdapter::dump: appending Points_Owner to netCDF file failed: " + *newFileName);
431 int_ptr = &mesh->Points->Owner[0];
432 if (! (ids->put(int_ptr, num_Points)) )
433 throw DataException("Error - MeshAdapter::dump: copy Points_Owner to netCDF buffer failed: " + *newFileName);
434
435 // Points_Color
436 if (! ( ids = dataFile.add_var("Points_Color", ncInt, ncdims[6])) )
437 throw DataException("Error - MeshAdapter::dump: appending Points_Color to netCDF file failed: " + *newFileName);
438 int_ptr = &mesh->Points->Color[0];
439 if (! (ids->put(int_ptr, num_Points)) )
440 throw DataException("Error - MeshAdapter::dump: copy Points_Color to netCDF buffer failed: " + *newFileName);
441
442 // Points_Nodes
443 // mesh->Nodes->Id[mesh->Points->Nodes[INDEX2(0,i,1)]]
444 if (! ( ids = dataFile.add_var("Points_Nodes", ncInt, ncdims[6]) ) )
445 throw DataException("Error - MeshAdapter::dump: appending Points_Nodes to netCDF file failed: " + *newFileName);
446 if (! (ids->put(&(mesh->Points->Nodes[0]), num_Points)) )
447 throw DataException("Error - MeshAdapter::dump: copy Points_Nodes to netCDF buffer failed: " + *newFileName);
448
449 }
450
451 // // // // // TagMap // // // // //
452
453 if (num_Tags>0) {
454
455 // Temp storage to gather node IDs
456 int *Tags_keys = TMPMEMALLOC(num_Tags, int);
457 char name_temp[4096];
458
459 /* Copy tag data into temp arrays */
460 tag_map = mesh->TagMap;
461 if (tag_map) {
462 int i = 0;
463 while (tag_map) {
464 Tags_keys[i++] = tag_map->tag_key;
465 tag_map=tag_map->next;
466 }
467 }
468
469 // Tags_keys
470 if (! ( ids = dataFile.add_var("Tags_keys", ncInt, ncdims[10])) )
471 throw DataException("Error - MeshAdapter::dump: appending Tags_keys to netCDF file failed: " + *newFileName);
472 int_ptr = &Tags_keys[0];
473 if (! (ids->put(int_ptr, num_Tags)) )
474 throw DataException("Error - MeshAdapter::dump: copy Tags_keys to netCDF buffer failed: " + *newFileName);
475
476 // Tags_names_*
477 // This is an array of strings, it should be stored as an array but instead I have hacked in one attribute per string
478 // because the NetCDF manual doesn't tell how to do an array of strings
479 tag_map = mesh->TagMap;
480 if (tag_map) {
481 int i = 0;
482 while (tag_map) {
483 sprintf(name_temp, "Tags_name_%d", i);
484 if (!dataFile.add_att(name_temp, tag_map->name) )
485 throw DataException("Error - MeshAdapter::dump: appending Tags_names_ to NetCDF file failed: " + *newFileName);
486 tag_map=tag_map->next;
487 i++;
488 }
489 }
490
491 TMPMEMFREE(Tags_keys);
492
493 }
494
495
496 // NetCDF file is closed by destructor of NcFile object
497 #else
498 Finley_setError(IO_ERROR, "MeshAdapter::dump: not configured with NetCDF. Please contact your installation manager.");
499 #endif /* USE_NETCDF */
500 checkFinleyError();
501 }
502
503 string MeshAdapter::getDescription() const
504 {
505 return "FinleyMesh";
506 }
507
508 string MeshAdapter::functionSpaceTypeAsString(int functionSpaceType) const
509 {
510 FunctionSpaceNamesMapType::iterator loc;
511 loc=m_functionSpaceTypeNames.find(functionSpaceType);
512 if (loc==m_functionSpaceTypeNames.end()) {
513 return "Invalid function space type code.";
514 } else {
515 return loc->second;
516 }
517 }
518
519 bool MeshAdapter::isValidFunctionSpaceType(int functionSpaceType) const
520 {
521 FunctionSpaceNamesMapType::iterator loc;
522 loc=m_functionSpaceTypeNames.find(functionSpaceType);
523 return (loc!=m_functionSpaceTypeNames.end());
524 }
525
526 void MeshAdapter::setFunctionSpaceTypeNames()
527 {
528 m_functionSpaceTypeNames.insert
529 (FunctionSpaceNamesMapType::value_type(DegreesOfFreedom,"Finley_DegreesOfFreedom"));
530 m_functionSpaceTypeNames.insert
531 (FunctionSpaceNamesMapType::value_type(ReducedDegreesOfFreedom,"Finley_ReducedDegreesOfFreedom"));
532 m_functionSpaceTypeNames.insert
533 (FunctionSpaceNamesMapType::value_type(Nodes,"Finley_Nodes"));
534 m_functionSpaceTypeNames.insert
535 (FunctionSpaceNamesMapType::value_type(ReducedNodes,"Finley_Reduced_Nodes"));
536 m_functionSpaceTypeNames.insert
537 (FunctionSpaceNamesMapType::value_type(Elements,"Finley_Elements"));
538 m_functionSpaceTypeNames.insert
539 (FunctionSpaceNamesMapType::value_type(ReducedElements,"Finley_Reduced_Elements"));
540 m_functionSpaceTypeNames.insert
541 (FunctionSpaceNamesMapType::value_type(FaceElements,"Finley_Face_Elements"));
542 m_functionSpaceTypeNames.insert
543 (FunctionSpaceNamesMapType::value_type(ReducedFaceElements,"Finley_Reduced_Face_Elements"));
544 m_functionSpaceTypeNames.insert
545 (FunctionSpaceNamesMapType::value_type(Points,"Finley_Points"));
546 m_functionSpaceTypeNames.insert
547 (FunctionSpaceNamesMapType::value_type(ContactElementsZero,"Finley_Contact_Elements_0"));
548 m_functionSpaceTypeNames.insert
549 (FunctionSpaceNamesMapType::value_type(ReducedContactElementsZero,"Finley_Reduced_Contact_Elements_0"));
550 m_functionSpaceTypeNames.insert
551 (FunctionSpaceNamesMapType::value_type(ContactElementsOne,"Finley_Contact_Elements_1"));
552 m_functionSpaceTypeNames.insert
553 (FunctionSpaceNamesMapType::value_type(ReducedContactElementsOne,"Finley_Reduced_Contact_Elements_1"));
554 }
555
556 int MeshAdapter::getContinuousFunctionCode() const
557 {
558 return Nodes;
559 }
560 int MeshAdapter::getReducedContinuousFunctionCode() const
561 {
562 return ReducedNodes;
563 }
564
565 int MeshAdapter::getFunctionCode() const
566 {
567 return Elements;
568 }
569 int MeshAdapter::getReducedFunctionCode() const
570 {
571 return ReducedElements;
572 }
573
574 int MeshAdapter::getFunctionOnBoundaryCode() const
575 {
576 return FaceElements;
577 }
578 int MeshAdapter::getReducedFunctionOnBoundaryCode() const
579 {
580 return ReducedFaceElements;
581 }
582
583 int MeshAdapter::getFunctionOnContactZeroCode() const
584 {
585 return ContactElementsZero;
586 }
587 int MeshAdapter::getReducedFunctionOnContactZeroCode() const
588 {
589 return ReducedContactElementsZero;
590 }
591
592 int MeshAdapter::getFunctionOnContactOneCode() const
593 {
594 return ContactElementsOne;
595 }
596 int MeshAdapter::getReducedFunctionOnContactOneCode() const
597 {
598 return ReducedContactElementsOne;
599 }
600
601 int MeshAdapter::getSolutionCode() const
602 {
603 return DegreesOfFreedom;
604 }
605
606 int MeshAdapter::getReducedSolutionCode() const
607 {
608 return ReducedDegreesOfFreedom;
609 }
610
611 int MeshAdapter::getDiracDeltaFunctionCode() const
612 {
613 return Points;
614 }
615
616 //
617 // return the spatial dimension of the Mesh:
618 //
619 int MeshAdapter::getDim() const
620 {
621 int numDim=Finley_Mesh_getDim(m_finleyMesh.get());
622 checkFinleyError();
623 return numDim;
624 }
625
626 //
627 // Return the number of data points summed across all MPI processes
628 //
629 int MeshAdapter::getNumDataPointsGlobal() const
630 {
631 return Finley_NodeFile_getGlobalNumNodes(m_finleyMesh.get()->Nodes);
632 }
633
634 //
635 // return the number of data points per sample and the number of samples
636 // needed to represent data on a parts of the mesh.
637 //
638 pair<int,int> MeshAdapter::getDataShape(int functionSpaceCode) const
639 {
640 int numDataPointsPerSample=0;
641 int numSamples=0;
642 Finley_Mesh* mesh=m_finleyMesh.get();
643 switch (functionSpaceCode) {
644 case(Nodes):
645 numDataPointsPerSample=1;
646 numSamples=Finley_NodeFile_getNumNodes(mesh->Nodes);
647 break;
648 case(ReducedNodes):
649 numDataPointsPerSample=1;
650 numSamples=Finley_NodeFile_getNumReducedNodes(mesh->Nodes);
651 break;
652 case(Elements):
653 if (mesh->Elements!=NULL) {
654 numSamples=mesh->Elements->numElements;
655 numDataPointsPerSample=mesh->Elements->ReferenceElement->numQuadNodes;
656 }
657 break;
658 case(ReducedElements):
659 if (mesh->Elements!=NULL) {
660 numSamples=mesh->Elements->numElements;
661 numDataPointsPerSample=mesh->Elements->ReferenceElementReducedOrder->numQuadNodes;
662 }
663 break;
664 case(FaceElements):
665 if (mesh->FaceElements!=NULL) {
666 numDataPointsPerSample=mesh->FaceElements->ReferenceElement->numQuadNodes;
667 numSamples=mesh->FaceElements->numElements;
668 }
669 break;
670 case(ReducedFaceElements):
671 if (mesh->FaceElements!=NULL) {
672 numDataPointsPerSample=mesh->FaceElements->ReferenceElementReducedOrder->numQuadNodes;
673 numSamples=mesh->FaceElements->numElements;
674 }
675 break;
676 case(Points):
677 if (mesh->Points!=NULL) {
678 numDataPointsPerSample=1;
679 numSamples=mesh->Points->numElements;
680 }
681 break;
682 case(ContactElementsZero):
683 if (mesh->ContactElements!=NULL) {
684 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
685 numSamples=mesh->ContactElements->numElements;
686 }
687 break;
688 case(ReducedContactElementsZero):
689 if (mesh->ContactElements!=NULL) {
690 numDataPointsPerSample=mesh->ContactElements->ReferenceElementReducedOrder->numQuadNodes;
691 numSamples=mesh->ContactElements->numElements;
692 }
693 break;
694 case(ContactElementsOne):
695 if (mesh->ContactElements!=NULL) {
696 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
697 numSamples=mesh->ContactElements->numElements;
698 }
699 break;
700 case(ReducedContactElementsOne):
701 if (mesh->ContactElements!=NULL) {
702 numDataPointsPerSample=mesh->ContactElements->ReferenceElementReducedOrder->numQuadNodes;
703 numSamples=mesh->ContactElements->numElements;
704 }
705 break;
706 case(DegreesOfFreedom):
707 if (mesh->Nodes!=NULL) {
708 numDataPointsPerSample=1;
709 numSamples=Finley_NodeFile_getNumDegreesOfFreedom(mesh->Nodes);
710 }
711 break;
712 case(ReducedDegreesOfFreedom):
713 if (mesh->Nodes!=NULL) {
714 numDataPointsPerSample=1;
715 numSamples=Finley_NodeFile_getNumReducedDegreesOfFreedom(mesh->Nodes);
716 }
717 break;
718 default:
719 stringstream temp;
720 temp << "Error - Invalid function space type: " << functionSpaceCode << " for domain: " << getDescription();
721 throw FinleyAdapterException(temp.str());
722 break;
723 }
724 return pair<int,int>(numDataPointsPerSample,numSamples);
725 }
726
727 //
728 // adds linear PDE of second order into a given stiffness matrix and right hand side:
729 //
730 void MeshAdapter::addPDEToSystem(
731 SystemMatrixAdapter& mat, escript::Data& rhs,
732 const escript::Data& A, const escript::Data& B, const escript::Data& C,const escript::Data& D,const escript::Data& X,const escript::Data& Y,
733 const escript::Data& d, const escript::Data& y,
734 const escript::Data& d_contact,const escript::Data& y_contact) const
735 {
736 escriptDataC _rhs=rhs.getDataC();
737 escriptDataC _A =A.getDataC();
738 escriptDataC _B=B.getDataC();
739 escriptDataC _C=C.getDataC();
740 escriptDataC _D=D.getDataC();
741 escriptDataC _X=X.getDataC();
742 escriptDataC _Y=Y.getDataC();
743 escriptDataC _d=d.getDataC();
744 escriptDataC _y=y.getDataC();
745 escriptDataC _d_contact=d_contact.getDataC();
746 escriptDataC _y_contact=y_contact.getDataC();
747
748 Finley_Mesh* mesh=m_finleyMesh.get();
749
750 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(), &_rhs, &_A, &_B, &_C, &_D, &_X, &_Y );
751 checkFinleyError();
752
753 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, mat.getPaso_SystemMatrix(), &_rhs, 0, 0, 0, &_d, 0, &_y );
754 checkFinleyError();
755
756 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, mat.getPaso_SystemMatrix(), &_rhs , 0, 0, 0, &_d_contact, 0, &_y_contact );
757 checkFinleyError();
758 }
759
760 void MeshAdapter::addPDEToLumpedSystem(
761 escript::Data& mat,
762 const escript::Data& D,
763 const escript::Data& d) const
764 {
765 escriptDataC _mat=mat.getDataC();
766 escriptDataC _D=D.getDataC();
767 escriptDataC _d=d.getDataC();
768
769 Finley_Mesh* mesh=m_finleyMesh.get();
770
771 Finley_Assemble_LumpedSystem(mesh->Nodes,mesh->Elements,&_mat, &_D);
772 Finley_Assemble_LumpedSystem(mesh->Nodes,mesh->FaceElements,&_mat, &_d);
773
774 checkFinleyError();
775 }
776
777
778 //
779 // adds linear PDE of second order into the right hand side only
780 //
781 void MeshAdapter::addPDEToRHS( escript::Data& rhs, const escript::Data& X,const escript::Data& Y, const escript::Data& y, const escript::Data& y_contact) const
782 {
783 Finley_Mesh* mesh=m_finleyMesh.get();
784
785 escriptDataC _rhs=rhs.getDataC();
786 escriptDataC _X=X.getDataC();
787 escriptDataC _Y=Y.getDataC();
788 escriptDataC _y=y.getDataC();
789 escriptDataC _y_contact=y_contact.getDataC();
790
791 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements, 0, &_rhs, 0, 0, 0, 0, &_X, &_Y );
792 checkFinleyError();
793
794 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, 0, &_rhs, 0, 0, 0, 0, 0, &_y );
795 checkFinleyError();
796
797 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, 0, &_rhs , 0, 0, 0, 0, 0, &_y_contact );
798 checkFinleyError();
799 }
800 //
801 // adds PDE of second order into a transport problem
802 //
803 void MeshAdapter::addPDEToTransportProblem(
804 TransportProblemAdapter& tp, escript::Data& source, const escript::Data& M,
805 const escript::Data& A, const escript::Data& B, const escript::Data& C,
806 const escript::Data& D,const escript::Data& X,const escript::Data& Y,
807 const escript::Data& d, const escript::Data& y,
808 const escript::Data& d_contact,const escript::Data& y_contact) const
809 {
810 DataTypes::ShapeType shape;
811 source.expand();
812 escriptDataC _source=source.getDataC();
813 escriptDataC _M=M.getDataC();
814 escriptDataC _A=A.getDataC();
815 escriptDataC _B=B.getDataC();
816 escriptDataC _C=C.getDataC();
817 escriptDataC _D=D.getDataC();
818 escriptDataC _X=X.getDataC();
819 escriptDataC _Y=Y.getDataC();
820 escriptDataC _d=d.getDataC();
821 escriptDataC _y=y.getDataC();
822 escriptDataC _d_contact=d_contact.getDataC();
823 escriptDataC _y_contact=y_contact.getDataC();
824
825 Finley_Mesh* mesh=m_finleyMesh.get();
826 Paso_FCTransportProblem* _tp = tp.getPaso_FCTransportProblem();
827
828 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->mass_matrix, &_source, 0, 0, 0, &_M, 0, 0 );
829 checkFinleyError();
830
831 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->transport_matrix, &_source, &_A, &_B, &_C, &_D, &_X, &_Y );
832 checkFinleyError();
833
834 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, _tp->transport_matrix, &_source, 0, 0, 0, &_d, 0, &_y );
835 checkFinleyError();
836
837 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, _tp->transport_matrix, &_source , 0, 0, 0, &_d_contact, 0, &_y_contact );
838 checkFinleyError();
839 }
840
841 //
842 // interpolates data between different function spaces:
843 //
844 void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
845 {
846 const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());
847 const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
848 if (inDomain!=*this)
849 throw FinleyAdapterException("Error - Illegal domain of interpolant.");
850 if (targetDomain!=*this)
851 throw FinleyAdapterException("Error - Illegal domain of interpolation target.");
852
853 Finley_Mesh* mesh=m_finleyMesh.get();
854 escriptDataC _target=target.getDataC();
855 escriptDataC _in=in.getDataC();
856 switch(in.getFunctionSpace().getTypeCode()) {
857 case(Nodes):
858 switch(target.getFunctionSpace().getTypeCode()) {
859 case(Nodes):
860 case(ReducedNodes):
861 case(DegreesOfFreedom):
862 case(ReducedDegreesOfFreedom):
863 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
864 break;
865 case(Elements):
866 case(ReducedElements):
867 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
868 break;
869 case(FaceElements):
870 case(ReducedFaceElements):
871 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
872 break;
873 case(Points):
874 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
875 break;
876 case(ContactElementsZero):
877 case(ReducedContactElementsZero):
878 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
879 break;
880 case(ContactElementsOne):
881 case(ReducedContactElementsOne):
882 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
883 break;
884 default:
885 stringstream temp;
886 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
887 throw FinleyAdapterException(temp.str());
888 break;
889 }
890 break;
891 case(ReducedNodes):
892 switch(target.getFunctionSpace().getTypeCode()) {
893 case(Nodes):
894 case(ReducedNodes):
895 case(DegreesOfFreedom):
896 case(ReducedDegreesOfFreedom):
897 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
898 break;
899 case(Elements):
900 case(ReducedElements):
901 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
902 break;
903 case(FaceElements):
904 case(ReducedFaceElements):
905 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
906 break;
907 case(Points):
908 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
909 break;
910 case(ContactElementsZero):
911 case(ReducedContactElementsZero):
912 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
913 break;
914 case(ContactElementsOne):
915 case(ReducedContactElementsOne):
916 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
917 break;
918 default:
919 stringstream temp;
920 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
921 throw FinleyAdapterException(temp.str());
922 break;
923 }
924 break;
925 case(Elements):
926 if (target.getFunctionSpace().getTypeCode()==Elements) {
927 Finley_Assemble_CopyElementData(mesh->Elements,&_target,&_in);
928 } else if (target.getFunctionSpace().getTypeCode()==ReducedElements) {
929 Finley_Assemble_AverageElementData(mesh->Elements,&_target,&_in);
930 } else {
931 throw FinleyAdapterException("Error - No interpolation with data on elements possible.");
932 }
933 break;
934 case(ReducedElements):
935 if (target.getFunctionSpace().getTypeCode()==ReducedElements) {
936 Finley_Assemble_CopyElementData(mesh->Elements,&_target,&_in);
937 } else {
938 throw FinleyAdapterException("Error - No interpolation with data on elements with reduced integration order possible.");
939 }
940 break;
941 case(FaceElements):
942 if (target.getFunctionSpace().getTypeCode()==FaceElements) {
943 Finley_Assemble_CopyElementData(mesh->FaceElements,&_target,&_in);
944 } else if (target.getFunctionSpace().getTypeCode()==ReducedFaceElements) {
945 Finley_Assemble_AverageElementData(mesh->FaceElements,&_target,&_in);
946 } else {
947 throw FinleyAdapterException("Error - No interpolation with data on face elements possible.");
948 }
949 break;
950 case(ReducedFaceElements):
951 if (target.getFunctionSpace().getTypeCode()==ReducedFaceElements) {
952 Finley_Assemble_CopyElementData(mesh->FaceElements,&_target,&_in);
953 } else {
954 throw FinleyAdapterException("Error - No interpolation with data on face elements with reduced integration order possible.");
955 }
956 break;
957 case(Points):
958 if (target.getFunctionSpace().getTypeCode()==Points) {
959 Finley_Assemble_CopyElementData(mesh->Points,&_target,&_in);
960 } else {
961 throw FinleyAdapterException("Error - No interpolation with data on points possible.");
962 }
963 break;
964 case(ContactElementsZero):
965 case(ContactElementsOne):
966 if (target.getFunctionSpace().getTypeCode()==ContactElementsZero || target.getFunctionSpace().getTypeCode()==ContactElementsOne) {
967 Finley_Assemble_CopyElementData(mesh->ContactElements,&_target,&_in);
968 } else if (target.getFunctionSpace().getTypeCode()==ReducedContactElementsZero || target.getFunctionSpace().getTypeCode()==ReducedContactElementsOne) {
969 Finley_Assemble_AverageElementData(mesh->ContactElements,&_target,&_in);
970 } else {
971 throw FinleyAdapterException("Error - No interpolation with data on contact elements possible.");
972 }
973 break;
974 case(ReducedContactElementsZero):
975 case(ReducedContactElementsOne):
976 if (target.getFunctionSpace().getTypeCode()==ReducedContactElementsZero || target.getFunctionSpace().getTypeCode()==ReducedContactElementsOne) {
977 Finley_Assemble_CopyElementData(mesh->ContactElements,&_target,&_in);
978 } else {
979 throw FinleyAdapterException("Error - No interpolation with data on contact elements with reduced integration order possible.");
980 }
981 break;
982 case(DegreesOfFreedom):
983 switch(target.getFunctionSpace().getTypeCode()) {
984 case(ReducedDegreesOfFreedom):
985 case(DegreesOfFreedom):
986 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
987 break;
988
989 case(Nodes):
990 case(ReducedNodes):
991 if (getMPISize()>1) {
992 escript::Data temp=escript::Data(in);
993 temp.expand();
994 escriptDataC _in2 = temp.getDataC();
995 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in2);
996 } else {
997 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
998 }
999 break;
1000 case(Elements):
1001 case(ReducedElements):
1002 if (getMPISize()>1) {
1003 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1004 escriptDataC _in2 = temp.getDataC();
1005 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);
1006 } else {
1007 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
1008 }
1009 break;
1010 case(FaceElements):
1011 case(ReducedFaceElements):
1012 if (getMPISize()>1) {
1013 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1014 escriptDataC _in2 = temp.getDataC();
1015 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1016
1017 } else {
1018 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
1019 }
1020 break;
1021 case(Points):
1022 if (getMPISize()>1) {
1023 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1024 escriptDataC _in2 = temp.getDataC();
1025 } else {
1026 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
1027 }
1028 break;
1029 case(ContactElementsZero):
1030 case(ContactElementsOne):
1031 case(ReducedContactElementsZero):
1032 case(ReducedContactElementsOne):
1033 if (getMPISize()>1) {
1034 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1035 escriptDataC _in2 = temp.getDataC();
1036 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1037 } else {
1038 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
1039 }
1040 break;
1041 default:
1042 stringstream temp;
1043 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
1044 throw FinleyAdapterException(temp.str());
1045 break;
1046 }
1047 break;
1048 case(ReducedDegreesOfFreedom):
1049 switch(target.getFunctionSpace().getTypeCode()) {
1050 case(Nodes):
1051 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to mesh nodes.");
1052 break;
1053 case(ReducedNodes):
1054 if (getMPISize()>1) {
1055 escript::Data temp=escript::Data(in);
1056 temp.expand();
1057 escriptDataC _in2 = temp.getDataC();
1058 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in2);
1059 } else {
1060 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1061 }
1062 break;
1063 case(DegreesOfFreedom):
1064 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");
1065 break;
1066 case(ReducedDegreesOfFreedom):
1067 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1068 break;
1069 case(Elements):
1070 case(ReducedElements):
1071 if (getMPISize()>1) {
1072 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1073 escriptDataC _in2 = temp.getDataC();
1074 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);
1075 } else {
1076 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
1077 }
1078 break;
1079 case(FaceElements):
1080 case(ReducedFaceElements):
1081 if (getMPISize()>1) {
1082 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1083 escriptDataC _in2 = temp.getDataC();
1084 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1085 } else {
1086 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
1087 }
1088 break;
1089 case(Points):
1090 if (getMPISize()>1) {
1091 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1092 escriptDataC _in2 = temp.getDataC();
1093 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in2,&_target);
1094 } else {
1095 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
1096 }
1097 break;
1098 case(ContactElementsZero):
1099 case(ContactElementsOne):
1100 case(ReducedContactElementsZero):
1101 case(ReducedContactElementsOne):
1102 if (getMPISize()>1) {
1103 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1104 escriptDataC _in2 = temp.getDataC();
1105 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1106 } else {
1107 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
1108 }
1109 break;
1110 default:
1111 stringstream temp;
1112 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
1113 throw FinleyAdapterException(temp.str());
1114 break;
1115 }
1116 break;
1117 default:
1118 stringstream temp;
1119 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type %d" << in.getFunctionSpace().getTypeCode();
1120 throw FinleyAdapterException(temp.str());
1121 break;
1122 }
1123 checkFinleyError();
1124 }
1125
1126 //
1127 // copies the locations of sample points into x:
1128 //
1129 void MeshAdapter::setToX(escript::Data& arg) const
1130 {
1131 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
1132 if (argDomain!=*this)
1133 throw FinleyAdapterException("Error - Illegal domain of data point locations");
1134 Finley_Mesh* mesh=m_finleyMesh.get();
1135 // in case of values node coordinates we can do the job directly:
1136 if (arg.getFunctionSpace().getTypeCode()==Nodes) {
1137 escriptDataC _arg=arg.getDataC();
1138 Finley_Assemble_NodeCoordinates(mesh->Nodes,&_arg);
1139 } else {
1140 escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
1141 escriptDataC _tmp_data=tmp_data.getDataC();
1142 Finley_Assemble_NodeCoordinates(mesh->Nodes,&_tmp_data);
1143 // this is then interpolated onto arg:
1144 interpolateOnDomain(arg,tmp_data);
1145 }
1146 checkFinleyError();
1147 }
1148
1149 //
1150 // return the normal vectors at the location of data points as a Data object:
1151 //
1152 void MeshAdapter::setToNormal(escript::Data& normal) const
1153 {
1154 const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());
1155 if (normalDomain!=*this)
1156 throw FinleyAdapterException("Error - Illegal domain of normal locations");
1157 Finley_Mesh* mesh=m_finleyMesh.get();
1158 escriptDataC _normal=normal.getDataC();
1159 switch(normal.getFunctionSpace().getTypeCode()) {
1160 case(Nodes):
1161 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for nodes");
1162 break;
1163 case(ReducedNodes):
1164 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced nodes");
1165 break;
1166 case(Elements):
1167 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements");
1168 break;
1169 case(ReducedElements):
1170 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements with reduced integration order");
1171 break;
1172 case (FaceElements):
1173 Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&_normal);
1174 break;
1175 case (ReducedFaceElements):
1176 Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&_normal);
1177 break;
1178 case(Points):
1179 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for point elements");
1180 break;
1181 case (ContactElementsOne):
1182 case (ContactElementsZero):
1183 Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&_normal);
1184 break;
1185 case (ReducedContactElementsOne):
1186 case (ReducedContactElementsZero):
1187 Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&_normal);
1188 break;
1189 case(DegreesOfFreedom):
1190 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for degrees of freedom.");
1191 break;
1192 case(ReducedDegreesOfFreedom):
1193 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced degrees of freedom.");
1194 break;
1195 default:
1196 stringstream temp;
1197 temp << "Error - Normal Vectors: Finley does not know anything about function space type " << normal.getFunctionSpace().getTypeCode();
1198 throw FinleyAdapterException(temp.str());
1199 break;
1200 }
1201 checkFinleyError();
1202 }
1203
1204 //
1205 // interpolates data to other domain:
1206 //
1207 void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
1208 {
1209 const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
1210 if (targetDomain!=*this)
1211 throw FinleyAdapterException("Error - Illegal domain of interpolation target");
1212
1213 throw FinleyAdapterException("Error - Finley does not allow interpolation across domains yet.");
1214 }
1215
1216 //
1217 // calculates the integral of a function defined of arg:
1218 //
1219 void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
1220 {
1221 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
1222 if (argDomain!=*this)
1223 throw FinleyAdapterException("Error - Illegal domain of integration kernel");
1224
1225 double blocktimer_start = blocktimer_time();
1226 Finley_Mesh* mesh=m_finleyMesh.get();
1227 escriptDataC _temp;
1228 escript::Data temp;
1229 escriptDataC _arg=arg.getDataC();
1230 switch(arg.getFunctionSpace().getTypeCode()) {
1231 case(Nodes):
1232 temp=escript::Data( arg, function(asAbstractContinuousDomain()) );
1233 _temp=temp.getDataC();
1234 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1235 break;
1236 case(ReducedNodes):
1237 temp=escript::Data( arg, function(asAbstractContinuousDomain()) );
1238 _temp=temp.getDataC();
1239 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1240 break;
1241 case(Elements):
1242 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_arg,&integrals[0]);
1243 break;
1244 case(ReducedElements):
1245 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_arg,&integrals[0]);
1246 break;
1247 case(FaceElements):
1248 Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&_arg,&integrals[0]);
1249 break;
1250 case(ReducedFaceElements):
1251 Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&_arg,&integrals[0]);
1252 break;
1253 case(Points):
1254 throw FinleyAdapterException("Error - Integral of data on points is not supported.");
1255 break;
1256 case(ContactElementsZero):
1257 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1258 break;
1259 case(ReducedContactElementsZero):
1260 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1261 break;
1262 case(ContactElementsOne):
1263 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1264 break;
1265 case(ReducedContactElementsOne):
1266 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1267 break;
1268 case(DegreesOfFreedom):
1269 temp=escript::Data( arg, function(asAbstractContinuousDomain()) );
1270 _temp=temp.getDataC();
1271 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1272 break;
1273 case(ReducedDegreesOfFreedom):
1274 temp=escript::Data( arg, function(asAbstractContinuousDomain()) );
1275 _temp=temp.getDataC();
1276 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1277 break;
1278 default:
1279 stringstream temp;
1280 temp << "Error - Integrals: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
1281 throw FinleyAdapterException(temp.str());
1282 break;
1283 }
1284 checkFinleyError();
1285 blocktimer_increment("integrate()", blocktimer_start);
1286 }
1287
1288 //
1289 // calculates the gradient of arg:
1290 //
1291 void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
1292 {
1293 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
1294 if (argDomain!=*this)
1295 throw FinleyAdapterException("Error - Illegal domain of gradient argument");
1296 const MeshAdapter& gradDomain=dynamic_cast<const MeshAdapter&>(grad.getFunctionSpace().getDomain());
1297 if (gradDomain!=*this)
1298 throw FinleyAdapterException("Error - Illegal domain of gradient");
1299
1300 Finley_Mesh* mesh=m_finleyMesh.get();
1301 escriptDataC _grad=grad.getDataC();
1302 escriptDataC nodeDataC;
1303 escript::Data temp;
1304 if (getMPISize()>1) {
1305 if( arg.getFunctionSpace().getTypeCode() == DegreesOfFreedom ) {
1306 temp=escript::Data( arg, continuousFunction(asAbstractContinuousDomain()) );
1307 nodeDataC = temp.getDataC();
1308 } else if( arg.getFunctionSpace().getTypeCode() == ReducedDegreesOfFreedom ) {
1309 temp=escript::Data( arg, reducedContinuousFunction(asAbstractContinuousDomain()) );
1310 nodeDataC = temp.getDataC();
1311 } else {
1312 nodeDataC = arg.getDataC();
1313 }
1314 } else {
1315 nodeDataC = arg.getDataC();
1316 }
1317 switch(grad.getFunctionSpace().getTypeCode()) {
1318 case(Nodes):
1319 throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
1320 break;
1321 case(ReducedNodes):
1322 throw FinleyAdapterException("Error - Gradient at reduced nodes is not supported.");
1323 break;
1324 case(Elements):
1325 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&_grad,&nodeDataC);
1326 break;
1327 case(ReducedElements):
1328 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&_grad,&nodeDataC);
1329 break;
1330 case(FaceElements):
1331 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&_grad,&nodeDataC);
1332 break;
1333 case(ReducedFaceElements):
1334 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&_grad,&nodeDataC);
1335 break;
1336 case(Points):
1337 throw FinleyAdapterException("Error - Gradient at points is not supported.");
1338 break;
1339 case(ContactElementsZero):
1340 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1341 break;
1342 case(ReducedContactElementsZero):
1343 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1344 break;
1345 case(ContactElementsOne):
1346 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1347 break;
1348 case(ReducedContactElementsOne):
1349 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1350 break;
1351 case(DegreesOfFreedom):
1352 throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
1353 break;
1354 case(ReducedDegreesOfFreedom):
1355 throw FinleyAdapterException("Error - Gradient at reduced degrees of freedom is not supported.");
1356 break;
1357 default:
1358 stringstream temp;
1359 temp << "Error - Gradient: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
1360 throw FinleyAdapterException(temp.str());
1361 break;
1362 }
1363 checkFinleyError();
1364 }
1365
1366 //
1367 // returns the size of elements:
1368 //
1369 void MeshAdapter::setToSize(escript::Data& size) const
1370 {
1371 Finley_Mesh* mesh=m_finleyMesh.get();
1372 escriptDataC tmp=size.getDataC();
1373 switch(size.getFunctionSpace().getTypeCode()) {
1374 case(Nodes):
1375 throw FinleyAdapterException("Error - Size of nodes is not supported.");
1376 break;
1377 case(ReducedNodes):
1378 throw FinleyAdapterException("Error - Size of reduced nodes is not supported.");
1379 break;
1380 case(Elements):
1381 Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
1382 break;
1383 case(ReducedElements):
1384 Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
1385 break;
1386 case(FaceElements):
1387 Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
1388 break;
1389 case(ReducedFaceElements):
1390 Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
1391 break;
1392 case(Points):
1393 throw FinleyAdapterException("Error - Size of point elements is not supported.");
1394 break;
1395 case(ContactElementsZero):
1396 case(ContactElementsOne):
1397 Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
1398 break;
1399 case(ReducedContactElementsZero):
1400 case(ReducedContactElementsOne):
1401 Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
1402 break;
1403 case(DegreesOfFreedom):
1404 throw FinleyAdapterException("Error - Size of degrees of freedom is not supported.");
1405 break;
1406 case(ReducedDegreesOfFreedom):
1407 throw FinleyAdapterException("Error - Size of reduced degrees of freedom is not supported.");
1408 break;
1409 default:
1410 stringstream temp;
1411 temp << "Error - Element size: Finley does not know anything about function space type " << size.getFunctionSpace().getTypeCode();
1412 throw FinleyAdapterException(temp.str());
1413 break;
1414 }
1415 checkFinleyError();
1416 }
1417
1418 // sets the location of nodes:
1419 void MeshAdapter::setNewX(const escript::Data& new_x)
1420 {
1421 Finley_Mesh* mesh=m_finleyMesh.get();
1422 escriptDataC tmp;
1423 const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());
1424 if (newDomain!=*this)
1425 throw FinleyAdapterException("Error - Illegal domain of new point locations");
1426 tmp = new_x.getDataC();
1427 Finley_Mesh_setCoordinates(mesh,&tmp);
1428 checkFinleyError();
1429 }
1430
1431 // saves a data array in openDX format:
1432 void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const
1433 {
1434 unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1435 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1436 /* win32 refactor */
1437 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
1438 for(int i=0;i<num_data;i++)
1439 {
1440 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
1441 }
1442
1443 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
1444 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
1445
1446 boost::python::list keys=arg.keys();
1447 for (int i=0;i<num_data;++i) {
1448 std::string n=boost::python::extract<std::string>(keys[i]);
1449 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
1450 if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
1451 throw FinleyAdapterException("Error in saveVTK: Data must be defined on same Domain");
1452 data[i]=d.getDataC();
1453 ptr_data[i]= &(data[i]);
1454 if (n.length()>MAX_namelength-1) {
1455 strncpy(names[i],n.c_str(),MAX_namelength-1);
1456 names[i][MAX_namelength-1]='\0';
1457 } else {
1458 strcpy(names[i],n.c_str());
1459 }
1460 }
1461 Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,names,ptr_data);
1462 checkFinleyError();
1463
1464 /* win32 refactor */
1465 TMPMEMFREE(data);
1466 TMPMEMFREE(ptr_data);
1467 for(int i=0;i<num_data;i++)
1468 {
1469 TMPMEMFREE(names[i]);
1470 }
1471 TMPMEMFREE(names);
1472
1473 return;
1474 }
1475
1476 // saves a data array in openVTK format:
1477 void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const
1478 {
1479 unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1480 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1481 /* win32 refactor */
1482 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
1483 for(int i=0;i<num_data;i++)
1484 {
1485 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
1486 }
1487
1488 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
1489 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
1490
1491 boost::python::list keys=arg.keys();
1492 for (int i=0;i<num_data;++i) {
1493 std::string n=boost::python::extract<std::string>(keys[i]);
1494 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
1495 if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
1496 throw FinleyAdapterException("Error in saveVTK: Data must be defined on same Domain");
1497 data[i]=d.getDataC();
1498 ptr_data[i]=&(data[i]);
1499 if (n.length()>MAX_namelength-1) {
1500 strncpy(names[i],n.c_str(),MAX_namelength-1);
1501 names[i][MAX_namelength-1]='\0';
1502 } else {
1503 strcpy(names[i],n.c_str());
1504 }
1505 }
1506 Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,names,ptr_data);
1507
1508 checkFinleyError();
1509 /* win32 refactor */
1510 TMPMEMFREE(data);
1511 TMPMEMFREE(ptr_data);
1512 for(int i=0;i<num_data;i++)
1513 {
1514 TMPMEMFREE(names[i]);
1515 }
1516 TMPMEMFREE(names);
1517
1518 return;
1519 }
1520
1521
1522 // creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
1523 SystemMatrixAdapter MeshAdapter::newSystemMatrix(
1524 const int row_blocksize,
1525 const escript::FunctionSpace& row_functionspace,
1526 const int column_blocksize,
1527 const escript::FunctionSpace& column_functionspace,
1528 const int type) const
1529 {
1530 int reduceRowOrder=0;
1531 int reduceColOrder=0;
1532 // is the domain right?
1533 const MeshAdapter& row_domain=dynamic_cast<const MeshAdapter&>(row_functionspace.getDomain());
1534 if (row_domain!=*this)
1535 throw FinleyAdapterException("Error - domain of row function space does not match the domain of matrix generator.");
1536 const MeshAdapter& col_domain=dynamic_cast<const MeshAdapter&>(column_functionspace.getDomain());
1537 if (col_domain!=*this)
1538 throw FinleyAdapterException("Error - domain of columnn function space does not match the domain of matrix generator.");
1539 // is the function space type right
1540 if (row_functionspace.getTypeCode()==DegreesOfFreedom) {
1541 reduceRowOrder=0;
1542 } else if (row_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1543 reduceRowOrder=1;
1544 } else {
1545 throw FinleyAdapterException("Error - illegal function space type for system matrix rows.");
1546 }
1547 if (column_functionspace.getTypeCode()==DegreesOfFreedom) {
1548 reduceColOrder=0;
1549 } else if (column_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1550 reduceColOrder=1;
1551 } else {
1552 throw FinleyAdapterException("Error - illegal function space type for system matrix columns.");
1553 }
1554 // generate matrix:
1555
1556 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);
1557 checkFinleyError();
1558 Paso_SystemMatrix* fsystemMatrix;
1559 int trilinos = 0;
1560 if (trilinos) {
1561 #ifdef TRILINOS
1562 /* Allocation Epetra_VrbMatrix here */
1563 #endif
1564 }
1565 else {
1566 fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);
1567 }
1568 checkPasoError();
1569 Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1570 return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
1571 }
1572 // creates a TransportProblemAdapter
1573 TransportProblemAdapter MeshAdapter::newTransportProblem(
1574 const double theta,
1575 const int blocksize,
1576 const escript::FunctionSpace& functionspace,
1577 const int type) const
1578 {
1579 int reduceOrder=0;
1580 // is the domain right?
1581 const MeshAdapter& domain=dynamic_cast<const MeshAdapter&>(functionspace.getDomain());
1582 if (domain!=*this)
1583 throw FinleyAdapterException("Error - domain of function space does not match the domain of transport problem generator.");
1584 // is the function space type right
1585 if (functionspace.getTypeCode()==DegreesOfFreedom) {
1586 reduceOrder=0;
1587 } else if (functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1588 reduceOrder=1;
1589 } else {
1590 throw FinleyAdapterException("Error - illegal function space type for system matrix rows.");
1591 }
1592 // generate matrix:
1593
1594 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);
1595 checkFinleyError();
1596 Paso_FCTransportProblem* transportProblem;
1597 transportProblem=Paso_FCTransportProblem_alloc(theta,fsystemMatrixPattern,blocksize);
1598 checkPasoError();
1599 Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1600 return TransportProblemAdapter(transportProblem,theta,blocksize,functionspace);
1601 }
1602
1603 //
1604 // vtkObject MeshAdapter::createVtkObject() const
1605 // TODO:
1606 //
1607
1608 //
1609 // returns true if data at the atom_type is considered as being cell centered:
1610 bool MeshAdapter::isCellOriented(int functionSpaceCode) const
1611 {
1612 switch(functionSpaceCode) {
1613 case(Nodes):
1614 case(DegreesOfFreedom):
1615 case(ReducedDegreesOfFreedom):
1616 return false;
1617 break;
1618 case(Elements):
1619 case(FaceElements):
1620 case(Points):
1621 case(ContactElementsZero):
1622 case(ContactElementsOne):
1623 case(ReducedElements):
1624 case(ReducedFaceElements):
1625 case(ReducedContactElementsZero):
1626 case(ReducedContactElementsOne):
1627 return true;
1628 break;
1629 default:
1630 stringstream temp;
1631 temp << "Error - Cell: Finley does not know anything about function space type " << functionSpaceCode;
1632 throw FinleyAdapterException(temp.str());
1633 break;
1634 }
1635 checkFinleyError();
1636 return false;
1637 }
1638
1639 bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const
1640 {
1641 switch(functionSpaceType_source) {
1642 case(Nodes):
1643 switch(functionSpaceType_target) {
1644 case(Nodes):
1645 case(ReducedNodes):
1646 case(ReducedDegreesOfFreedom):
1647 case(DegreesOfFreedom):
1648 case(Elements):
1649 case(ReducedElements):
1650 case(FaceElements):
1651 case(ReducedFaceElements):
1652 case(Points):
1653 case(ContactElementsZero):
1654 case(ReducedContactElementsZero):
1655 case(ContactElementsOne):
1656 case(ReducedContactElementsOne):
1657 return true;
1658 default:
1659 stringstream temp;
1660 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1661 throw FinleyAdapterException(temp.str());
1662 }
1663 break;
1664 case(ReducedNodes):
1665 switch(functionSpaceType_target) {
1666 case(ReducedNodes):
1667 case(ReducedDegreesOfFreedom):
1668 case(Elements):
1669 case(ReducedElements):
1670 case(FaceElements):
1671 case(ReducedFaceElements):
1672 case(Points):
1673 case(ContactElementsZero):
1674 case(ReducedContactElementsZero):
1675 case(ContactElementsOne):
1676 case(ReducedContactElementsOne):
1677 return true;
1678 case(Nodes):
1679 case(DegreesOfFreedom):
1680 return false;
1681 default:
1682 stringstream temp;
1683 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1684 throw FinleyAdapterException(temp.str());
1685 }
1686 break;
1687 case(Elements):
1688 if (functionSpaceType_target==Elements) {
1689 return true;
1690 } else if (functionSpaceType_target==ReducedElements) {
1691 return true;
1692 } else {
1693 return false;
1694 }
1695 case(ReducedElements):
1696 if (functionSpaceType_target==ReducedElements) {
1697 return true;
1698 } else {
1699 return false;
1700 }
1701 case(FaceElements):
1702 if (functionSpaceType_target==FaceElements) {
1703 return true;
1704 } else if (functionSpaceType_target==ReducedFaceElements) {
1705 return true;
1706 } else {
1707 return false;
1708 }
1709 case(ReducedFaceElements):
1710 if (functionSpaceType_target==ReducedFaceElements) {
1711 return true;
1712 } else {
1713 return false;
1714 }
1715 case(Points):
1716 if (functionSpaceType_target==Points) {
1717 return true;
1718 } else {
1719 return false;
1720 }
1721 case(ContactElementsZero):
1722 case(ContactElementsOne):
1723 if (functionSpaceType_target==ContactElementsZero || functionSpaceType_target==ContactElementsOne) {
1724 return true;
1725 } else if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1726 return true;
1727 } else {
1728 return false;
1729 }
1730 case(ReducedContactElementsZero):
1731 case(ReducedContactElementsOne):
1732 if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1733 return true;
1734 } else {
1735 return false;
1736 }
1737 case(DegreesOfFreedom):
1738 switch(functionSpaceType_target) {
1739 case(ReducedDegreesOfFreedom):
1740 case(DegreesOfFreedom):
1741 case(Nodes):
1742 case(ReducedNodes):
1743 case(Elements):
1744 case(ReducedElements):
1745 case(Points):
1746 case(FaceElements):
1747 case(ReducedFaceElements):
1748 case(ContactElementsZero):
1749 case(ReducedContactElementsZero):
1750 case(ContactElementsOne):
1751 case(ReducedContactElementsOne):
1752 return true;
1753 default:
1754 stringstream temp;
1755 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1756 throw FinleyAdapterException(temp.str());
1757 }
1758 break;
1759 case(ReducedDegreesOfFreedom):
1760 switch(functionSpaceType_target) {
1761 case(ReducedDegreesOfFreedom):
1762 case(ReducedNodes):
1763 case(Elements):
1764 case(ReducedElements):
1765 case(FaceElements):
1766 case(ReducedFaceElements):
1767 case(Points):
1768 case(ContactElementsZero):
1769 case(ReducedContactElementsZero):
1770 case(ContactElementsOne):
1771 case(ReducedContactElementsOne):
1772 return true;
1773 case(Nodes):
1774 case(DegreesOfFreedom):
1775 return false;
1776 default:
1777 stringstream temp;
1778 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1779 throw FinleyAdapterException(temp.str());
1780 }
1781 break;
1782 default:
1783 stringstream temp;
1784 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;
1785 throw FinleyAdapterException(temp.str());
1786 break;
1787 }
1788 checkFinleyError();
1789 return false;
1790 }
1791
1792 bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const
1793 {
1794 return false;
1795 }
1796
1797 bool MeshAdapter::operator==(const AbstractDomain& other) const
1798 {
1799 const MeshAdapter* temp=dynamic_cast<const MeshAdapter*>(&other);
1800 if (temp!=0) {
1801 return (m_finleyMesh==temp->m_finleyMesh);
1802 } else {
1803 return false;
1804 }
1805 }
1806
1807 bool MeshAdapter::operator!=(const AbstractDomain& other) const
1808 {
1809 return !(operator==(other));
1810 }
1811
1812 int MeshAdapter::getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const
1813 {
1814 int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
1815 checkPasoError();
1816 return out;
1817 }
1818
1819 escript::Data MeshAdapter::getX() const
1820 {
1821 return continuousFunction(asAbstractContinuousDomain()).getX();
1822 }
1823
1824 escript::Data MeshAdapter::getNormal() const
1825 {
1826 return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
1827 }
1828
1829 escript::Data MeshAdapter::getSize() const
1830 {
1831 return function(asAbstractContinuousDomain()).getSize();
1832 }
1833
1834 int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1835 {
1836 int *out = NULL;
1837 Finley_Mesh* mesh=m_finleyMesh.get();
1838 switch (functionSpaceType) {
1839 case(Nodes):
1840 out=mesh->Nodes->Id;
1841 break;
1842 case(ReducedNodes):
1843 out=mesh->Nodes->reducedNodesId;
1844 break;
1845 case(Elements):
1846 out=mesh->Elements->Id;
1847 break;
1848 case(ReducedElements):
1849 out=mesh->Elements->Id;
1850 break;
1851 case(FaceElements):
1852 out=mesh->FaceElements->Id;
1853 break;
1854 case(ReducedFaceElements):
1855 out=mesh->FaceElements->Id;
1856 break;
1857 case(Points):
1858 out=mesh->Points->Id;
1859 break;
1860 case(ContactElementsZero):
1861 out=mesh->ContactElements->Id;
1862 break;
1863 case(ReducedContactElementsZero):
1864 out=mesh->ContactElements->Id;
1865 break;
1866 case(ContactElementsOne):
1867 out=mesh->ContactElements->Id;
1868 break;
1869 case(ReducedContactElementsOne):
1870 out=mesh->ContactElements->Id;
1871 break;
1872 case(DegreesOfFreedom):
1873 out=mesh->Nodes->degreesOfFreedomId;
1874 break;
1875 case(ReducedDegreesOfFreedom):
1876 out=mesh->Nodes->reducedDegreesOfFreedomId;
1877 break;
1878 default:
1879 stringstream temp;
1880 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1881 throw FinleyAdapterException(temp.str());
1882 break;
1883 }
1884 return out;
1885 }
1886 int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
1887 {
1888 int out=0;
1889 Finley_Mesh* mesh=m_finleyMesh.get();
1890 switch (functionSpaceType) {
1891 case(Nodes):
1892 out=mesh->Nodes->Tag[sampleNo];
1893 break;
1894 case(ReducedNodes):
1895 throw FinleyAdapterException(" Error - ReducedNodes does not support tags.");
1896 break;
1897 case(Elements):
1898 out=mesh->Elements->Tag[sampleNo];
1899 break;
1900 case(ReducedElements):
1901 out=mesh->Elements->Tag[sampleNo];
1902 break;
1903 case(FaceElements):
1904 out=mesh->FaceElements->Tag[sampleNo];
1905 break;
1906 case(ReducedFaceElements):
1907 out=mesh->FaceElements->Tag[sampleNo];
1908 break;
1909 case(Points):
1910 out=mesh->Points->Tag[sampleNo];
1911 break;
1912 case(ContactElementsZero):
1913 out=mesh->ContactElements->Tag[sampleNo];
1914 break;
1915 case(ReducedContactElementsZero):
1916 out=mesh->ContactElements->Tag[sampleNo];
1917 break;
1918 case(ContactElementsOne):
1919 out=mesh->ContactElements->Tag[sampleNo];
1920 break;
1921 case(ReducedContactElementsOne):
1922 out=mesh->ContactElements->Tag[sampleNo];
1923 break;
1924 case(DegreesOfFreedom):
1925 throw FinleyAdapterException(" Error - DegreesOfFreedom does not support tags.");
1926 break;
1927 case(ReducedDegreesOfFreedom):
1928 throw FinleyAdapterException(" Error - ReducedDegreesOfFreedom does not support tags.");
1929 break;
1930 default:
1931 stringstream temp;
1932 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1933 throw FinleyAdapterException(temp.str());
1934 break;
1935 }
1936 return out;
1937 }
1938
1939
1940 void MeshAdapter::setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const
1941 {
1942 Finley_Mesh* mesh=m_finleyMesh.get();
1943 escriptDataC tmp=mask.getDataC();
1944 switch(functionSpaceType) {
1945 case(Nodes):
1946 Finley_NodeFile_setTags(mesh->Nodes,newTag,&tmp);
1947 break;
1948 case(ReducedNodes):
1949 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
1950 break;
1951 case(DegreesOfFreedom):
1952 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
1953 break;
1954 case(ReducedDegreesOfFreedom):
1955 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
1956 break;
1957 case(Elements):
1958 Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1959 break;
1960 case(ReducedElements):
1961 Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1962 break;
1963 case(FaceElements):
1964 Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
1965 break;
1966 case(ReducedFaceElements):
1967 Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
1968 break;
1969 case(Points):
1970 Finley_ElementFile_setTags(mesh->Points,newTag,&tmp);
1971 break;
1972 case(ContactElementsZero):
1973 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1974 break;
1975 case(ReducedContactElementsZero):
1976 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1977 break;
1978 case(ContactElementsOne):
1979 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1980 break;
1981 case(ReducedContactElementsOne):
1982 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1983 break;
1984 default:
1985 stringstream temp;
1986 temp << "Error - Finley does not know anything about function space type " << functionSpaceType;
1987 throw FinleyAdapterException(temp.str());
1988 }
1989 checkFinleyError();
1990 return;
1991 }
1992
1993 void MeshAdapter::setTagMap(const std::string& name, int tag)
1994 {
1995 Finley_Mesh* mesh=m_finleyMesh.get();
1996 Finley_Mesh_addTagMap(mesh, name.c_str(),tag);
1997 checkPasoError();
1998 // throwStandardException("MeshAdapter::set TagMap is not implemented.");
1999 }
2000
2001 int MeshAdapter::getTag(const std::string& name) const
2002 {
2003 Finley_Mesh* mesh=m_finleyMesh.get();
2004 int tag=0;
2005 tag=Finley_Mesh_getTag(mesh, name.c_str());
2006 checkPasoError();
2007 // throwStandardException("MeshAdapter::getTag is not implemented.");
2008 return tag;
2009 }
2010
2011 bool MeshAdapter::isValidTagName(const std::string& name) const
2012 {
2013 Finley_Mesh* mesh=m_finleyMesh.get();
2014 return Finley_Mesh_isValidTagName(mesh,name.c_str());
2015 }
2016
2017 std::string MeshAdapter::showTagNames() const
2018 {
2019 stringstream temp;
2020 Finley_Mesh* mesh=m_finleyMesh.get();
2021 Finley_TagMap* tag_map=mesh->TagMap;
2022 while (tag_map) {
2023 temp << tag_map->name;
2024 tag_map=tag_map->next;
2025 if (tag_map) temp << ", ";
2026 }
2027 return temp.str();
2028 }
2029
2030 int MeshAdapter::getNumberOfTagsInUse(int functionSpaceCode) const
2031 {
2032 Finley_Mesh* mesh=m_finleyMesh.get();
2033 dim_t numTags=0;
2034 switch(functionSpaceCode) {
2035 case(Nodes):
2036 numTags=mesh->Nodes->numTagsInUse;
2037 break;
2038 case(ReducedNodes):
2039 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
2040 break;
2041 case(DegreesOfFreedom):
2042 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
2043 break;
2044 case(ReducedDegreesOfFreedom):
2045 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
2046 break;
2047 case(Elements):
2048 case(ReducedElements):
2049 numTags=mesh->Elements->numTagsInUse;
2050 break;
2051 case(FaceElements):
2052 case(ReducedFaceElements):
2053 numTags=mesh->FaceElements->numTagsInUse;
2054 break;
2055 case(Points):
2056 numTags=mesh->Points->numTagsInUse;
2057 break;
2058 case(ContactElementsZero):
2059 case(ReducedContactElementsZero):
2060 case(ContactElementsOne):
2061 case(ReducedContactElementsOne):
2062 numTags=mesh->ContactElements->numTagsInUse;
2063 break;
2064 default:
2065 stringstream temp;
2066 temp << "Error - Finley does not know anything about function space type " << functionSpaceCode;
2067 throw FinleyAdapterException(temp.str());
2068 }
2069 return numTags;
2070 }
2071 int* MeshAdapter::borrowListOfTagsInUse(int functionSpaceCode) const
2072 {
2073 Finley_Mesh* mesh=m_finleyMesh.get();
2074 index_t* tags=NULL;
2075 switch(functionSpaceCode) {
2076 case(Nodes):
2077 tags=mesh->Nodes->tagsInUse;
2078 break;
2079 case(ReducedNodes):
2080 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
2081 break;
2082 case(DegreesOfFreedom):
2083 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
2084 break;
2085 case(ReducedDegreesOfFreedom):
2086 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
2087 break;
2088 case(Elements):
2089 case(ReducedElements):
2090 tags=mesh->Elements->tagsInUse;
2091 break;
2092 case(FaceElements):
2093 case(ReducedFaceElements):
2094 tags=mesh->FaceElements->tagsInUse;
2095 break;
2096 case(Points):
2097 tags=mesh->Points->tagsInUse;
2098 break;
2099 case(ContactElementsZero):
2100 case(ReducedContactElementsZero):
2101 case(ContactElementsOne):
2102 case(ReducedContactElementsOne):
2103 tags=mesh->ContactElements->tagsInUse;
2104 break;
2105 default:
2106 stringstream temp;
2107 temp << "Error - Finley does not know anything about function space type " << functionSpaceCode;
2108 throw FinleyAdapterException(temp.str());
2109 }
2110 return tags;
2111 }
2112
2113
2114 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26