/[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 4687 - (show annotations)
Wed Feb 19 00:03:29 2014 UTC (5 years, 2 months ago) by jfenwick
Original Path: trunk/finley/src/CPPAdapter/MeshAdapter.cpp
File size: 91259 byte(s)
Remove randomFill python method from ripley domains.
All random data objects (for all domain types) should be generated 
using esys.escript.RandomData()

The only filtered random we have is gaussian on ripley but
it is triggered by passing the tuple as the last arg of RandomData().

While the interface is a bit more complicated (in that you always need
 to pass in shape and functionspace) it does mean we have a 
common interface for all domains. 

Removed randomFill from DataExpanded.
The reasoning behind this is to force domains to call the util function
themselves and enforce whatever consistancy requirements they have.

Added version of blocktools to deal with 2D case in Ripley.
Use blocktools for the 2D transfers [This was cleaner than modifying the
previous implementation to deal with variable shaped points].

Note that under MPI, ripley can not generate random data (even unfiltered)
if any of its per rank dimensions is <4 elements on any side.

Unit tests for these calls are in but some extra checks still needed.



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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26