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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1059 - (hide 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 jgs 102 // $Id$
2 jgs 82 /*
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 jgs 472
16 ksteube 817 #ifdef PASO_MPI
17     #include <mpi.h>
18     #endif
19 jgs 203 #include "MeshAdapter.h"
20 jgs 149
21 robwdcock 682 #include "escript/Data.h"
22     #include "escript/DataFactory.h"
23 jgs 480
24 jgs 82 using namespace std;
25     using namespace escript;
26    
27     namespace finley {
28 jgs 149
29 jgs 82 //
30 jgs 149 // define the static constants
31 jgs 82 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 gross 1059 const int MeshAdapter::ReducedElements=FINLEY_REDUCED_ELEMENTS;
37 jgs 82 const int MeshAdapter::FaceElements=FINLEY_FACE_ELEMENTS;
38 gross 1059 const int MeshAdapter::ReducedFaceElements=FINLEY_REDUCED_FACE_ELEMENTS;
39 jgs 82 const int MeshAdapter::Points=FINLEY_POINTS;
40     const int MeshAdapter::ContactElementsZero=FINLEY_CONTACT_ELEMENTS_1;
41 gross 1059 const int MeshAdapter::ReducedContactElementsZero=FINLEY_REDUCED_CONTACT_ELEMENTS_1;
42 jgs 82 const int MeshAdapter::ContactElementsOne=FINLEY_CONTACT_ELEMENTS_2;
43 gross 1059 const int MeshAdapter::ReducedContactElementsOne=FINLEY_REDUCED_CONTACT_ELEMENTS_2;
44 jgs 82
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 jgs 149
54 jgs 82 //
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 woo409 757 char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
81 jgs 82 strcpy(fName,fileName.c_str());
82     Finley_Mesh_write(m_finleyMesh.get(),fName);
83     checkFinleyError();
84 woo409 757 TMPMEMFREE(fName);
85 jgs 82 }
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 gross 1059 (FunctionSpaceNamesMapType::value_type(ReducedElements,"Finley_Reduced_Elements"));
122     m_functionSpaceTypeNames.insert
123 jgs 82 (FunctionSpaceNamesMapType::value_type(FaceElements,"Finley_Face_Elements"));
124     m_functionSpaceTypeNames.insert
125 gross 1059 (FunctionSpaceNamesMapType::value_type(ReducedFaceElements,"Finley_Reduced_Face_Elements"));
126     m_functionSpaceTypeNames.insert
127 jgs 82 (FunctionSpaceNamesMapType::value_type(Points,"Finley_Points"));
128     m_functionSpaceTypeNames.insert
129     (FunctionSpaceNamesMapType::value_type(ContactElementsZero,"Finley_Contact_Elements_0"));
130     m_functionSpaceTypeNames.insert
131 gross 1059 (FunctionSpaceNamesMapType::value_type(ReducedContactElementsZero,"Finley_Reduced_Contact_Elements_0"));
132     m_functionSpaceTypeNames.insert
133 jgs 82 (FunctionSpaceNamesMapType::value_type(ContactElementsOne,"Finley_Contact_Elements_1"));
134 gross 1059 m_functionSpaceTypeNames.insert
135     (FunctionSpaceNamesMapType::value_type(ReducedContactElementsOne,"Finley_Reduced_Contact_Elements_1"));
136 jgs 82 }
137    
138     int MeshAdapter::getContinuousFunctionCode() const
139     {
140     return Nodes;
141     }
142 jgs 149
143 jgs 82 int MeshAdapter::getFunctionCode() const
144     {
145     return Elements;
146     }
147 gross 1059 int MeshAdapter::getReducedFunctionCode() const
148     {
149     return ReducedElements;
150     }
151 jgs 149
152 jgs 82 int MeshAdapter::getFunctionOnBoundaryCode() const
153     {
154     return FaceElements;
155     }
156 gross 1059 int MeshAdapter::getReducedFunctionOnBoundaryCode() const
157     {
158     return ReducedFaceElements;
159     }
160 jgs 149
161 jgs 82 int MeshAdapter::getFunctionOnContactZeroCode() const
162     {
163     return ContactElementsZero;
164     }
165 gross 1059 int MeshAdapter::getReducedFunctionOnContactZeroCode() const
166     {
167     return ReducedContactElementsZero;
168     }
169 jgs 149
170 jgs 82 int MeshAdapter::getFunctionOnContactOneCode() const
171     {
172     return ContactElementsOne;
173     }
174 gross 1059 int MeshAdapter::getReducedFunctionOnContactOneCode() const
175     {
176     return ReducedContactElementsOne;
177     }
178 jgs 82
179     int MeshAdapter::getSolutionCode() const
180     {
181     return DegreesOfFreedom;
182     }
183 jgs 149
184 jgs 82 int MeshAdapter::getReducedSolutionCode() const
185     {
186     return ReducedDegreesOfFreedom;
187     }
188 jgs 149
189 jgs 82 int MeshAdapter::getDiracDeltaFunctionCode() const
190     {
191     return Points;
192     }
193 jgs 149
194 jgs 82 //
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 jgs 149
204 jgs 82 //
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 bcumming 751 #ifndef PASO_MPI
252 jgs 82 numSamples=mesh->Nodes->numDegreesOfFreedom;
253 bcumming 751 #else
254     numSamples=mesh->Nodes->degreeOfFreedomDistribution->numLocal;
255     #endif
256 jgs 82 }
257     break;
258     case(ReducedDegreesOfFreedom):
259     if (mesh->Nodes!=NULL) {
260     numDataPointsPerSample=1;
261 bcumming 751 #ifndef PASO_MPI
262 jgs 82 numSamples=mesh->Nodes->reducedNumDegreesOfFreedom;
263 bcumming 751 #else
264     numSamples=mesh->Nodes->reducedDegreeOfFreedomDistribution->numLocal;
265     #endif
266 jgs 82 }
267     break;
268     default:
269     stringstream temp;
270 jgs 150 temp << "Error - Invalid function space type: " << functionSpaceCode << " for domain: " << getDescription();
271 jgs 82 throw FinleyAdapterException(temp.str());
272     break;
273     }
274     return pair<int,int>(numDataPointsPerSample,numSamples);
275     }
276 jgs 149
277 jgs 82 //
278     // adds linear PDE of second order into a given stiffness matrix and right hand side:
279     //
280     void MeshAdapter::addPDEToSystem(
281 jgs 480 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 jgs 82 {
286 gross 798 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 bcumming 751
298 jgs 82 Finley_Mesh* mesh=m_finleyMesh.get();
299 gross 798
300 bcumming 751 Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(), &_rhs, &_A, &_B, &_C, &_D, &_X, &_Y );
301 jgs 82 checkFinleyError();
302 bcumming 751
303 gross 798 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, mat.getPaso_SystemMatrix(), &_rhs, 0, 0, 0, &_d, 0, &_y );
304 jgs 82 checkFinleyError();
305 bcumming 751
306 gross 798 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, mat.getPaso_SystemMatrix(), &_rhs , 0, 0, 0, &_d_contact, 0, &_y_contact );
307 jgs 82 checkFinleyError();
308     }
309 jgs 149
310 jgs 82 //
311 jgs 102 // adds linear PDE of second order into the right hand side only
312     //
313 bcumming 751 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 jgs 102 {
315     Finley_Mesh* mesh=m_finleyMesh.get();
316 jgs 147
317 gross 798 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 jgs 102 checkFinleyError();
325 jgs 147
326 gross 798 Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, 0, &_rhs, 0, 0, 0, 0, 0, &_y );
327     checkFinleyError();
328 jgs 147
329 gross 798 Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, 0, &_rhs , 0, 0, 0, 0, 0, &_y_contact );
330 jgs 102 checkFinleyError();
331     }
332 jgs 149
333 jgs 102 //
334 jgs 82 // interpolates data between different function spaces:
335     //
336 jgs 480 void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
337 jgs 82 {
338     const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());
339     const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
340 bcumming 751 if (inDomain!=*this)
341     throw FinleyAdapterException("Error - Illegal domain of interpolant.");
342 jgs 82 if (targetDomain!=*this)
343 bcumming 751 throw FinleyAdapterException("Error - Illegal domain of interpolation target.");
344 jgs 82
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 jgs 150 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 jgs 82 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 jgs 150 throw FinleyAdapterException("Error - No interpolation with data on elements possible.");
381 jgs 82 }
382     break;
383     case(FaceElements):
384     if (target.getFunctionSpace().getTypeCode()==FaceElements) {
385     Finley_Assemble_CopyElementData(mesh->FaceElements,&(target.getDataC()),&(in.getDataC()));
386     } else {
387 jgs 150 throw FinleyAdapterException("Error - No interpolation with data on face elements possible.");
388 jgs 82 break;
389     }
390     case(Points):
391     if (target.getFunctionSpace().getTypeCode()==Points) {
392     Finley_Assemble_CopyElementData(mesh->Points,&(target.getDataC()),&(in.getDataC()));
393     } else {
394 jgs 150 throw FinleyAdapterException("Error - No interpolation with data on points possible.");
395 jgs 82 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 jgs 150 throw FinleyAdapterException("Error - No interpolation with data on contact elements possible.");
404 jgs 82 break;
405     }
406     break;
407 bcumming 751 case(DegreesOfFreedom):
408 jgs 82 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 bcumming 751 #ifndef PASO_MPI
415 jgs 82 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 bcumming 751 #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 jgs 82 default:
457 jgs 150 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 jgs 82 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 jgs 147 case(Nodes):
482 jgs 150 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to mesh nodes.");
483 jgs 147 break;
484 jgs 82 case(DegreesOfFreedom):
485 jgs 150 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");
486 jgs 82 break;
487     default:
488 jgs 150 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 jgs 82 break;
492     }
493     break;
494     default:
495 jgs 150 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 jgs 82 break;
499     }
500     checkFinleyError();
501     }
502    
503     //
504     // copies the locations of sample points into x:
505     //
506 jgs 480 void MeshAdapter::setToX(escript::Data& arg) const
507 jgs 82 {
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 jgs 480 escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
517 jgs 82 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 jgs 149
524 jgs 82 //
525     // return the normal vectors at the location of data points as a Data object:
526     //
527 jgs 480 void MeshAdapter::setToNormal(escript::Data& normal) const
528 jgs 82 {
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 jgs 150 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for nodes");
536 jgs 82 break;
537     case(Elements):
538 jgs 150 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements");
539 jgs 82 break;
540     case (FaceElements):
541     Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&(normal.getDataC()));
542     break;
543     case(Points):
544 jgs 150 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for point elements");
545 jgs 82 break;
546     case (ContactElementsOne):
547     case (ContactElementsZero):
548     Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&(normal.getDataC()));
549     break;
550     case(DegreesOfFreedom):
551 jgs 150 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for degrees of freedom.");
552 jgs 82 break;
553     case(ReducedDegreesOfFreedom):
554 jgs 150 throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced degrees of freedom.");
555 jgs 82 break;
556     default:
557 jgs 150 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 jgs 82 break;
561     }
562     checkFinleyError();
563     }
564 jgs 149
565 jgs 82 //
566     // interpolates data to other domain:
567     //
568 jgs 480 void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
569 jgs 82 {
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 jgs 150 throw FinleyAdapterException("Error - Finley does not allow interpolation across domains yet.");
575 jgs 82 }
576 jgs 149
577 jgs 82 //
578     // calculates the integral of a function defined of arg:
579     //
580 jgs 480 void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
581 jgs 82 {
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 jgs 150 throw FinleyAdapterException("Error - Integral of data on nodes is not supported.");
590 jgs 82 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 jgs 150 throw FinleyAdapterException("Error - Integral of data on points is not supported.");
599 jgs 82 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 jgs 150 throw FinleyAdapterException("Error - Integral of data on degrees of freedom is not supported.");
608 jgs 82 break;
609     case(ReducedDegreesOfFreedom):
610 jgs 150 throw FinleyAdapterException("Error - Integral of data on reduced degrees of freedom is not supported.");
611 jgs 82 break;
612     default:
613 jgs 150 stringstream temp;
614     temp << "Error - Integrals: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
615     throw FinleyAdapterException(temp.str());
616 jgs 82 break;
617     }
618     checkFinleyError();
619     }
620 jgs 149
621 jgs 82 //
622     // calculates the gradient of arg:
623     //
624 jgs 480 void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
625 jgs 82 {
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 bcumming 751 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 jgs 82 switch(grad.getFunctionSpace().getTypeCode()) {
648     case(Nodes):
649 jgs 150 throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
650 jgs 82 break;
651     case(Elements):
652 bcumming 751 Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&nodeDataC);
653 jgs 82 break;
654     case(FaceElements):
655 bcumming 751 Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&nodeDataC);
656 jgs 82 break;
657     case(Points):
658 jgs 150 throw FinleyAdapterException("Error - Gradient at points is not supported.");
659 jgs 82 break;
660     case(ContactElementsZero):
661 bcumming 751 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
662 jgs 82 break;
663     case(ContactElementsOne):
664 bcumming 751 Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
665 jgs 82 break;
666     case(DegreesOfFreedom):
667 jgs 150 throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
668 jgs 82 break;
669     case(ReducedDegreesOfFreedom):
670 jgs 150 throw FinleyAdapterException("Error - Gradient at reduced degrees of freedom is not supported.");
671 jgs 82 break;
672     default:
673 jgs 150 stringstream temp;
674     temp << "Error - Gradient: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
675     throw FinleyAdapterException(temp.str());
676 jgs 82 break;
677     }
678     checkFinleyError();
679     }
680 jgs 149
681 jgs 82 //
682     // returns the size of elements:
683     //
684 jgs 480 void MeshAdapter::setToSize(escript::Data& size) const
685 jgs 82 {
686     Finley_Mesh* mesh=m_finleyMesh.get();
687     escriptDataC tmp=size.getDataC();
688     switch(size.getFunctionSpace().getTypeCode()) {
689     case(Nodes):
690 jgs 150 throw FinleyAdapterException("Error - Size of nodes is not supported.");
691 jgs 82 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 jgs 150 throw FinleyAdapterException("Error - Size of point elements is not supported.");
700 jgs 82 break;
701     case(ContactElementsZero):
702     case(ContactElementsOne):
703     Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
704     break;
705     case(DegreesOfFreedom):
706 jgs 150 throw FinleyAdapterException("Error - Size of degrees of freedom is not supported.");
707 jgs 82 break;
708     case(ReducedDegreesOfFreedom):
709 jgs 150 throw FinleyAdapterException("Error - Size of reduced degrees of freedom is not supported.");
710 jgs 82 break;
711     default:
712 jgs 150 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 jgs 82 break;
716     }
717     checkFinleyError();
718     }
719 jgs 149
720 jgs 82 // sets the location of nodes:
721 jgs 480 void MeshAdapter::setNewX(const escript::Data& new_x)
722 jgs 82 {
723     Finley_Mesh* mesh=m_finleyMesh.get();
724 bcumming 751 escriptDataC tmp;
725 jgs 82 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 bcumming 751 tmp = new_x.getDataC();
729     Finley_Mesh_setCoordinates(mesh,&tmp);
730 jgs 82 checkFinleyError();
731     }
732 jgs 149
733 jgs 82 // saves a data array in openDX format:
734 jgs 153 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 woo409 757 /* 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 jgs 153
745 woo409 757 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 jgs 153 boost::python::list keys=arg.keys();
750     for (int i=0;i<num_data;++i) {
751 jgs 480 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
752 jgs 153 if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
753 dhawcroft 793 throw FinleyAdapterException("Error in saveDX: Data must be defined on same Domain");
754 jgs 153 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 woo409 757
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 jgs 153 return;
779 jgs 82 }
780 jgs 149
781 jgs 110 // saves a data array in openVTK format:
782 jgs 153 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 woo409 757 /* 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 jgs 153
793 woo409 757 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 jgs 153 boost::python::list keys=arg.keys();
798     for (int i=0;i<num_data;++i) {
799 jgs 480 escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
800 jgs 153 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 dhawcroft 793 #ifndef PASO_MPI
814 jgs 153 Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
815 dhawcroft 793 #else
816     Finley_Mesh_saveVTK_MPIO(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
817     #endif
818    
819     checkFinleyError();
820 woo409 757 /* 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 jgs 153 return;
831 jgs 110 }
832 jgs 153
833    
834 jgs 82 // 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 jgs 102 const int type) const
841 jgs 82 {
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 jgs 102
868 jgs 150 Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);
869 jgs 82 checkFinleyError();
870 ksteube 971 Paso_SystemMatrix* fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);
871 jgs 150 checkPasoError();
872     Paso_SystemMatrixPattern_dealloc(fsystemMatrixPattern);
873 jgs 82 return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
874     }
875 jgs 149
876 jgs 82 //
877     // vtkObject MeshAdapter::createVtkObject() const
878     // TODO:
879     //
880 jgs 149
881 jgs 82 //
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 jgs 150 stringstream temp;
900     temp << "Error - Cell: Finley does not know anything about function space type " << functionSpaceCode;
901     throw FinleyAdapterException(temp.str());
902 jgs 82 break;
903     }
904     checkFinleyError();
905     return false;
906     }
907 jgs 149
908 jgs 82 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 jgs 150 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 jgs 82 }
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 jgs 150 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 jgs 82 }
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 jgs 147 case(Nodes):
980 jgs 82 case(DegreesOfFreedom):
981     return false;
982     default:
983 jgs 150 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 jgs 82 }
987     break;
988     default:
989 jgs 150 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 jgs 82 break;
993     }
994     checkFinleyError();
995     return false;
996     }
997 jgs 149
998 jgs 82 bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const
999     {
1000     return false;
1001     }
1002 jgs 104
1003 jgs 121 bool MeshAdapter::operator==(const AbstractDomain& other) const
1004 jgs 82 {
1005 jgs 121 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 jgs 82 }
1012    
1013 jgs 121 bool MeshAdapter::operator!=(const AbstractDomain& other) const
1014 jgs 104 {
1015 jgs 121 return !(operator==(other));
1016 jgs 104 }
1017    
1018 jgs 150 int MeshAdapter::getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const
1019 jgs 102 {
1020 jgs 150 int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
1021     checkPasoError();
1022 jgs 102 return out;
1023     }
1024 jgs 149
1025 jgs 480 escript::Data MeshAdapter::getX() const
1026 jgs 102 {
1027     return continuousFunction(asAbstractContinuousDomain()).getX();
1028     }
1029 jgs 149
1030 jgs 480 escript::Data MeshAdapter::getNormal() const
1031 jgs 102 {
1032     return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
1033     }
1034 jgs 149
1035 jgs 480 escript::Data MeshAdapter::getSize() const
1036 jgs 102 {
1037     return function(asAbstractContinuousDomain()).getSize();
1038     }
1039    
1040 jgs 110 int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
1041     {
1042 gross 547 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 jgs 110 }
1075 gross 964 int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1076 jgs 110 {
1077 gross 964 int *out=0,i;
1078 gross 547 Finley_Mesh* mesh=m_finleyMesh.get();
1079     switch (functionSpaceType) {
1080     case(Nodes):
1081     if (mesh->Nodes!=NULL) {
1082 gross 964 out=mesh->Nodes->Id;
1083 gross 547 break;
1084     }
1085     case(Elements):
1086 gross 964 out=mesh->Elements->Id;
1087 gross 547 break;
1088     case(FaceElements):
1089 gross 964 out=mesh->FaceElements->Id;
1090 gross 547 break;
1091     case(Points):
1092 gross 964 out=mesh->Points->Id;
1093 gross 547 break;
1094     case(ContactElementsZero):
1095 gross 964 out=mesh->ContactElements->Id;
1096 gross 547 break;
1097     case(ContactElementsOne):
1098 gross 964 out=mesh->ContactElements->Id;
1099 gross 547 break;
1100     case(DegreesOfFreedom):
1101 gross 964 out=mesh->Nodes->degreeOfFreedomId;
1102 gross 547 break;
1103     case(ReducedDegreesOfFreedom):
1104 gross 964 out=mesh->Nodes->reducedDegreeOfFreedomId;
1105 gross 547 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 jgs 110 }
1114    
1115 gross 767 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 gross 1044 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 gross 767
1161 gross 1044 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 jgs 82 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26