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

Annotation of /branches/diaplayground/finley/src/CPPAdapter/MeshAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 757 - (hide annotations)
Mon Jun 26 13:12:56 2006 UTC (12 years, 9 months ago) by woo409
Original Path: trunk/finley/src/CPPAdapter/MeshAdapter.cpp
File size: 41135 byte(s)
+ Merge of intelc_win32 branch (revision 741:755) with trunk. Tested on iVEC altix (run_tests and py_tests all pass)

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