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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1059 - (show annotations)
Fri Mar 23 11:01:53 2007 UTC (12 years, 3 months ago) by gross
File size: 44072 byte(s)
first steps toward reduced element integration order. The escript bit is done but the finley part still needs work.
1 // $Id$
2 /*
3 ******************************************************************************
4 * *
5 * COPYRIGHT ACcESS 2004 - All Rights Reserved *
6 * *
7 * This software is the property of ACcESS. No part of this code *
8 * may be copied in any form or by any means without the expressed written *
9 * consent of ACcESS. Copying, use or modification of this software *
10 * by any unauthorised person is illegal unless that person has a software *
11 * license agreement with ACcESS. *
12 * *
13 ******************************************************************************
14 */
15
16 #ifdef PASO_MPI
17 #include <mpi.h>
18 #endif
19 #include "MeshAdapter.h"
20
21 #include "escript/Data.h"
22 #include "escript/DataFactory.h"
23
24 using namespace std;
25 using namespace escript;
26
27 namespace finley {
28
29 //
30 // define the static constants
31 MeshAdapter::FunctionSpaceNamesMapType MeshAdapter::m_functionSpaceTypeNames;
32 const int MeshAdapter::DegreesOfFreedom=FINLEY_DEGREES_OF_FREEDOM;
33 const int MeshAdapter::ReducedDegreesOfFreedom=FINLEY_REDUCED_DEGREES_OF_FREEDOM;
34 const int MeshAdapter::Nodes=FINLEY_NODES;
35 const int MeshAdapter::Elements=FINLEY_ELEMENTS;
36 const int MeshAdapter::ReducedElements=FINLEY_REDUCED_ELEMENTS;
37 const int MeshAdapter::FaceElements=FINLEY_FACE_ELEMENTS;
38 const int MeshAdapter::ReducedFaceElements=FINLEY_REDUCED_FACE_ELEMENTS;
39 const int MeshAdapter::Points=FINLEY_POINTS;
40 const int MeshAdapter::ContactElementsZero=FINLEY_CONTACT_ELEMENTS_1;
41 const int MeshAdapter::ReducedContactElementsZero=FINLEY_REDUCED_CONTACT_ELEMENTS_1;
42 const int MeshAdapter::ContactElementsOne=FINLEY_CONTACT_ELEMENTS_2;
43 const int MeshAdapter::ReducedContactElementsOne=FINLEY_REDUCED_CONTACT_ELEMENTS_2;
44
45 MeshAdapter::MeshAdapter(Finley_Mesh* finleyMesh)
46 {
47 setFunctionSpaceTypeNames();
48 //
49 // need to use a null_deleter as Finley_Mesh_dealloc deletes the pointer
50 // for us.
51 m_finleyMesh.reset(finleyMesh,null_deleter());
52 }
53
54 //
55 // The copy constructor should just increment the use count
56 MeshAdapter::MeshAdapter(const MeshAdapter& in):
57 m_finleyMesh(in.m_finleyMesh)
58 {
59 setFunctionSpaceTypeNames();
60 }
61
62 MeshAdapter::~MeshAdapter()
63 {
64 //
65 // I hope the case for the pointer being zero has been taken care of.
66 // cout << "In MeshAdapter destructor." << endl;
67 if (m_finleyMesh.unique()) {
68 // cout << "Calling dealloc." << endl;
69 Finley_Mesh_dealloc(m_finleyMesh.get());
70 // cout << "Finished dealloc." << endl;
71 }
72 }
73
74 Finley_Mesh* MeshAdapter::getFinley_Mesh() const {
75 return m_finleyMesh.get();
76 }
77
78 void MeshAdapter::write(const std::string& fileName) const
79 {
80 char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
81 strcpy(fName,fileName.c_str());
82 Finley_Mesh_write(m_finleyMesh.get(),fName);
83 checkFinleyError();
84 TMPMEMFREE(fName);
85 }
86
87 string MeshAdapter::getDescription() const
88 {
89 return "FinleyMesh";
90 }
91
92 string MeshAdapter::functionSpaceTypeAsString(int functionSpaceType) const
93 {
94 FunctionSpaceNamesMapType::iterator loc;
95 loc=m_functionSpaceTypeNames.find(functionSpaceType);
96 if (loc==m_functionSpaceTypeNames.end()) {
97 return "Invalid function space type code.";
98 } else {
99 return loc->second;
100 }
101 }
102
103 bool MeshAdapter::isValidFunctionSpaceType(int functionSpaceType) const
104 {
105 FunctionSpaceNamesMapType::iterator loc;
106 loc=m_functionSpaceTypeNames.find(functionSpaceType);
107 return (loc!=m_functionSpaceTypeNames.end());
108 }
109
110 void MeshAdapter::setFunctionSpaceTypeNames()
111 {
112 m_functionSpaceTypeNames.insert
113 (FunctionSpaceNamesMapType::value_type(DegreesOfFreedom,"Finley_DegreesOfFreedom"));
114 m_functionSpaceTypeNames.insert
115 (FunctionSpaceNamesMapType::value_type(ReducedDegreesOfFreedom,"Finley_ReducedDegreesOfFreedom"));
116 m_functionSpaceTypeNames.insert
117 (FunctionSpaceNamesMapType::value_type(Nodes,"Finley_Nodes"));
118 m_functionSpaceTypeNames.insert
119 (FunctionSpaceNamesMapType::value_type(Elements,"Finley_Elements"));
120 m_functionSpaceTypeNames.insert
121 (FunctionSpaceNamesMapType::value_type(ReducedElements,"Finley_Reduced_Elements"));
122 m_functionSpaceTypeNames.insert
123 (FunctionSpaceNamesMapType::value_type(FaceElements,"Finley_Face_Elements"));
124 m_functionSpaceTypeNames.insert
125 (FunctionSpaceNamesMapType::value_type(ReducedFaceElements,"Finley_Reduced_Face_Elements"));
126 m_functionSpaceTypeNames.insert
127 (FunctionSpaceNamesMapType::value_type(Points,"Finley_Points"));
128 m_functionSpaceTypeNames.insert
129 (FunctionSpaceNamesMapType::value_type(ContactElementsZero,"Finley_Contact_Elements_0"));
130 m_functionSpaceTypeNames.insert
131 (FunctionSpaceNamesMapType::value_type(ReducedContactElementsZero,"Finley_Reduced_Contact_Elements_0"));
132 m_functionSpaceTypeNames.insert
133 (FunctionSpaceNamesMapType::value_type(ContactElementsOne,"Finley_Contact_Elements_1"));
134 m_functionSpaceTypeNames.insert
135 (FunctionSpaceNamesMapType::value_type(ReducedContactElementsOne,"Finley_Reduced_Contact_Elements_1"));
136 }
137
138 int MeshAdapter::getContinuousFunctionCode() const
139 {
140 return Nodes;
141 }
142
143 int MeshAdapter::getFunctionCode() const
144 {
145 return Elements;
146 }
147 int MeshAdapter::getReducedFunctionCode() const
148 {
149 return ReducedElements;
150 }
151
152 int MeshAdapter::getFunctionOnBoundaryCode() const
153 {
154 return FaceElements;
155 }
156 int MeshAdapter::getReducedFunctionOnBoundaryCode() const
157 {
158 return ReducedFaceElements;
159 }
160
161 int MeshAdapter::getFunctionOnContactZeroCode() const
162 {
163 return ContactElementsZero;
164 }
165 int MeshAdapter::getReducedFunctionOnContactZeroCode() const
166 {
167 return ReducedContactElementsZero;
168 }
169
170 int MeshAdapter::getFunctionOnContactOneCode() const
171 {
172 return ContactElementsOne;
173 }
174 int MeshAdapter::getReducedFunctionOnContactOneCode() const
175 {
176 return ReducedContactElementsOne;
177 }
178
179 int MeshAdapter::getSolutionCode() const
180 {
181 return DegreesOfFreedom;
182 }
183
184 int MeshAdapter::getReducedSolutionCode() const
185 {
186 return ReducedDegreesOfFreedom;
187 }
188
189 int MeshAdapter::getDiracDeltaFunctionCode() const
190 {
191 return Points;
192 }
193
194 //
195 // return the spatial dimension of the Mesh:
196 //
197 int MeshAdapter::getDim() const
198 {
199 int numDim=Finley_Mesh_getDim(m_finleyMesh.get());
200 checkFinleyError();
201 return numDim;
202 }
203
204 //
205 // return the number of data points per sample and the number of samples
206 // needed to represent data on a parts of the mesh.
207 //
208 pair<int,int> MeshAdapter::getDataShape(int functionSpaceCode) const
209 {
210 int numDataPointsPerSample=0;
211 int numSamples=0;
212 Finley_Mesh* mesh=m_finleyMesh.get();
213 switch (functionSpaceCode) {
214 case(Nodes):
215 numDataPointsPerSample=1;
216 if (mesh->Nodes!=NULL) numSamples=mesh->Nodes->numNodes;
217 break;
218 case(Elements):
219 if (mesh->Elements!=NULL) {
220 numSamples=mesh->Elements->numElements;
221 numDataPointsPerSample=mesh->Elements->ReferenceElement->numQuadNodes;
222 }
223 break;
224 case(FaceElements):
225 if (mesh->FaceElements!=NULL) {
226 numDataPointsPerSample=mesh->FaceElements->ReferenceElement->numQuadNodes;
227 numSamples=mesh->FaceElements->numElements;
228 }
229 break;
230 case(Points):
231 if (mesh->Points!=NULL) {
232 numDataPointsPerSample=1;
233 numSamples=mesh->Points->numElements;
234 }
235 break;
236 case(ContactElementsZero):
237 if (mesh->ContactElements!=NULL) {
238 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
239 numSamples=mesh->ContactElements->numElements;
240 }
241 break;
242 case(ContactElementsOne):
243 if (mesh->ContactElements!=NULL) {
244 numDataPointsPerSample=mesh->ContactElements->ReferenceElement->numQuadNodes;
245 numSamples=mesh->ContactElements->numElements;
246 }
247 break;
248 case(DegreesOfFreedom):
249 if (mesh->Nodes!=NULL) {
250 numDataPointsPerSample=1;
251 #ifndef PASO_MPI
252 numSamples=mesh->Nodes->numDegreesOfFreedom;
253 #else
254 numSamples=mesh->Nodes->degreeOfFreedomDistribution->numLocal;
255 #endif
256 }
257 break;
258 case(ReducedDegreesOfFreedom):
259 if (mesh->Nodes!=NULL) {
260 numDataPointsPerSample=1;
261 #ifndef PASO_MPI
262 numSamples=mesh->Nodes->reducedNumDegreesOfFreedom;
263 #else
264 numSamples=mesh->Nodes->reducedDegreeOfFreedomDistribution->numLocal;
265 #endif
266 }
267 break;
268 default:
269 stringstream temp;
270 temp << "Error - Invalid function space type: " << functionSpaceCode << " for domain: " << getDescription();
271 throw FinleyAdapterException(temp.str());
272 break;
273 }
274 return pair<int,int>(numDataPointsPerSample,numSamples);
275 }
276
277 //
278 // adds linear PDE of second order into a given stiffness matrix and right hand side:
279 //
280 void MeshAdapter::addPDEToSystem(
281 SystemMatrixAdapter& mat, escript::Data& rhs,
282 const escript::Data& A, const escript::Data& B, const escript::Data& C,const escript::Data& D,const escript::Data& X,const escript::Data& Y,
283 const escript::Data& d, const escript::Data& y,
284 const escript::Data& d_contact,const escript::Data& y_contact) const
285 {
286 escriptDataC _rhs=rhs.getDataC();
287 escriptDataC _A =A.getDataC();
288 escriptDataC _B=B.getDataC();
289 escriptDataC _C=C.getDataC();
290 escriptDataC _D=D.getDataC();
291 escriptDataC _X=X.getDataC();
292 escriptDataC _Y=Y.getDataC();
293 escriptDataC _d=d.getDataC();
294 escriptDataC _y=y.getDataC();
295 escriptDataC _d_contact=d_contact.getDataC();
296 escriptDataC _y_contact=y_contact.getDataC();
297
298 Finley_Mesh* mesh=m_finleyMesh.get();
299
300 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(), &_rhs, &_A, &_B, &_C, &_D, &_X, &_Y );
301 checkFinleyError();
302
303 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, mat.getPaso_SystemMatrix(), &_rhs, 0, 0, 0, &_d, 0, &_y );
304 checkFinleyError();
305
306 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, mat.getPaso_SystemMatrix(), &_rhs , 0, 0, 0, &_d_contact, 0, &_y_contact );
307 checkFinleyError();
308 }
309
310 //
311 // adds linear PDE of second order into the right hand side only
312 //
313 void MeshAdapter::addPDEToRHS( escript::Data& rhs, const escript::Data& X,const escript::Data& Y, const escript::Data& y, const escript::Data& y_contact) const
314 {
315 Finley_Mesh* mesh=m_finleyMesh.get();
316
317 escriptDataC _rhs=rhs.getDataC();
318 escriptDataC _X=X.getDataC();
319 escriptDataC _Y=Y.getDataC();
320 escriptDataC _y=y.getDataC();
321 escriptDataC _y_contact=y_contact.getDataC();
322
323 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements, 0, &_rhs, 0, 0, 0, 0, &_X, &_Y );
324 checkFinleyError();
325
326 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, 0, &_rhs, 0, 0, 0, 0, 0, &_y );
327 checkFinleyError();
328
329 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, 0, &_rhs , 0, 0, 0, 0, 0, &_y_contact );
330 checkFinleyError();
331 }
332
333 //
334 // interpolates data between different function spaces:
335 //
336 void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
337 {
338 const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());
339 const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
340 if (inDomain!=*this)
341 throw FinleyAdapterException("Error - Illegal domain of interpolant.");
342 if (targetDomain!=*this)
343 throw FinleyAdapterException("Error - Illegal domain of interpolation target.");
344
345 Finley_Mesh* mesh=m_finleyMesh.get();
346 switch(in.getFunctionSpace().getTypeCode()) {
347 case(Nodes):
348 switch(target.getFunctionSpace().getTypeCode()) {
349 case(Nodes):
350 case(ReducedDegreesOfFreedom):
351 case(DegreesOfFreedom):
352 Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));
353 break;
354 case(Elements):
355 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(in.getDataC()),&(target.getDataC()));
356 break;
357 case(FaceElements):
358 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&(in.getDataC()),&(target.getDataC()));
359 break;
360 case(Points):
361 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&(in.getDataC()),&(target.getDataC()));
362 break;
363 case(ContactElementsZero):
364 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
365 break;
366 case(ContactElementsOne):
367 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
368 break;
369 default:
370 stringstream temp;
371 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
372 throw FinleyAdapterException(temp.str());
373 break;
374 }
375 break;
376 case(Elements):
377 if (target.getFunctionSpace().getTypeCode()==Elements) {
378 Finley_Assemble_CopyElementData(mesh->Elements,&(target.getDataC()),&(in.getDataC()));
379 } else {
380 throw FinleyAdapterException("Error - No interpolation with data on elements possible.");
381 }
382 break;
383 case(FaceElements):
384 if (target.getFunctionSpace().getTypeCode()==FaceElements) {
385 Finley_Assemble_CopyElementData(mesh->FaceElements,&(target.getDataC()),&(in.getDataC()));
386 } else {
387 throw FinleyAdapterException("Error - No interpolation with data on face elements possible.");
388 break;
389 }
390 case(Points):
391 if (target.getFunctionSpace().getTypeCode()==Points) {
392 Finley_Assemble_CopyElementData(mesh->Points,&(target.getDataC()),&(in.getDataC()));
393 } else {
394 throw FinleyAdapterException("Error - No interpolation with data on points possible.");
395 break;
396 }
397 break;
398 case(ContactElementsZero):
399 case(ContactElementsOne):
400 if (target.getFunctionSpace().getTypeCode()==ContactElementsZero || target.getFunctionSpace().getTypeCode()==ContactElementsOne) {
401 Finley_Assemble_CopyElementData(mesh->ContactElements,&(target.getDataC()),&(in.getDataC()));
402 } else {
403 throw FinleyAdapterException("Error - No interpolation with data on contact elements possible.");
404 break;
405 }
406 break;
407 case(DegreesOfFreedom):
408 switch(target.getFunctionSpace().getTypeCode()) {
409 case(ReducedDegreesOfFreedom):
410 case(DegreesOfFreedom):
411 case(Nodes):
412 Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));
413 break;
414 #ifndef PASO_MPI
415 case(Elements):
416 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(in.getDataC()),&(target.getDataC()));
417 break;
418 case(FaceElements):
419 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&(in.getDataC()),&(target.getDataC()));
420 break;
421 case(Points):
422 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&(in.getDataC()),&(target.getDataC()));
423 break;
424 case(ContactElementsZero):
425 case(ContactElementsOne):
426 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
427 break;
428 #else
429 /* need to copy Degrees of freedom data to nodal data so that the external values are available */
430 case(Elements):
431 {
432 escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
433 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(nodeTemp.getDataC()),&(target.getDataC()));
434 break;
435 }
436 case(FaceElements):
437 {
438 escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
439 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&(nodeTemp.getDataC()),&(target.getDataC()));
440 break;
441 }
442 case(Points):
443 {
444 escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
445 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&(nodeTemp.getDataC()),&(target.getDataC()));
446 break;
447 }
448 case(ContactElementsZero):
449 case(ContactElementsOne):
450 {
451 escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
452 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(nodeTemp.getDataC()),&(target.getDataC()));
453 break;
454 }
455 #endif
456 default:
457 stringstream temp;
458 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
459 throw FinleyAdapterException(temp.str());
460 break;
461 }
462 break;
463 case(ReducedDegreesOfFreedom):
464 switch(target.getFunctionSpace().getTypeCode()) {
465 case(ReducedDegreesOfFreedom):
466 Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));
467 break;
468 case(Elements):
469 Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(in.getDataC()),&(target.getDataC()));
470 break;
471 case(FaceElements):
472 Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&(in.getDataC()),&(target.getDataC()));
473 break;
474 case(Points):
475 Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&(in.getDataC()),&(target.getDataC()));
476 break;
477 case(ContactElementsZero):
478 case(ContactElementsOne):
479 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
480 break;
481 case(Nodes):
482 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to mesh nodes.");
483 break;
484 case(DegreesOfFreedom):
485 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");
486 break;
487 default:
488 stringstream temp;
489 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
490 throw FinleyAdapterException(temp.str());
491 break;
492 }
493 break;
494 default:
495 stringstream temp;
496 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type %d" << in.getFunctionSpace().getTypeCode();
497 throw FinleyAdapterException(temp.str());
498 break;
499 }
500 checkFinleyError();
501 }
502
503 //
504 // copies the locations of sample points into x:
505 //
506 void MeshAdapter::setToX(escript::Data& arg) const
507 {
508 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
509 if (argDomain!=*this)
510 throw FinleyAdapterException("Error - Illegal domain of data point locations");
511 Finley_Mesh* mesh=m_finleyMesh.get();
512 // in case of values node coordinates we can do the job directly:
513 if (arg.getFunctionSpace().getTypeCode()==Nodes) {
514 Finley_Assemble_NodeCoordinates(mesh->Nodes,&(arg.getDataC()));
515 } else {
516 escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
517 Finley_Assemble_NodeCoordinates(mesh->Nodes,&(tmp_data.getDataC()));
518 // this is then interpolated onto arg:
519 interpolateOnDomain(arg,tmp_data);
520 }
521 checkFinleyError();
522 }
523
524 //
525 // return the normal vectors at the location of data points as a Data object:
526 //
527 void MeshAdapter::setToNormal(escript::Data& normal) const
528 {
529 const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());
530 if (normalDomain!=*this)
531 throw FinleyAdapterException("Error - Illegal domain of normal locations");
532 Finley_Mesh* mesh=m_finleyMesh.get();
533 switch(normal.getFunctionSpace().getTypeCode()) {
534 case(Nodes):
535 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for nodes");
536 break;
537 case(Elements):
538 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements");
539 break;
540 case (FaceElements):
541 Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&(normal.getDataC()));
542 break;
543 case(Points):
544 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for point elements");
545 break;
546 case (ContactElementsOne):
547 case (ContactElementsZero):
548 Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&(normal.getDataC()));
549 break;
550 case(DegreesOfFreedom):
551 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for degrees of freedom.");
552 break;
553 case(ReducedDegreesOfFreedom):
554 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced degrees of freedom.");
555 break;
556 default:
557 stringstream temp;
558 temp << "Error - Normal Vectors: Finley does not know anything about function space type " << normal.getFunctionSpace().getTypeCode();
559 throw FinleyAdapterException(temp.str());
560 break;
561 }
562 checkFinleyError();
563 }
564
565 //
566 // interpolates data to other domain:
567 //
568 void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
569 {
570 const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
571 if (targetDomain!=*this)
572 throw FinleyAdapterException("Error - Illegal domain of interpolation target");
573
574 throw FinleyAdapterException("Error - Finley does not allow interpolation across domains yet.");
575 }
576
577 //
578 // calculates the integral of a function defined of arg:
579 //
580 void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
581 {
582 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
583 if (argDomain!=*this)
584 throw FinleyAdapterException("Error - Illegal domain of integration kernel");
585
586 Finley_Mesh* mesh=m_finleyMesh.get();
587 switch(arg.getFunctionSpace().getTypeCode()) {
588 case(Nodes):
589 throw FinleyAdapterException("Error - Integral of data on nodes is not supported.");
590 break;
591 case(Elements):
592 Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&(arg.getDataC()),&integrals[0]);
593 break;
594 case(FaceElements):
595 Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&(arg.getDataC()),&integrals[0]);
596 break;
597 case(Points):
598 throw FinleyAdapterException("Error - Integral of data on points is not supported.");
599 break;
600 case(ContactElementsZero):
601 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);
602 break;
603 case(ContactElementsOne):
604 Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);
605 break;
606 case(DegreesOfFreedom):
607 throw FinleyAdapterException("Error - Integral of data on degrees of freedom is not supported.");
608 break;
609 case(ReducedDegreesOfFreedom):
610 throw FinleyAdapterException("Error - Integral of data on reduced degrees of freedom is not supported.");
611 break;
612 default:
613 stringstream temp;
614 temp << "Error - Integrals: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
615 throw FinleyAdapterException(temp.str());
616 break;
617 }
618 checkFinleyError();
619 }
620
621 //
622 // calculates the gradient of arg:
623 //
624 void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
625 {
626 const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
627 if (argDomain!=*this)
628 throw FinleyAdapterException("Error - Illegal domain of gradient argument");
629 const MeshAdapter& gradDomain=dynamic_cast<const MeshAdapter&>(grad.getFunctionSpace().getDomain());
630 if (gradDomain!=*this)
631 throw FinleyAdapterException("Error - Illegal domain of gradient");
632
633 Finley_Mesh* mesh=m_finleyMesh.get();
634 escriptDataC nodeDataC;
635 #ifdef PASO_MPI
636 escript::Data nodeTemp( arg, continuousFunction(asAbstractContinuousDomain()) );
637 if( arg.getFunctionSpace().getTypeCode() != Nodes )
638 {
639 Finley_Assemble_CopyNodalData(mesh->Nodes,&(nodeTemp.getDataC()),&(arg.getDataC()));
640 nodeDataC = nodeTemp.getDataC();
641 }
642 else
643 nodeDataC = arg.getDataC();
644 #else
645 nodeDataC = arg.getDataC();
646 #endif
647 switch(grad.getFunctionSpace().getTypeCode()) {
648 case(Nodes):
649 throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
650 break;
651 case(Elements):
652 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&nodeDataC);
653 break;
654 case(FaceElements):
655 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&nodeDataC);
656 break;
657 case(Points):
658 throw FinleyAdapterException("Error - Gradient at points is not supported.");
659 break;
660 case(ContactElementsZero):
661 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
662 break;
663 case(ContactElementsOne):
664 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
665 break;
666 case(DegreesOfFreedom):
667 throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
668 break;
669 case(ReducedDegreesOfFreedom):
670 throw FinleyAdapterException("Error - Gradient at reduced degrees of freedom is not supported.");
671 break;
672 default:
673 stringstream temp;
674 temp << "Error - Gradient: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
675 throw FinleyAdapterException(temp.str());
676 break;
677 }
678 checkFinleyError();
679 }
680
681 //
682 // returns the size of elements:
683 //
684 void MeshAdapter::setToSize(escript::Data& size) const
685 {
686 Finley_Mesh* mesh=m_finleyMesh.get();
687 escriptDataC tmp=size.getDataC();
688 switch(size.getFunctionSpace().getTypeCode()) {
689 case(Nodes):
690 throw FinleyAdapterException("Error - Size of nodes is not supported.");
691 break;
692 case(Elements):
693 Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
694 break;
695 case(FaceElements):
696 Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
697 break;
698 case(Points):
699 throw FinleyAdapterException("Error - Size of point elements is not supported.");
700 break;
701 case(ContactElementsZero):
702 case(ContactElementsOne):
703 Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
704 break;
705 case(DegreesOfFreedom):
706 throw FinleyAdapterException("Error - Size of degrees of freedom is not supported.");
707 break;
708 case(ReducedDegreesOfFreedom):
709 throw FinleyAdapterException("Error - Size of reduced degrees of freedom is not supported.");
710 break;
711 default:
712 stringstream temp;
713 temp << "Error - Element size: Finley does not know anything about function space type " << size.getFunctionSpace().getTypeCode();
714 throw FinleyAdapterException(temp.str());
715 break;
716 }
717 checkFinleyError();
718 }
719
720 // sets the location of nodes:
721 void MeshAdapter::setNewX(const escript::Data& new_x)
722 {
723 Finley_Mesh* mesh=m_finleyMesh.get();
724 escriptDataC tmp;
725 const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());
726 if (newDomain!=*this)
727 throw FinleyAdapterException("Error - Illegal domain of new point locations");
728 tmp = new_x.getDataC();
729 Finley_Mesh_setCoordinates(mesh,&tmp);
730 checkFinleyError();
731 }
732
733 // saves a data array in openDX format:
734 void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const
735 {
736 int MAX_namelength=256;
737 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
738 /* win32 refactor */
739 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
740 for(int i=0;i<num_data;i++)
741 {
742 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
743 }
744
745 char* *c_names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
746 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
747 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
748
749 boost::python::list keys=arg.keys();
750 for (int i=0;i<num_data;++i) {
751 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
752 if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
753 throw FinleyAdapterException("Error in saveDX: Data must be defined on same Domain");
754 data[i]=d.getDataC();
755 ptr_data[i]=&(data[i]);
756 std::string n=boost::python::extract<std::string>(keys[i]);
757 c_names[i]=&(names[i][0]);
758 if (n.length()>MAX_namelength-1) {
759 strncpy(c_names[i],n.c_str(),MAX_namelength-1);
760 c_names[i][MAX_namelength-1]='\0';
761 } else {
762 strcpy(c_names[i],n.c_str());
763 }
764 }
765 Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
766 checkFinleyError();
767
768 /* win32 refactor */
769 TMPMEMFREE(c_names);
770 TMPMEMFREE(data);
771 TMPMEMFREE(ptr_data);
772 for(int i=0;i<num_data;i++)
773 {
774 TMPMEMFREE(names[i]);
775 }
776 TMPMEMFREE(names);
777
778 return;
779 }
780
781 // saves a data array in openVTK format:
782 void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const
783 {
784 int MAX_namelength=256;
785 const int num_data=boost::python::extract<int>(arg.attr("__len__")());
786 /* win32 refactor */
787 char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
788 for(int i=0;i<num_data;i++)
789 {
790 names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
791 }
792
793 char* *c_names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
794 escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
795 escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
796
797 boost::python::list keys=arg.keys();
798 for (int i=0;i<num_data;++i) {
799 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
800 if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
801 throw FinleyAdapterException("Error in saveVTK: Data must be defined on same Domain");
802 data[i]=d.getDataC();
803 ptr_data[i]=&(data[i]);
804 std::string n=boost::python::extract<std::string>(keys[i]);
805 c_names[i]=&(names[i][0]);
806 if (n.length()>MAX_namelength-1) {
807 strncpy(c_names[i],n.c_str(),MAX_namelength-1);
808 c_names[i][MAX_namelength-1]='\0';
809 } else {
810 strcpy(c_names[i],n.c_str());
811 }
812 }
813 #ifndef PASO_MPI
814 Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
815 #else
816 Finley_Mesh_saveVTK_MPIO(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
817 #endif
818
819 checkFinleyError();
820 /* win32 refactor */
821 TMPMEMFREE(c_names);
822 TMPMEMFREE(data);
823 TMPMEMFREE(ptr_data);
824 for(int i=0;i<num_data;i++)
825 {
826 TMPMEMFREE(names[i]);
827 }
828 TMPMEMFREE(names);
829
830 return;
831 }
832
833
834 // creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
835 SystemMatrixAdapter MeshAdapter::newSystemMatrix(
836 const int row_blocksize,
837 const escript::FunctionSpace& row_functionspace,
838 const int column_blocksize,
839 const escript::FunctionSpace& column_functionspace,
840 const int type) const
841 {
842 int reduceRowOrder=0;
843 int reduceColOrder=0;
844 // is the domain right?
845 const MeshAdapter& row_domain=dynamic_cast<const MeshAdapter&>(row_functionspace.getDomain());
846 if (row_domain!=*this)
847 throw FinleyAdapterException("Error - domain of row function space does not match the domain of matrix generator.");
848 const MeshAdapter& col_domain=dynamic_cast<const MeshAdapter&>(column_functionspace.getDomain());
849 if (col_domain!=*this)
850 throw FinleyAdapterException("Error - domain of columnn function space does not match the domain of matrix generator.");
851 // is the function space type right
852 if (row_functionspace.getTypeCode()==DegreesOfFreedom) {
853 reduceRowOrder=0;
854 } else if (row_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
855 reduceRowOrder=1;
856 } else {
857 throw FinleyAdapterException("Error - illegal function space type for system matrix rows.");
858 }
859 if (column_functionspace.getTypeCode()==DegreesOfFreedom) {
860 reduceColOrder=0;
861 } else if (column_functionspace.getTypeCode()==ReducedDegreesOfFreedom) {
862 reduceColOrder=1;
863 } else {
864 throw FinleyAdapterException("Error - illegal function space type for system matrix columns.");
865 }
866 // generate matrix:
867
868 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);
869 checkFinleyError();
870 Paso_SystemMatrix* fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);
871 checkPasoError();
872 Paso_SystemMatrixPattern_dealloc(fsystemMatrixPattern);
873 return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
874 }
875
876 //
877 // vtkObject MeshAdapter::createVtkObject() const
878 // TODO:
879 //
880
881 //
882 // returns true if data at the atom_type is considered as being cell centered:
883 bool MeshAdapter::isCellOriented(int functionSpaceCode) const
884 {
885 switch(functionSpaceCode) {
886 case(Nodes):
887 case(DegreesOfFreedom):
888 case(ReducedDegreesOfFreedom):
889 return false;
890 break;
891 case(Elements):
892 case(FaceElements):
893 case(Points):
894 case(ContactElementsZero):
895 case(ContactElementsOne):
896 return true;
897 break;
898 default:
899 stringstream temp;
900 temp << "Error - Cell: Finley does not know anything about function space type " << functionSpaceCode;
901 throw FinleyAdapterException(temp.str());
902 break;
903 }
904 checkFinleyError();
905 return false;
906 }
907
908 bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const
909 {
910 switch(functionSpaceType_source) {
911 case(Nodes):
912 switch(functionSpaceType_target) {
913 case(Nodes):
914 case(ReducedDegreesOfFreedom):
915 case(DegreesOfFreedom):
916 case(Elements):
917 case(FaceElements):
918 case(Points):
919 case(ContactElementsZero):
920 case(ContactElementsOne):
921 return true;
922 default:
923 stringstream temp;
924 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
925 throw FinleyAdapterException(temp.str());
926 }
927 break;
928 case(Elements):
929 if (functionSpaceType_target==Elements) {
930 return true;
931 } else {
932 return false;
933 }
934 case(FaceElements):
935 if (functionSpaceType_target==FaceElements) {
936 return true;
937 } else {
938 return false;
939 }
940 case(Points):
941 if (functionSpaceType_target==Points) {
942 return true;
943 } else {
944 return false;
945 }
946 case(ContactElementsZero):
947 case(ContactElementsOne):
948 if (functionSpaceType_target==ContactElementsZero || functionSpaceType_target==ContactElementsOne) {
949 return true;
950 } else {
951 return false;
952 }
953 case(DegreesOfFreedom):
954 switch(functionSpaceType_target) {
955 case(ReducedDegreesOfFreedom):
956 case(DegreesOfFreedom):
957 case(Nodes):
958 case(Elements):
959 case(FaceElements):
960 case(Points):
961 case(ContactElementsZero):
962 case(ContactElementsOne):
963 return true;
964 default:
965 stringstream temp;
966 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
967 throw FinleyAdapterException(temp.str());
968 }
969 break;
970 case(ReducedDegreesOfFreedom):
971 switch(functionSpaceType_target) {
972 case(ReducedDegreesOfFreedom):
973 case(Elements):
974 case(FaceElements):
975 case(Points):
976 case(ContactElementsZero):
977 case(ContactElementsOne):
978 return true;
979 case(Nodes):
980 case(DegreesOfFreedom):
981 return false;
982 default:
983 stringstream temp;
984 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
985 throw FinleyAdapterException(temp.str());
986 }
987 break;
988 default:
989 stringstream temp;
990 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;
991 throw FinleyAdapterException(temp.str());
992 break;
993 }
994 checkFinleyError();
995 return false;
996 }
997
998 bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const
999 {
1000 return false;
1001 }
1002
1003 bool MeshAdapter::operator==(const AbstractDomain& other) const
1004 {
1005 const MeshAdapter* temp=dynamic_cast<const MeshAdapter*>(&other);
1006 if (temp!=0) {
1007 return (m_finleyMesh==temp->m_finleyMesh);
1008 } else {
1009 return false;
1010 }
1011 }
1012
1013 bool MeshAdapter::operator!=(const AbstractDomain& other) const
1014 {
1015 return !(operator==(other));
1016 }
1017
1018 int MeshAdapter::getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const
1019 {
1020 int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
1021 checkPasoError();
1022 return out;
1023 }
1024
1025 escript::Data MeshAdapter::getX() const
1026 {
1027 return continuousFunction(asAbstractContinuousDomain()).getX();
1028 }
1029
1030 escript::Data MeshAdapter::getNormal() const
1031 {
1032 return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
1033 }
1034
1035 escript::Data MeshAdapter::getSize() const
1036 {
1037 return function(asAbstractContinuousDomain()).getSize();
1038 }
1039
1040 int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
1041 {
1042 int out=0;
1043 Finley_Mesh* mesh=m_finleyMesh.get();
1044 switch (functionSpaceType) {
1045 case(Nodes):
1046 out=mesh->Nodes->Tag[sampleNo];
1047 break;
1048 case(Elements):
1049 out=mesh->Elements->Tag[sampleNo];
1050 break;
1051 case(FaceElements):
1052 out=mesh->FaceElements->Tag[sampleNo];
1053 break;
1054 case(Points):
1055 out=mesh->Points->Tag[sampleNo];
1056 break;
1057 case(ContactElementsZero):
1058 out=mesh->ContactElements->Tag[sampleNo];
1059 break;
1060 case(ContactElementsOne):
1061 out=mesh->ContactElements->Tag[sampleNo];
1062 break;
1063 case(DegreesOfFreedom):
1064 break;
1065 case(ReducedDegreesOfFreedom):
1066 break;
1067 default:
1068 stringstream temp;
1069 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1070 throw FinleyAdapterException(temp.str());
1071 break;
1072 }
1073 return out;
1074 }
1075 int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1076 {
1077 int *out=0,i;
1078 Finley_Mesh* mesh=m_finleyMesh.get();
1079 switch (functionSpaceType) {
1080 case(Nodes):
1081 if (mesh->Nodes!=NULL) {
1082 out=mesh->Nodes->Id;
1083 break;
1084 }
1085 case(Elements):
1086 out=mesh->Elements->Id;
1087 break;
1088 case(FaceElements):
1089 out=mesh->FaceElements->Id;
1090 break;
1091 case(Points):
1092 out=mesh->Points->Id;
1093 break;
1094 case(ContactElementsZero):
1095 out=mesh->ContactElements->Id;
1096 break;
1097 case(ContactElementsOne):
1098 out=mesh->ContactElements->Id;
1099 break;
1100 case(DegreesOfFreedom):
1101 out=mesh->Nodes->degreeOfFreedomId;
1102 break;
1103 case(ReducedDegreesOfFreedom):
1104 out=mesh->Nodes->reducedDegreeOfFreedomId;
1105 break;
1106 default:
1107 stringstream temp;
1108 temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1109 throw FinleyAdapterException(temp.str());
1110 break;
1111 }
1112 return out;
1113 }
1114
1115 void MeshAdapter::setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const
1116 {
1117 Finley_Mesh* mesh=m_finleyMesh.get();
1118 escriptDataC tmp=mask.getDataC();
1119 switch(functionSpaceType) {
1120 case(Nodes):
1121 Finley_NodeFile_setTags(mesh->Nodes,newTag,&tmp);
1122 break;
1123 case(DegreesOfFreedom):
1124 throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
1125 break;
1126 case(ReducedDegreesOfFreedom):
1127 throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
1128 break;
1129 case(Elements):
1130 Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1131 break;
1132 case(FaceElements):
1133 Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
1134 break;
1135 case(Points):
1136 Finley_ElementFile_setTags(mesh->Points,newTag,&tmp);
1137 break;
1138 case(ContactElementsZero):
1139 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1140 break;
1141 case(ContactElementsOne):
1142 Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1143 break;
1144 default:
1145 stringstream temp;
1146 temp << "Error - Finley does not know anything about function space type " << functionSpaceType;
1147 throw FinleyAdapterException(temp.str());
1148 }
1149 checkFinleyError();
1150 return;
1151 }
1152
1153 void MeshAdapter::setTagMap(const std::string& name, int tag)
1154 {
1155 Finley_Mesh* mesh=m_finleyMesh.get();
1156 Finley_Mesh_addTagMap(mesh, name.c_str(),tag);
1157 checkPasoError();
1158 // throwStandardException("MeshAdapter::set TagMap is not implemented.");
1159 }
1160
1161 int MeshAdapter::getTag(const std::string& name) const
1162 {
1163 Finley_Mesh* mesh=m_finleyMesh.get();
1164 int tag=0;
1165 tag=Finley_Mesh_getTag(mesh, name.c_str());
1166 checkPasoError();
1167 // throwStandardException("MeshAdapter::getTag is not implemented.");
1168 return tag;
1169 }
1170
1171 bool MeshAdapter::isValidTagName(const std::string& name) const
1172 {
1173 Finley_Mesh* mesh=m_finleyMesh.get();
1174 return Finley_Mesh_isValidTagName(mesh,name.c_str());
1175 }
1176
1177 std::string MeshAdapter::showTagNames() const
1178 {
1179 stringstream temp;
1180 Finley_Mesh* mesh=m_finleyMesh.get();
1181 Finley_TagMap* tag_map=mesh->TagMap;
1182 while (tag_map) {
1183 temp << tag_map->name;
1184 tag_map=tag_map->next;
1185 if (tag_map) temp << ", ";
1186 }
1187 return temp.str();
1188 }
1189
1190 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26