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

Contents of /branches/diaplayground/finley/src/CPPAdapter/MeshAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1802 - (show annotations)
Tue Sep 23 01:03:29 2008 UTC (10 years, 6 months ago) by jfenwick
Original Path: trunk/finley/src/CPPAdapter/MeshAdapter.cpp
File size: 84417 byte(s)
Added canTag methods to FunctionSpace and AbstractDomain (and its 
offspring).
This checks to see if the domain supports tags for the given type of 
function space.

Constructors for DataTagged now throw exceptions if you attempt to make 
a DataTagged with a FunctionSpace which does not support tags.
To allow the default constructor to work, NullDomain has a single 
functioncode which "supports" tagging.

Fixed a bug in DataTagged::toString and DataTypes::pointToString.

Added FunctionSpace::getListOfTagsSTL.

algorithm(DataTagged, BinaryFunction) in DataAlgorithm now only 
processes tags known to be in use.
This fixes mantis issue #0000186.

Added comment to Data.h intro warning about holding references if the 
underlying DataAbstract changes.

_python_ unit tests have been updated to test TaggedData with invalid 
FunctionSpaces and to give the correct answers to Lsup etc.


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26