/[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 2087 - (show annotations)
Mon Nov 24 04:51:30 2008 UTC (10 years, 4 months ago) by jfenwick
Original Path: trunk/finley/src/CPPAdapter/MeshAdapter.cpp
File size: 85375 byte(s)
Addressing mantis issue #221.
Interpolation.. and probeInterpolation.. now "work" for the NullDomain.
Work means throw a descriptive exception if you try to move into or out 
of the NullDomain.
The bad_cast exception related to this has been fixed.


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 int numDim=Finley_Mesh_getDim(m_finleyMesh.get());
650 checkFinleyError();
651 return numDim;
652 }
653
654 //
655 // Return the number of data points summed across all MPI processes
656 //
657 int MeshAdapter::getNumDataPointsGlobal() const
658 {
659 return Finley_NodeFile_getGlobalNumNodes(m_finleyMesh.get()->Nodes);
660 }
661
662 //
663 // return the number of data points per sample and the number of samples
664 // needed to represent data on a parts of the mesh.
665 //
666 pair<int,int> MeshAdapter::getDataShape(int functionSpaceCode) const
667 {
668 int numDataPointsPerSample=0;
669 int numSamples=0;
670 Finley_Mesh* mesh=m_finleyMesh.get();
671 switch (functionSpaceCode) {
672 case(Nodes):
673 numDataPointsPerSample=1;
674 numSamples=Finley_NodeFile_getNumNodes(mesh->Nodes);
675 break;
676 case(ReducedNodes):
677 numDataPointsPerSample=1;
678 numSamples=Finley_NodeFile_getNumReducedNodes(mesh->Nodes);
679 break;
680 case(Elements):
681 if (mesh->Elements!=NULL) {
682 numSamples=mesh->Elements->numElements;
683 numDataPointsPerSample=mesh->Elements->ReferenceElement->numQuadNodes;
684 }
685 break;
686 case(ReducedElements):
687 if (mesh->Elements!=NULL) {
688 numSamples=mesh->Elements->numElements;
689 numDataPointsPerSample=mesh->Elements->ReferenceElementReducedOrder->numQuadNodes;
690 }
691 break;
692 case(FaceElements):
693 if (mesh->FaceElements!=NULL) {
694 numDataPointsPerSample=mesh->FaceElements->ReferenceElement->numQuadNodes;
695 numSamples=mesh->FaceElements->numElements;
696 }
697 break;
698 case(ReducedFaceElements):
699 if (mesh->FaceElements!=NULL) {
700 numDataPointsPerSample=mesh->FaceElements->ReferenceElementReducedOrder->numQuadNodes;
701 numSamples=mesh->FaceElements->numElements;
702 }
703 break;
704 case(Points):
705 if (mesh->Points!=NULL) {
706 numDataPointsPerSample=1;
707 numSamples=mesh->Points->numElements;
708 }
709 break;
710 case(ContactElementsZero):
711 if (mesh->ContactElements!=NULL) {
712 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
713 numSamples=mesh->ContactElements->numElements;
714 }
715 break;
716 case(ReducedContactElementsZero):
717 if (mesh->ContactElements!=NULL) {
718 numDataPointsPerSample=mesh->ContactElements->ReferenceElementReducedOrder->numQuadNodes;
719 numSamples=mesh->ContactElements->numElements;
720 }
721 break;
722 case(ContactElementsOne):
723 if (mesh->ContactElements!=NULL) {
724 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
725 numSamples=mesh->ContactElements->numElements;
726 }
727 break;
728 case(ReducedContactElementsOne):
729 if (mesh->ContactElements!=NULL) {
730 numDataPointsPerSample=mesh->ContactElements->ReferenceElementReducedOrder->numQuadNodes;
731 numSamples=mesh->ContactElements->numElements;
732 }
733 break;
734 case(DegreesOfFreedom):
735 if (mesh->Nodes!=NULL) {
736 numDataPointsPerSample=1;
737 numSamples=Finley_NodeFile_getNumDegreesOfFreedom(mesh->Nodes);
738 }
739 break;
740 case(ReducedDegreesOfFreedom):
741 if (mesh->Nodes!=NULL) {
742 numDataPointsPerSample=1;
743 numSamples=Finley_NodeFile_getNumReducedDegreesOfFreedom(mesh->Nodes);
744 }
745 break;
746 default:
747 stringstream temp;
748 temp << "Error - Invalid function space type: " << functionSpaceCode << " for domain: " << getDescription();
749 throw FinleyAdapterException(temp.str());
750 break;
751 }
752 return pair<int,int>(numDataPointsPerSample,numSamples);
753 }
754
755 //
756 // adds linear PDE of second order into a given stiffness matrix and right hand side:
757 //
758 void MeshAdapter::addPDEToSystem(
759 SystemMatrixAdapter& mat, escript::Data& rhs,
760 const escript::Data& A, const escript::Data& B, const escript::Data& C,const escript::Data& D,const escript::Data& X,const escript::Data& Y,
761 const escript::Data& d, const escript::Data& y,
762 const escript::Data& d_contact,const escript::Data& y_contact) const
763 {
764 escriptDataC _rhs=rhs.getDataC();
765 escriptDataC _A =A.getDataC();
766 escriptDataC _B=B.getDataC();
767 escriptDataC _C=C.getDataC();
768 escriptDataC _D=D.getDataC();
769 escriptDataC _X=X.getDataC();
770 escriptDataC _Y=Y.getDataC();
771 escriptDataC _d=d.getDataC();
772 escriptDataC _y=y.getDataC();
773 escriptDataC _d_contact=d_contact.getDataC();
774 escriptDataC _y_contact=y_contact.getDataC();
775
776 Finley_Mesh* mesh=m_finleyMesh.get();
777
778 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(), &_rhs, &_A, &_B, &_C, &_D, &_X, &_Y );
779 checkFinleyError();
780
781 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, mat.getPaso_SystemMatrix(), &_rhs, 0, 0, 0, &_d, 0, &_y );
782 checkFinleyError();
783
784 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, mat.getPaso_SystemMatrix(), &_rhs , 0, 0, 0, &_d_contact, 0, &_y_contact );
785 checkFinleyError();
786 }
787
788 void MeshAdapter::addPDEToLumpedSystem(
789 escript::Data& mat,
790 const escript::Data& D,
791 const escript::Data& d) const
792 {
793 escriptDataC _mat=mat.getDataC();
794 escriptDataC _D=D.getDataC();
795 escriptDataC _d=d.getDataC();
796
797 Finley_Mesh* mesh=m_finleyMesh.get();
798
799 Finley_Assemble_LumpedSystem(mesh->Nodes,mesh->Elements,&_mat, &_D);
800 Finley_Assemble_LumpedSystem(mesh->Nodes,mesh->FaceElements,&_mat, &_d);
801
802 checkFinleyError();
803 }
804
805
806 //
807 // adds linear PDE of second order into the right hand side only
808 //
809 void MeshAdapter::addPDEToRHS( escript::Data& rhs, const escript::Data& X,const escript::Data& Y, const escript::Data& y, const escript::Data& y_contact) const
810 {
811 Finley_Mesh* mesh=m_finleyMesh.get();
812
813 escriptDataC _rhs=rhs.getDataC();
814 escriptDataC _X=X.getDataC();
815 escriptDataC _Y=Y.getDataC();
816 escriptDataC _y=y.getDataC();
817 escriptDataC _y_contact=y_contact.getDataC();
818
819 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements, 0, &_rhs, 0, 0, 0, 0, &_X, &_Y );
820 checkFinleyError();
821
822 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, 0, &_rhs, 0, 0, 0, 0, 0, &_y );
823 checkFinleyError();
824
825 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, 0, &_rhs , 0, 0, 0, 0, 0, &_y_contact );
826 checkFinleyError();
827 }
828 //
829 // adds PDE of second order into a transport problem
830 //
831 void MeshAdapter::addPDEToTransportProblem(
832 TransportProblemAdapter& tp, escript::Data& source, const escript::Data& M,
833 const escript::Data& A, const escript::Data& B, const escript::Data& C,
834 const escript::Data& D,const escript::Data& X,const escript::Data& Y,
835 const escript::Data& d, const escript::Data& y,
836 const escript::Data& d_contact,const escript::Data& y_contact) const
837 {
838 DataTypes::ShapeType shape;
839 source.expand();
840 escriptDataC _source=source.getDataC();
841 escriptDataC _M=M.getDataC();
842 escriptDataC _A=A.getDataC();
843 escriptDataC _B=B.getDataC();
844 escriptDataC _C=C.getDataC();
845 escriptDataC _D=D.getDataC();
846 escriptDataC _X=X.getDataC();
847 escriptDataC _Y=Y.getDataC();
848 escriptDataC _d=d.getDataC();
849 escriptDataC _y=y.getDataC();
850 escriptDataC _d_contact=d_contact.getDataC();
851 escriptDataC _y_contact=y_contact.getDataC();
852
853 Finley_Mesh* mesh=m_finleyMesh.get();
854 Paso_FCTransportProblem* _tp = tp.getPaso_FCTransportProblem();
855
856 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->mass_matrix, &_source, 0, 0, 0, &_M, 0, 0 );
857 checkFinleyError();
858
859 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->transport_matrix, &_source, &_A, &_B, &_C, &_D, &_X, &_Y );
860 checkFinleyError();
861
862 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, _tp->transport_matrix, &_source, 0, 0, 0, &_d, 0, &_y );
863 checkFinleyError();
864
865 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, _tp->transport_matrix, &_source , 0, 0, 0, &_d_contact, 0, &_y_contact );
866 checkFinleyError();
867 }
868
869 //
870 // interpolates data between different function spaces:
871 //
872 void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
873 {
874 const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(*(in.getFunctionSpace().getDomain()));
875 const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(*(target.getFunctionSpace().getDomain()));
876 if (inDomain!=*this)
877 throw FinleyAdapterException("Error - Illegal domain of interpolant.");
878 if (targetDomain!=*this)
879 throw FinleyAdapterException("Error - Illegal domain of interpolation target.");
880
881 Finley_Mesh* mesh=m_finleyMesh.get();
882 escriptDataC _target=target.getDataC();
883 escriptDataC _in=in.getDataC();
884 switch(in.getFunctionSpace().getTypeCode()) {
885 case(Nodes):
886 switch(target.getFunctionSpace().getTypeCode()) {
887 case(Nodes):
888 case(ReducedNodes):
889 case(DegreesOfFreedom):
890 case(ReducedDegreesOfFreedom):
891 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
892 break;
893 case(Elements):
894 case(ReducedElements):
895 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
896 break;
897 case(FaceElements):
898 case(ReducedFaceElements):
899 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
900 break;
901 case(Points):
902 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
903 break;
904 case(ContactElementsZero):
905 case(ReducedContactElementsZero):
906 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
907 break;
908 case(ContactElementsOne):
909 case(ReducedContactElementsOne):
910 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
911 break;
912 default:
913 stringstream temp;
914 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
915 throw FinleyAdapterException(temp.str());
916 break;
917 }
918 break;
919 case(ReducedNodes):
920 switch(target.getFunctionSpace().getTypeCode()) {
921 case(Nodes):
922 case(ReducedNodes):
923 case(DegreesOfFreedom):
924 case(ReducedDegreesOfFreedom):
925 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
926 break;
927 case(Elements):
928 case(ReducedElements):
929 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
930 break;
931 case(FaceElements):
932 case(ReducedFaceElements):
933 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
934 break;
935 case(Points):
936 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
937 break;
938 case(ContactElementsZero):
939 case(ReducedContactElementsZero):
940 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
941 break;
942 case(ContactElementsOne):
943 case(ReducedContactElementsOne):
944 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
945 break;
946 default:
947 stringstream temp;
948 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
949 throw FinleyAdapterException(temp.str());
950 break;
951 }
952 break;
953 case(Elements):
954 if (target.getFunctionSpace().getTypeCode()==Elements) {
955 Finley_Assemble_CopyElementData(mesh->Elements,&_target,&_in);
956 } else if (target.getFunctionSpace().getTypeCode()==ReducedElements) {
957 Finley_Assemble_AverageElementData(mesh->Elements,&_target,&_in);
958 } else {
959 throw FinleyAdapterException("Error - No interpolation with data on elements possible.");
960 }
961 break;
962 case(ReducedElements):
963 if (target.getFunctionSpace().getTypeCode()==ReducedElements) {
964 Finley_Assemble_CopyElementData(mesh->Elements,&_target,&_in);
965 } else {
966 throw FinleyAdapterException("Error - No interpolation with data on elements with reduced integration order possible.");
967 }
968 break;
969 case(FaceElements):
970 if (target.getFunctionSpace().getTypeCode()==FaceElements) {
971 Finley_Assemble_CopyElementData(mesh->FaceElements,&_target,&_in);
972 } else if (target.getFunctionSpace().getTypeCode()==ReducedFaceElements) {
973 Finley_Assemble_AverageElementData(mesh->FaceElements,&_target,&_in);
974 } else {
975 throw FinleyAdapterException("Error - No interpolation with data on face elements possible.");
976 }
977 break;
978 case(ReducedFaceElements):
979 if (target.getFunctionSpace().getTypeCode()==ReducedFaceElements) {
980 Finley_Assemble_CopyElementData(mesh->FaceElements,&_target,&_in);
981 } else {
982 throw FinleyAdapterException("Error - No interpolation with data on face elements with reduced integration order possible.");
983 }
984 break;
985 case(Points):
986 if (target.getFunctionSpace().getTypeCode()==Points) {
987 Finley_Assemble_CopyElementData(mesh->Points,&_target,&_in);
988 } else {
989 throw FinleyAdapterException("Error - No interpolation with data on points possible.");
990 }
991 break;
992 case(ContactElementsZero):
993 case(ContactElementsOne):
994 if (target.getFunctionSpace().getTypeCode()==ContactElementsZero || target.getFunctionSpace().getTypeCode()==ContactElementsOne) {
995 Finley_Assemble_CopyElementData(mesh->ContactElements,&_target,&_in);
996 } else if (target.getFunctionSpace().getTypeCode()==ReducedContactElementsZero || target.getFunctionSpace().getTypeCode()==ReducedContactElementsOne) {
997 Finley_Assemble_AverageElementData(mesh->ContactElements,&_target,&_in);
998 } else {
999 throw FinleyAdapterException("Error - No interpolation with data on contact elements possible.");
1000 }
1001 break;
1002 case(ReducedContactElementsZero):
1003 case(ReducedContactElementsOne):
1004 if (target.getFunctionSpace().getTypeCode()==ReducedContactElementsZero || target.getFunctionSpace().getTypeCode()==ReducedContactElementsOne) {
1005 Finley_Assemble_CopyElementData(mesh->ContactElements,&_target,&_in);
1006 } else {
1007 throw FinleyAdapterException("Error - No interpolation with data on contact elements with reduced integration order possible.");
1008 }
1009 break;
1010 case(DegreesOfFreedom):
1011 switch(target.getFunctionSpace().getTypeCode()) {
1012 case(ReducedDegreesOfFreedom):
1013 case(DegreesOfFreedom):
1014 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1015 break;
1016
1017 case(Nodes):
1018 case(ReducedNodes):
1019 if (getMPISize()>1) {
1020 escript::Data temp=escript::Data(in);
1021 temp.expand();
1022 escriptDataC _in2 = temp.getDataC();
1023 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in2);
1024 } else {
1025 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1026 }
1027 break;
1028 case(Elements):
1029 case(ReducedElements):
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->Elements,&_in2,&_target);
1034 } else {
1035 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
1036 }
1037 break;
1038 case(FaceElements):
1039 case(ReducedFaceElements):
1040 if (getMPISize()>1) {
1041 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1042 escriptDataC _in2 = temp.getDataC();
1043 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1044
1045 } else {
1046 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
1047 }
1048 break;
1049 case(Points):
1050 if (getMPISize()>1) {
1051 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1052 escriptDataC _in2 = temp.getDataC();
1053 } else {
1054 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
1055 }
1056 break;
1057 case(ContactElementsZero):
1058 case(ContactElementsOne):
1059 case(ReducedContactElementsZero):
1060 case(ReducedContactElementsOne):
1061 if (getMPISize()>1) {
1062 escript::Data temp=escript::Data( in, continuousFunction(asAbstractContinuousDomain()) );
1063 escriptDataC _in2 = temp.getDataC();
1064 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1065 } else {
1066 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
1067 }
1068 break;
1069 default:
1070 stringstream temp;
1071 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
1072 throw FinleyAdapterException(temp.str());
1073 break;
1074 }
1075 break;
1076 case(ReducedDegreesOfFreedom):
1077 switch(target.getFunctionSpace().getTypeCode()) {
1078 case(Nodes):
1079 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to mesh nodes.");
1080 break;
1081 case(ReducedNodes):
1082 if (getMPISize()>1) {
1083 escript::Data temp=escript::Data(in);
1084 temp.expand();
1085 escriptDataC _in2 = temp.getDataC();
1086 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in2);
1087 } else {
1088 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1089 }
1090 break;
1091 case(DegreesOfFreedom):
1092 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");
1093 break;
1094 case(ReducedDegreesOfFreedom):
1095 Finley_Assemble_CopyNodalData(mesh->Nodes,&_target,&_in);
1096 break;
1097 case(Elements):
1098 case(ReducedElements):
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->Elements,&_in2,&_target);
1103 } else {
1104 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in,&_target);
1105 }
1106 break;
1107 case(FaceElements):
1108 case(ReducedFaceElements):
1109 if (getMPISize()>1) {
1110 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1111 escriptDataC _in2 = temp.getDataC();
1112 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1113 } else {
1114 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in,&_target);
1115 }
1116 break;
1117 case(Points):
1118 if (getMPISize()>1) {
1119 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1120 escriptDataC _in2 = temp.getDataC();
1121 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in2,&_target);
1122 } else {
1123 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
1124 }
1125 break;
1126 case(ContactElementsZero):
1127 case(ContactElementsOne):
1128 case(ReducedContactElementsZero):
1129 case(ReducedContactElementsOne):
1130 if (getMPISize()>1) {
1131 escript::Data temp=escript::Data( in, reducedContinuousFunction(asAbstractContinuousDomain()) );
1132 escriptDataC _in2 = temp.getDataC();
1133 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1134 } else {
1135 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in,&_target);
1136 }
1137 break;
1138 default:
1139 stringstream temp;
1140 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
1141 throw FinleyAdapterException(temp.str());
1142 break;
1143 }
1144 break;
1145 default:
1146 stringstream temp;
1147 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type %d" << in.getFunctionSpace().getTypeCode();
1148 throw FinleyAdapterException(temp.str());
1149 break;
1150 }
1151 checkFinleyError();
1152 }
1153
1154 //
1155 // copies the locations of sample points into x:
1156 //
1157 void MeshAdapter::setToX(escript::Data& arg) const
1158 {
1159 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(*(arg.getFunctionSpace().getDomain()));
1160 if (argDomain!=*this)
1161 throw FinleyAdapterException("Error - Illegal domain of data point locations");
1162 Finley_Mesh* mesh=m_finleyMesh.get();
1163 // in case of values node coordinates we can do the job directly:
1164 if (arg.getFunctionSpace().getTypeCode()==Nodes) {
1165 escriptDataC _arg=arg.getDataC();
1166 Finley_Assemble_NodeCoordinates(mesh->Nodes,&_arg);
1167 } else {
1168 escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
1169 escriptDataC _tmp_data=tmp_data.getDataC();
1170 Finley_Assemble_NodeCoordinates(mesh->Nodes,&_tmp_data);
1171 // this is then interpolated onto arg:
1172 interpolateOnDomain(arg,tmp_data);
1173 }
1174 checkFinleyError();
1175 }
1176
1177 //
1178 // return the normal vectors at the location of data points as a Data object:
1179 //
1180 void MeshAdapter::setToNormal(escript::Data& normal) const
1181 {
1182 /* const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());*/
1183 const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(*(normal.getFunctionSpace().getDomain()));
1184 if (normalDomain!=*this)
1185 throw FinleyAdapterException("Error - Illegal domain of normal locations");
1186 Finley_Mesh* mesh=m_finleyMesh.get();
1187 escriptDataC _normal=normal.getDataC();
1188 switch(normal.getFunctionSpace().getTypeCode()) {
1189 case(Nodes):
1190 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for nodes");
1191 break;
1192 case(ReducedNodes):
1193 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced nodes");
1194 break;
1195 case(Elements):
1196 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements");
1197 break;
1198 case(ReducedElements):
1199 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements with reduced integration order");
1200 break;
1201 case (FaceElements):
1202 Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&_normal);
1203 break;
1204 case (ReducedFaceElements):
1205 Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&_normal);
1206 break;
1207 case(Points):
1208 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for point elements");
1209 break;
1210 case (ContactElementsOne):
1211 case (ContactElementsZero):
1212 Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&_normal);
1213 break;
1214 case (ReducedContactElementsOne):
1215 case (ReducedContactElementsZero):
1216 Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&_normal);
1217 break;
1218 case(DegreesOfFreedom):
1219 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for degrees of freedom.");
1220 break;
1221 case(ReducedDegreesOfFreedom):
1222 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced degrees of freedom.");
1223 break;
1224 default:
1225 stringstream temp;
1226 temp << "Error - Normal Vectors: Finley does not know anything about function space type " << normal.getFunctionSpace().getTypeCode();
1227 throw FinleyAdapterException(temp.str());
1228 break;
1229 }
1230 checkFinleyError();
1231 }
1232
1233 //
1234 // interpolates data to other domain:
1235 //
1236 void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
1237 {
1238 const_Domain_ptr targetDomain_p=target.getFunctionSpace().getDomain();
1239 const MeshAdapter* targetDomain=dynamic_cast<const MeshAdapter*>(targetDomain_p.get());
1240 if (targetDomain!=this)
1241 throw FinleyAdapterException("Error - Illegal domain of interpolation target");
1242
1243 throw FinleyAdapterException("Error - Finley does not allow interpolation across domains yet.");
1244 }
1245
1246 //
1247 // calculates the integral of a function defined of arg:
1248 //
1249 void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
1250 {
1251 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(*(arg.getFunctionSpace().getDomain()));
1252 if (argDomain!=*this)
1253 throw FinleyAdapterException("Error - Illegal domain of integration kernel");
1254
1255 double blocktimer_start = blocktimer_time();
1256 Finley_Mesh* mesh=m_finleyMesh.get();
1257 escriptDataC _temp;
1258 escript::Data temp;
1259 escriptDataC _arg=arg.getDataC();
1260 switch(arg.getFunctionSpace().getTypeCode()) {
1261 case(Nodes):
1262 temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
1263 _temp=temp.getDataC();
1264 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1265 break;
1266 case(ReducedNodes):
1267 temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
1268 _temp=temp.getDataC();
1269 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1270 break;
1271 case(Elements):
1272 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_arg,&integrals[0]);
1273 break;
1274 case(ReducedElements):
1275 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_arg,&integrals[0]);
1276 break;
1277 case(FaceElements):
1278 Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&_arg,&integrals[0]);
1279 break;
1280 case(ReducedFaceElements):
1281 Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&_arg,&integrals[0]);
1282 break;
1283 case(Points):
1284 throw FinleyAdapterException("Error - Integral of data on points is not supported.");
1285 break;
1286 case(ContactElementsZero):
1287 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1288 break;
1289 case(ReducedContactElementsZero):
1290 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1291 break;
1292 case(ContactElementsOne):
1293 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1294 break;
1295 case(ReducedContactElementsOne):
1296 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1297 break;
1298 case(DegreesOfFreedom):
1299 temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
1300 _temp=temp.getDataC();
1301 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1302 break;
1303 case(ReducedDegreesOfFreedom):
1304 temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
1305 _temp=temp.getDataC();
1306 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1307 break;
1308 default:
1309 stringstream temp;
1310 temp << "Error - Integrals: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
1311 throw FinleyAdapterException(temp.str());
1312 break;
1313 }
1314 checkFinleyError();
1315 blocktimer_increment("integrate()", blocktimer_start);
1316 }
1317
1318 //
1319 // calculates the gradient of arg:
1320 //
1321 void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
1322 {
1323 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(*(arg.getFunctionSpace().getDomain()));
1324 if (argDomain!=*this)
1325 throw FinleyAdapterException("Error - Illegal domain of gradient argument");
1326 const MeshAdapter& gradDomain=dynamic_cast<const MeshAdapter&>(*(grad.getFunctionSpace().getDomain()));
1327 if (gradDomain!=*this)
1328 throw FinleyAdapterException("Error - Illegal domain of gradient");
1329
1330 Finley_Mesh* mesh=m_finleyMesh.get();
1331 escriptDataC _grad=grad.getDataC();
1332 escriptDataC nodeDataC;
1333 escript::Data temp;
1334 if (getMPISize()>1) {
1335 if( arg.getFunctionSpace().getTypeCode() == DegreesOfFreedom ) {
1336 temp=escript::Data( arg, continuousFunction(asAbstractContinuousDomain()) );
1337 nodeDataC = temp.getDataC();
1338 } else if( arg.getFunctionSpace().getTypeCode() == ReducedDegreesOfFreedom ) {
1339 temp=escript::Data( arg, reducedContinuousFunction(asAbstractContinuousDomain()) );
1340 nodeDataC = temp.getDataC();
1341 } else {
1342 nodeDataC = arg.getDataC();
1343 }
1344 } else {
1345 nodeDataC = arg.getDataC();
1346 }
1347 switch(grad.getFunctionSpace().getTypeCode()) {
1348 case(Nodes):
1349 throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
1350 break;
1351 case(ReducedNodes):
1352 throw FinleyAdapterException("Error - Gradient at reduced nodes is not supported.");
1353 break;
1354 case(Elements):
1355 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&_grad,&nodeDataC);
1356 break;
1357 case(ReducedElements):
1358 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&_grad,&nodeDataC);
1359 break;
1360 case(FaceElements):
1361 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&_grad,&nodeDataC);
1362 break;
1363 case(ReducedFaceElements):
1364 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&_grad,&nodeDataC);
1365 break;
1366 case(Points):
1367 throw FinleyAdapterException("Error - Gradient at points is not supported.");
1368 break;
1369 case(ContactElementsZero):
1370 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1371 break;
1372 case(ReducedContactElementsZero):
1373 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1374 break;
1375 case(ContactElementsOne):
1376 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1377 break;
1378 case(ReducedContactElementsOne):
1379 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&_grad,&nodeDataC);
1380 break;
1381 case(DegreesOfFreedom):
1382 throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
1383 break;
1384 case(ReducedDegreesOfFreedom):
1385 throw FinleyAdapterException("Error - Gradient at reduced degrees of freedom is not supported.");
1386 break;
1387 default:
1388 stringstream temp;
1389 temp << "Error - Gradient: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
1390 throw FinleyAdapterException(temp.str());
1391 break;
1392 }
1393 checkFinleyError();
1394 }
1395
1396 //
1397 // returns the size of elements:
1398 //
1399 void MeshAdapter::setToSize(escript::Data& size) const
1400 {
1401 Finley_Mesh* mesh=m_finleyMesh.get();
1402 escriptDataC tmp=size.getDataC();
1403 switch(size.getFunctionSpace().getTypeCode()) {
1404 case(Nodes):
1405 throw FinleyAdapterException("Error - Size of nodes is not supported.");
1406 break;
1407 case(ReducedNodes):
1408 throw FinleyAdapterException("Error - Size of reduced nodes is not supported.");
1409 break;
1410 case(Elements):
1411 Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
1412 break;
1413 case(ReducedElements):
1414 Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
1415 break;
1416 case(FaceElements):
1417 Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
1418 break;
1419 case(ReducedFaceElements):
1420 Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
1421 break;
1422 case(Points):
1423 throw FinleyAdapterException("Error - Size of point elements is not supported.");
1424 break;
1425 case(ContactElementsZero):
1426 case(ContactElementsOne):
1427 Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
1428 break;
1429 case(ReducedContactElementsZero):
1430 case(ReducedContactElementsOne):
1431 Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
1432 break;
1433 case(DegreesOfFreedom):
1434 throw FinleyAdapterException("Error - Size of degrees of freedom is not supported.");
1435 break;
1436 case(ReducedDegreesOfFreedom):
1437 throw FinleyAdapterException("Error - Size of reduced degrees of freedom is not supported.");
1438 break;
1439 default:
1440 stringstream temp;
1441 temp << "Error - Element size: Finley does not know anything about function space type " << size.getFunctionSpace().getTypeCode();
1442 throw FinleyAdapterException(temp.str());
1443 break;
1444 }
1445 checkFinleyError();
1446 }
1447
1448 // sets the location of nodes:
1449 void MeshAdapter::setNewX(const escript::Data& new_x)
1450 {
1451 Finley_Mesh* mesh=m_finleyMesh.get();
1452 escriptDataC tmp;
1453 const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(*(new_x.getFunctionSpace().getDomain()));
1454 if (newDomain!=*this)
1455 throw FinleyAdapterException("Error - Illegal domain of new point locations");
1456 tmp = new_x.getDataC();
1457 Finley_Mesh_setCoordinates(mesh,&tmp);
1458 checkFinleyError();
1459 }
1460
1461 // saves a data array in openDX format:
1462 void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const
1463 {
1464 unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1465 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1466 /* win32 refactor */
1467 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
1468 for(int i=0;i<num_data;i++)
1469 {
1470 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
1471 }
1472
1473 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
1474 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
1475
1476 boost::python::list keys=arg.keys();
1477 for (int i=0;i<num_data;++i) {
1478 std::string n=boost::python::extract<std::string>(keys[i]);
1479 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
1480 if (dynamic_cast<const MeshAdapter&>(*(d.getFunctionSpace().getDomain())) !=*this)
1481 throw FinleyAdapterException("Error in saveVTK: Data must be defined on same Domain");
1482 data[i]=d.getDataC();
1483 ptr_data[i]= &(data[i]);
1484 if (n.length()>MAX_namelength-1) {
1485 strncpy(names[i],n.c_str(),MAX_namelength-1);
1486 names[i][MAX_namelength-1]='\0';
1487 } else {
1488 strcpy(names[i],n.c_str());
1489 }
1490 }
1491 Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,names,ptr_data);
1492 checkFinleyError();
1493
1494 /* win32 refactor */
1495 TMPMEMFREE(data);
1496 TMPMEMFREE(ptr_data);
1497 for(int i=0;i<num_data;i++)
1498 {
1499 TMPMEMFREE(names[i]);
1500 }
1501 TMPMEMFREE(names);
1502
1503 return;
1504 }
1505
1506 // saves a data array in openVTK format:
1507 void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const
1508 {
1509 unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1510 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1511 /* win32 refactor */
1512 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
1513 for(int i=0;i<num_data;i++)
1514 {
1515 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
1516 }
1517
1518 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
1519 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
1520
1521 boost::python::list keys=arg.keys();
1522 for (int i=0;i<num_data;++i) {
1523 std::string n=boost::python::extract<std::string>(keys[i]);
1524 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
1525 if (dynamic_cast<const MeshAdapter&>(*(d.getFunctionSpace().getDomain())) !=*this)
1526 throw FinleyAdapterException("Error in saveVTK: Data must be defined on same Domain");
1527 data[i]=d.getDataC();
1528 ptr_data[i]=&(data[i]);
1529 if (n.length()>MAX_namelength-1) {
1530 strncpy(names[i],n.c_str(),MAX_namelength-1);
1531 names[i][MAX_namelength-1]='\0';
1532 } else {
1533 strcpy(names[i],n.c_str());
1534 }
1535 }
1536 Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,names,ptr_data);
1537
1538 checkFinleyError();
1539 /* win32 refactor */
1540 TMPMEMFREE(data);
1541 TMPMEMFREE(ptr_data);
1542 for(int i=0;i<num_data;i++)
1543 {
1544 TMPMEMFREE(names[i]);
1545 }
1546 TMPMEMFREE(names);
1547
1548 return;
1549 }
1550
1551
1552 // creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
1553 SystemMatrixAdapter MeshAdapter::newSystemMatrix(
1554 const int row_blocksize,
1555 const escript::FunctionSpace& row_functionspace,
1556 const int column_blocksize,
1557 const escript::FunctionSpace& column_functionspace,
1558 const int type) const
1559 {
1560 int reduceRowOrder=0;
1561 int reduceColOrder=0;
1562 // is the domain right?
1563 const MeshAdapter& row_domain=dynamic_cast<const MeshAdapter&>(*(row_functionspace.getDomain()));
1564 if (row_domain!=*this)
1565 throw FinleyAdapterException("Error - domain of row function space does not match the domain of matrix generator.");
1566 const MeshAdapter& col_domain=dynamic_cast<const MeshAdapter&>(*(column_functionspace.getDomain()));
1567 if (col_domain!=*this)
1568 throw FinleyAdapterException("Error - domain of columnn function space does not match the domain of matrix generator.");
1569 // is the function space type right
1570 if (row_functionspace.getTypeCode()==DegreesOfFreedom) {
1571 reduceRowOrder=0;
1572 } else if (row_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1573 reduceRowOrder=1;
1574 } else {
1575 throw FinleyAdapterException("Error - illegal function space type for system matrix rows.");
1576 }
1577 if (column_functionspace.getTypeCode()==DegreesOfFreedom) {
1578 reduceColOrder=0;
1579 } else if (column_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1580 reduceColOrder=1;
1581 } else {
1582 throw FinleyAdapterException("Error - illegal function space type for system matrix columns.");
1583 }
1584 // generate matrix:
1585
1586 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);
1587 checkFinleyError();
1588 Paso_SystemMatrix* fsystemMatrix;
1589 int trilinos = 0;
1590 if (trilinos) {
1591 #ifdef TRILINOS
1592 /* Allocation Epetra_VrbMatrix here */
1593 #endif
1594 }
1595 else {
1596 fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);
1597 }
1598 checkPasoError();
1599 Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1600 return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
1601 }
1602 // creates a TransportProblemAdapter
1603 TransportProblemAdapter MeshAdapter::newTransportProblem(
1604 const double theta,
1605 const int blocksize,
1606 const escript::FunctionSpace& functionspace,
1607 const int type) const
1608 {
1609 int reduceOrder=0;
1610 // is the domain right?
1611 const MeshAdapter& domain=dynamic_cast<const MeshAdapter&>(*(functionspace.getDomain()));
1612 if (domain!=*this)
1613 throw FinleyAdapterException("Error - domain of function space does not match the domain of transport problem generator.");
1614 // is the function space type right
1615 if (functionspace.getTypeCode()==DegreesOfFreedom) {
1616 reduceOrder=0;
1617 } else if (functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
1618 reduceOrder=1;
1619 } else {
1620 throw FinleyAdapterException("Error - illegal function space type for system matrix rows.");
1621 }
1622 // generate matrix:
1623
1624 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);
1625 checkFinleyError();
1626 Paso_FCTransportProblem* transportProblem;
1627 transportProblem=Paso_FCTransportProblem_alloc(theta,fsystemMatrixPattern,blocksize);
1628 checkPasoError();
1629 Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1630 return TransportProblemAdapter(transportProblem,theta,blocksize,functionspace);
1631 }
1632
1633 //
1634 // vtkObject MeshAdapter::createVtkObject() const
1635 // TODO:
1636 //
1637
1638 //
1639 // returns true if data at the atom_type is considered as being cell centered:
1640 bool MeshAdapter::isCellOriented(int functionSpaceCode) const
1641 {
1642 switch(functionSpaceCode) {
1643 case(Nodes):
1644 case(DegreesOfFreedom):
1645 case(ReducedDegreesOfFreedom):
1646 return false;
1647 break;
1648 case(Elements):
1649 case(FaceElements):
1650 case(Points):
1651 case(ContactElementsZero):
1652 case(ContactElementsOne):
1653 case(ReducedElements):
1654 case(ReducedFaceElements):
1655 case(ReducedContactElementsZero):
1656 case(ReducedContactElementsOne):
1657 return true;
1658 break;
1659 default:
1660 stringstream temp;
1661 temp << "Error - Cell: Finley does not know anything about function space type " << functionSpaceCode;
1662 throw FinleyAdapterException(temp.str());
1663 break;
1664 }
1665 checkFinleyError();
1666 return false;
1667 }
1668
1669 bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const
1670 {
1671 switch(functionSpaceType_source) {
1672 case(Nodes):
1673 switch(functionSpaceType_target) {
1674 case(Nodes):
1675 case(ReducedNodes):
1676 case(ReducedDegreesOfFreedom):
1677 case(DegreesOfFreedom):
1678 case(Elements):
1679 case(ReducedElements):
1680 case(FaceElements):
1681 case(ReducedFaceElements):
1682 case(Points):
1683 case(ContactElementsZero):
1684 case(ReducedContactElementsZero):
1685 case(ContactElementsOne):
1686 case(ReducedContactElementsOne):
1687 return true;
1688 default:
1689 stringstream temp;
1690 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1691 throw FinleyAdapterException(temp.str());
1692 }
1693 break;
1694 case(ReducedNodes):
1695 switch(functionSpaceType_target) {
1696 case(ReducedNodes):
1697 case(ReducedDegreesOfFreedom):
1698 case(Elements):
1699 case(ReducedElements):
1700 case(FaceElements):
1701 case(ReducedFaceElements):
1702 case(Points):
1703 case(ContactElementsZero):
1704 case(ReducedContactElementsZero):
1705 case(ContactElementsOne):
1706 case(ReducedContactElementsOne):
1707 return true;
1708 case(Nodes):
1709 case(DegreesOfFreedom):
1710 return false;
1711 default:
1712 stringstream temp;
1713 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1714 throw FinleyAdapterException(temp.str());
1715 }
1716 break;
1717 case(Elements):
1718 if (functionSpaceType_target==Elements) {
1719 return true;
1720 } else if (functionSpaceType_target==ReducedElements) {
1721 return true;
1722 } else {
1723 return false;
1724 }
1725 case(ReducedElements):
1726 if (functionSpaceType_target==ReducedElements) {
1727 return true;
1728 } else {
1729 return false;
1730 }
1731 case(FaceElements):
1732 if (functionSpaceType_target==FaceElements) {
1733 return true;
1734 } else if (functionSpaceType_target==ReducedFaceElements) {
1735 return true;
1736 } else {
1737 return false;
1738 }
1739 case(ReducedFaceElements):
1740 if (functionSpaceType_target==ReducedFaceElements) {
1741 return true;
1742 } else {
1743 return false;
1744 }
1745 case(Points):
1746 if (functionSpaceType_target==Points) {
1747 return true;
1748 } else {
1749 return false;
1750 }
1751 case(ContactElementsZero):
1752 case(ContactElementsOne):
1753 if (functionSpaceType_target==ContactElementsZero || functionSpaceType_target==ContactElementsOne) {
1754 return true;
1755 } else if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1756 return true;
1757 } else {
1758 return false;
1759 }
1760 case(ReducedContactElementsZero):
1761 case(ReducedContactElementsOne):
1762 if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1763 return true;
1764 } else {
1765 return false;
1766 }
1767 case(DegreesOfFreedom):
1768 switch(functionSpaceType_target) {
1769 case(ReducedDegreesOfFreedom):
1770 case(DegreesOfFreedom):
1771 case(Nodes):
1772 case(ReducedNodes):
1773 case(Elements):
1774 case(ReducedElements):
1775 case(Points):
1776 case(FaceElements):
1777 case(ReducedFaceElements):
1778 case(ContactElementsZero):
1779 case(ReducedContactElementsZero):
1780 case(ContactElementsOne):
1781 case(ReducedContactElementsOne):
1782 return true;
1783 default:
1784 stringstream temp;
1785 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1786 throw FinleyAdapterException(temp.str());
1787 }
1788 break;
1789 case(ReducedDegreesOfFreedom):
1790 switch(functionSpaceType_target) {
1791 case(ReducedDegreesOfFreedom):
1792 case(ReducedNodes):
1793 case(Elements):
1794 case(ReducedElements):
1795 case(FaceElements):
1796 case(ReducedFaceElements):
1797 case(Points):
1798 case(ContactElementsZero):
1799 case(ReducedContactElementsZero):
1800 case(ContactElementsOne):
1801 case(ReducedContactElementsOne):
1802 return true;
1803 case(Nodes):
1804 case(DegreesOfFreedom):
1805 return false;
1806 default:
1807 stringstream temp;
1808 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1809 throw FinleyAdapterException(temp.str());
1810 }
1811 break;
1812 default:
1813 stringstream temp;
1814 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;
1815 throw FinleyAdapterException(temp.str());
1816 break;
1817 }
1818 checkFinleyError();
1819 return false;
1820 }
1821
1822 bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const
1823 {
1824 return false;
1825 }
1826
1827 bool MeshAdapter::operator==(const AbstractDomain& other) const
1828 {
1829 const MeshAdapter* temp=dynamic_cast<const MeshAdapter*>(&other);
1830 if (temp!=0) {
1831 return (m_finleyMesh==temp->m_finleyMesh);
1832 } else {
1833 return false;
1834 }
1835 }
1836
1837 bool MeshAdapter::operator!=(const AbstractDomain& other) const
1838 {
1839 return !(operator==(other));
1840 }
1841
1842 int MeshAdapter::getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const
1843 {
1844 int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(preconditioner), SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
1845 checkPasoError();
1846 return out;
1847 }
1848 int MeshAdapter::getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const
1849 {
1850 int out=Paso_FCTransportProblem_getTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(preconditioner), SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
1851 checkPasoError();
1852 return out;
1853 }
1854
1855 escript::Data MeshAdapter::getX() const
1856 {
1857 return continuousFunction(asAbstractContinuousDomain()).getX();
1858 }
1859
1860 escript::Data MeshAdapter::getNormal() const
1861 {
1862 return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
1863 }
1864
1865 escript::Data MeshAdapter::getSize() const
1866 {
1867 return escript::function(asAbstractContinuousDomain()).getSize();
1868 }
1869
1870 int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1871 {
1872 int *out = NULL;
1873 Finley_Mesh* mesh=m_finleyMesh.get();
1874 switch (functionSpaceType) {
1875 case(Nodes):
1876 out=mesh->Nodes->Id;
1877 break;
1878 case(ReducedNodes):
1879 out=mesh->Nodes->reducedNodesId;
1880 break;
1881 case(Elements):
1882 out=mesh->Elements->Id;
1883 break;
1884 case(ReducedElements):
1885 out=mesh->Elements->Id;
1886 break;
1887 case(FaceElements):
1888 out=mesh->FaceElements->Id;
1889 break;
1890 case(ReducedFaceElements):
1891 out=mesh->FaceElements->Id;
1892 break;
1893 case(Points):
1894 out=mesh->Points->Id;
1895 break;
1896 case(ContactElementsZero):
1897 out=mesh->ContactElements->Id;
1898 break;
1899 case(ReducedContactElementsZero):
1900 out=mesh->ContactElements->Id;
1901 break;
1902 case(ContactElementsOne):
1903 out=mesh->ContactElements->Id;
1904 break;
1905 case(ReducedContactElementsOne):
1906 out=mesh->ContactElements->Id;
1907 break;
1908 case(DegreesOfFreedom):
1909 out=mesh->Nodes->degreesOfFreedomId;
1910 break;
1911 case(ReducedDegreesOfFreedom):
1912 out=mesh->Nodes->reducedDegreesOfFreedomId;
1913 break;
1914 default:
1915 stringstream temp;
1916 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1917 throw FinleyAdapterException(temp.str());
1918 break;
1919 }
1920 return out;
1921 }
1922 int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
1923 {
1924 int out=0;
1925 Finley_Mesh* mesh=m_finleyMesh.get();
1926 switch (functionSpaceType) {
1927 case(Nodes):
1928 out=mesh->Nodes->Tag[sampleNo];
1929 break;
1930 case(ReducedNodes):
1931 throw FinleyAdapterException(" Error - ReducedNodes does not support tags.");
1932 break;
1933 case(Elements):
1934 out=mesh->Elements->Tag[sampleNo];
1935 break;
1936 case(ReducedElements):
1937 out=mesh->Elements->Tag[sampleNo];
1938 break;
1939 case(FaceElements):
1940 out=mesh->FaceElements->Tag[sampleNo];
1941 break;
1942 case(ReducedFaceElements):
1943 out=mesh->FaceElements->Tag[sampleNo];
1944 break;
1945 case(Points):
1946 out=mesh->Points->Tag[sampleNo];
1947 break;
1948 case(ContactElementsZero):
1949 out=mesh->ContactElements->Tag[sampleNo];
1950 break;
1951 case(ReducedContactElementsZero):
1952 out=mesh->ContactElements->Tag[sampleNo];
1953 break;
1954 case(ContactElementsOne):
1955 out=mesh->ContactElements->Tag[sampleNo];
1956 break;
1957 case(ReducedContactElementsOne):
1958 out=mesh->ContactElements->Tag[sampleNo];
1959 break;
1960 case(DegreesOfFreedom):
1961 throw FinleyAdapterException(" Error - DegreesOfFreedom does not support tags.");
1962 break;
1963 case(ReducedDegreesOfFreedom):
1964 throw FinleyAdapterException(" Error - ReducedDegreesOfFreedom does not support tags.");
1965 break;
1966 default:
1967 stringstream temp;
1968 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1969 throw FinleyAdapterException(temp.str());
1970 break;
1971 }
1972 return out;
1973 }
1974
1975
1976 void MeshAdapter::setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const
1977 {
1978 Finley_Mesh* mesh=m_finleyMesh.get();
1979 escriptDataC tmp=mask.getDataC();
1980 switch(functionSpaceType) {
1981 case(Nodes):
1982 Finley_NodeFile_setTags(mesh->Nodes,newTag,&tmp);
1983 break;
1984 case(ReducedNodes):
1985 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
1986 break;
1987 case(DegreesOfFreedom):
1988 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
1989 break;
1990 case(ReducedDegreesOfFreedom):
1991 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
1992 break;
1993 case(Elements):
1994 Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1995 break;
1996 case(ReducedElements):
1997 Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1998 break;
1999 case(FaceElements):
2000 Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
2001 break;
2002 case(ReducedFaceElements):
2003 Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
2004 break;
2005 case(Points):
2006 Finley_ElementFile_setTags(mesh->Points,newTag,&tmp);
2007 break;
2008 case(ContactElementsZero):
2009 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
2010 break;
2011 case(ReducedContactElementsZero):
2012 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
2013 break;
2014 case(ContactElementsOne):
2015 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
2016 break;
2017 case(ReducedContactElementsOne):
2018 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
2019 break;
2020 default:
2021 stringstream temp;
2022 temp << "Error - Finley does not know anything about function space type " << functionSpaceType;
2023 throw FinleyAdapterException(temp.str());
2024 }
2025 checkFinleyError();
2026 return;
2027 }
2028
2029 void MeshAdapter::setTagMap(const std::string& name, int tag)
2030 {
2031 Finley_Mesh* mesh=m_finleyMesh.get();
2032 Finley_Mesh_addTagMap(mesh, name.c_str(),tag);
2033 checkPasoError();
2034 // throwStandardException("MeshAdapter::set TagMap is not implemented.");
2035 }
2036
2037 int MeshAdapter::getTag(const std::string& name) const
2038 {
2039 Finley_Mesh* mesh=m_finleyMesh.get();
2040 int tag=0;
2041 tag=Finley_Mesh_getTag(mesh, name.c_str());
2042 checkPasoError();
2043 // throwStandardException("MeshAdapter::getTag is not implemented.");
2044 return tag;
2045 }
2046
2047 bool MeshAdapter::isValidTagName(const std::string& name) const
2048 {
2049 Finley_Mesh* mesh=m_finleyMesh.get();
2050 return Finley_Mesh_isValidTagName(mesh,name.c_str());
2051 }
2052
2053 std::string MeshAdapter::showTagNames() const
2054 {
2055 stringstream temp;
2056 Finley_Mesh* mesh=m_finleyMesh.get();
2057 Finley_TagMap* tag_map=mesh->TagMap;
2058 while (tag_map) {
2059 temp << tag_map->name;
2060 tag_map=tag_map->next;
2061 if (tag_map) temp << ", ";
2062 }
2063 return temp.str();
2064 }
2065
2066 int MeshAdapter::getNumberOfTagsInUse(int functionSpaceCode) const
2067 {
2068 Finley_Mesh* mesh=m_finleyMesh.get();
2069 dim_t numTags=0;
2070 switch(functionSpaceCode) {
2071 case(Nodes):
2072 numTags=mesh->Nodes->numTagsInUse;
2073 break;
2074 case(ReducedNodes):
2075 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
2076 break;
2077 case(DegreesOfFreedom):
2078 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
2079 break;
2080 case(ReducedDegreesOfFreedom):
2081 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
2082 break;
2083 case(Elements):
2084 case(ReducedElements):
2085 numTags=mesh->Elements->numTagsInUse;
2086 break;
2087 case(FaceElements):
2088 case(ReducedFaceElements):
2089 numTags=mesh->FaceElements->numTagsInUse;
2090 break;
2091 case(Points):
2092 numTags=mesh->Points->numTagsInUse;
2093 break;
2094 case(ContactElementsZero):
2095 case(ReducedContactElementsZero):
2096 case(ContactElementsOne):
2097 case(ReducedContactElementsOne):
2098 numTags=mesh->ContactElements->numTagsInUse;
2099 break;
2100 default:
2101 stringstream temp;
2102 temp << "Error - Finley does not know anything about function space type " << functionSpaceCode;
2103 throw FinleyAdapterException(temp.str());
2104 }
2105 return numTags;
2106 }
2107 int* MeshAdapter::borrowListOfTagsInUse(int functionSpaceCode) const
2108 {
2109 Finley_Mesh* mesh=m_finleyMesh.get();
2110 index_t* tags=NULL;
2111 switch(functionSpaceCode) {
2112 case(Nodes):
2113 tags=mesh->Nodes->tagsInUse;
2114 break;
2115 case(ReducedNodes):
2116 throw FinleyAdapterException("Error - ReducedNodes does not support tags");
2117 break;
2118 case(DegreesOfFreedom):
2119 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
2120 break;
2121 case(ReducedDegreesOfFreedom):
2122 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
2123 break;
2124 case(Elements):
2125 case(ReducedElements):
2126 tags=mesh->Elements->tagsInUse;
2127 break;
2128 case(FaceElements):
2129 case(ReducedFaceElements):
2130 tags=mesh->FaceElements->tagsInUse;
2131 break;
2132 case(Points):
2133 tags=mesh->Points->tagsInUse;
2134 break;
2135 case(ContactElementsZero):
2136 case(ReducedContactElementsZero):
2137 case(ContactElementsOne):
2138 case(ReducedContactElementsOne):
2139 tags=mesh->ContactElements->tagsInUse;
2140 break;
2141 default:
2142 stringstream temp;
2143 temp << "Error - Finley does not know anything about function space type " << functionSpaceCode;
2144 throw FinleyAdapterException(temp.str());
2145 }
2146 return tags;
2147 }
2148
2149
2150 bool MeshAdapter::canTag(int functionSpaceCode) const
2151 {
2152 switch(functionSpaceCode) {
2153 case(Nodes):
2154 case(Elements):
2155 case(ReducedElements):
2156 case(FaceElements):
2157 case(ReducedFaceElements):
2158 case(Points):
2159 case(ContactElementsZero):
2160 case(ReducedContactElementsZero):
2161 case(ContactElementsOne):
2162 case(ReducedContactElementsOne):
2163 return true;
2164 case(ReducedNodes):
2165 case(DegreesOfFreedom):
2166 case(ReducedDegreesOfFreedom):
2167 return false;
2168 default:
2169 return false;
2170 }
2171 }
2172
2173
2174 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26