/[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 2149 - (show annotations)
Wed Dec 10 05:08:17 2008 UTC (10 years, 4 months ago) by caltinay
File size: 85399 byte(s)
Fixed a typo in error message.

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26