/[escript]/trunk/dataexporter/src/ElementData.cpp
ViewVC logotype

Annotation of /trunk/dataexporter/src/ElementData.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2201 - (hide annotations)
Fri Jan 9 00:45:44 2009 UTC (12 years, 3 months ago) by caltinay
Original Path: trunk/tools/libescriptreader/src/escriptreader/ElementData.cpp
File size: 23177 byte(s)
Avoid NULL reference.

1 caltinay 2183
2     /*******************************************************
3     *
4     * Copyright (c) 2003-2008 by University of Queensland
5     * Earth Systems Science Computational Center (ESSCC)
6     * http://www.uq.edu.au/esscc
7     *
8     * Primary Business: Queensland, Australia
9     * Licensed under the Open Software License version 3.0
10     * http://www.opensource.org/licenses/osl-3.0.php
11     *
12     *******************************************************/
13    
14     //
15     // ElementData.cpp
16     //
17     #include <escriptreader/ElementData.h>
18     #include <netcdf.hh>
19     #if HAVE_SILO
20     #include <silo.h>
21     #endif
22    
23     using namespace std;
24    
25     // The following arrays contain indices to convert unsupported element
26     // types into supported ones (lines, triangles, quads, hexahedrons)
27    
28     static const size_t line3indices[2*2] = {
29     0, 2,
30     2, 1
31     };
32     static const size_t tri6indices[4*3] = {
33     0, 3, 5,
34     3, 1, 4,
35     5, 4, 2,
36     3, 4, 5
37     };
38     static const size_t rec8indices[6*3] = {
39     0, 4, 7,
40     4, 1, 5,
41     5, 2, 6,
42     6, 3, 7,
43     4, 6, 7,
44     4, 5, 6
45     };
46     static const size_t rec9indices[4*4] = {
47     0, 4, 8, 7,
48     4, 1, 5, 8,
49     7, 8, 6, 3,
50     8, 5, 2, 6
51     };
52     static const size_t tet10indices[8*4] = {
53     6, 0, 7, 4,
54     2, 6, 9, 5,
55     9, 7, 3, 8,
56     4, 1, 8, 5,
57     6, 7, 9, 8,
58     5, 6, 9, 8,
59     6, 5, 4, 8,
60     7, 6, 4, 8
61     };
62     static const size_t hex20indices[36*3] = {
63     0, 8, 12, 8, 1, 13, 13, 5, 16,
64     16, 4, 12, 8, 13, 16, 8, 16, 12,
65     1, 9, 13, 9, 2, 14, 14, 6, 17,
66     17, 5, 13, 9, 14, 17, 9, 17, 13,
67     2, 10, 14, 10, 3, 15, 15, 7, 18,
68     18, 14, 6, 10, 15, 18, 10, 18, 14,
69     3, 11, 15, 11, 0, 12, 12, 4, 19,
70     19, 7, 15, 11, 12, 19, 11, 19, 15,
71     4, 16, 19, 16, 5, 17, 17, 6, 18,
72     18, 7, 19, 16, 17, 18, 16, 18, 19,
73     3, 10, 11, 10, 2, 9, 9, 1, 8,
74     8, 0, 11, 10, 9, 8, 10, 8, 11
75     };
76     static const size_t hex27indices[8*8] = {
77     0, 8, 20, 11, 12, 21, 26, 24,
78     8, 1, 9, 20, 21, 13, 22, 26,
79     11, 20, 10, 3, 24, 26, 23, 15,
80     20, 9, 2, 10, 26, 22, 14, 23,
81     12, 21, 26, 24, 4, 16, 25, 19,
82     21, 13, 22, 26, 16, 5, 17, 25,
83     24, 26, 23, 15, 19, 25, 18, 7,
84     26, 22, 14, 23, 25, 17, 6, 18
85     };
86    
87 caltinay 2187 namespace EscriptReader {
88    
89 caltinay 2183 //
90     // Constructor
91     //
92     ElementData::ElementData(const string& elementName, const Mesh* mainMesh)
93     : name(elementName), count(0), reducedCount(0), numGhostElements(0),
94     numReducedGhostElements(0), fullMesh(NULL), reducedMesh(NULL),
95     originalMesh(mainMesh), fullMeshIsOriginalMesh(false)
96     {
97     }
98    
99     //
100     // Copy constructor
101     //
102     ElementData::ElementData(const ElementData& e)
103     {
104     name = e.name;
105     count = e.count;
106     reducedCount = e.reducedCount;
107     numGhostElements = e.numGhostElements;
108     numReducedGhostElements = e.numReducedGhostElements;
109     numDims = e.numDims;
110     type = e.type;
111     reducedType = e.reducedType;
112     nodesPerElement = e.nodesPerElement;
113     reducedNodesPerElement = e.reducedNodesPerElement;
114    
115     originalMesh = e.originalMesh;
116     fullMeshIsOriginalMesh = e.fullMeshIsOriginalMesh;
117     if (fullMeshIsOriginalMesh)
118     fullMesh = const_cast<Mesh*>(originalMesh);
119     else if (e.fullMesh)
120     fullMesh = new Mesh(*e.fullMesh);
121     else
122     fullMesh = NULL;
123    
124     if (e.reducedMesh)
125     reducedMesh = new Mesh(*e.reducedMesh);
126     else
127     reducedMesh = NULL;
128    
129     nodes = e.nodes;
130     reducedNodes = e.reducedNodes;
131     ID = e.ID;
132     color = e.color;
133     tag = e.tag;
134     owner = e.owner;
135     reducedOwner = e.reducedOwner;
136     indexArray = e.indexArray;
137     reducedIndexArray = e.reducedIndexArray;
138     ID2idx = e.ID2idx;
139     }
140    
141     //
142     // Destructor
143     //
144     ElementData::~ElementData()
145     {
146     if (reducedMesh)
147     delete reducedMesh;
148     if (fullMesh && !fullMeshIsOriginalMesh)
149     delete fullMesh;
150     }
151    
152     StringVec ElementData::getMeshNames() const
153     {
154     StringVec res;
155     if (fullMesh && !fullMeshIsOriginalMesh)
156     res.push_back(fullMesh->getName());
157     if (reducedMesh)
158     res.push_back(reducedMesh->getName());
159     return res;
160     }
161    
162     StringVec ElementData::getVarNames() const
163     {
164     StringVec res;
165     if (count > 0) {
166     res.push_back(name + string("_Color"));
167     res.push_back(name + string("_Id"));
168     res.push_back(name + string("_Owner"));
169     res.push_back(name + string("_Tag"));
170     }
171     return res;
172     }
173    
174     const IntVec& ElementData::getVarDataByName(const string varName) const
175     {
176     if (varName == name+string("_Color"))
177     return color;
178     else if (varName == name+string("_Id"))
179     return ID;
180     else if (varName == name+string("_Owner")) {
181     if (reducedCount > 0)
182     return reducedOwner;
183     else
184     return owner;
185     } else if (varName == name+string("_Tag"))
186     return tag;
187     else
188 caltinay 2201 throw "Invalid variable name";
189 caltinay 2183 }
190    
191     //
192     // Reads element data from given NetCDF file
193     //
194     bool ElementData::readFromNc(NcFile* ncfile)
195     {
196     NcAtt* att;
197     NcVar* var;
198    
199     string num_str("num_");
200     num_str += name;
201    
202     att = ncfile->get_att(num_str.c_str());
203     count = att->as_int(0);
204    
205     // Only attempt to read data if there are any elements.
206     // Having no elements is not an error.
207     if (count == 0)
208     return true;
209    
210     att = ncfile->get_att((num_str + string("_numNodes")).c_str());
211     nodesPerElement = att->as_int(0);
212    
213     nodes.insert(nodes.end(), count*nodesPerElement, 0);
214     var = ncfile->get_var((name + string("_Nodes")).c_str());
215     var->get(&nodes[0], count, nodesPerElement);
216    
217     color.insert(color.end(), count, 0);
218     var = ncfile->get_var((name + string("_Color")).c_str());
219     var->get(&color[0], count);
220    
221     ID.insert(ID.end(), count, 0);
222     var = ncfile->get_var((name + string("_Id")).c_str());
223     var->get(&ID[0], count);
224    
225     owner.insert(owner.end(), count, 0);
226     var = ncfile->get_var((name + string("_Owner")).c_str());
227     var->get(&owner[0], count);
228    
229     tag.insert(tag.end(), count, 0);
230     var = ncfile->get_var((name + string("_Tag")).c_str());
231     var->get(&tag[0], count);
232    
233     att = ncfile->get_att((name + string("_TypeId")).c_str());
234     FinleyElementInfo f = getFinleyTypeInfo((ElementTypeId)att->as_int(0));
235     type = f.elementType;
236     reducedType = f.reducedElementType;
237    
238     // build elementID->index map
239     buildIndexMap();
240    
241     if (f.elementFactor > 1 || f.reducedElementSize != nodesPerElement)
242     buildReducedElements(f);
243    
244     buildMeshes();
245    
246     return true;
247     }
248    
249     //
250     //
251     //
252     void ElementData::reorderArray(IntVec& v, const IntVec& idx,
253     int elementsPerIndex)
254     {
255     IntVec newArray(v.size());
256     IntVec::iterator arrIt = newArray.begin();
257     IntVec::const_iterator idxIt;
258     if (elementsPerIndex == 1) {
259     for (idxIt=idx.begin(); idxIt!=idx.end(); idxIt++) {
260     *arrIt++ = v[*idxIt];
261     }
262     } else {
263     for (idxIt=idx.begin(); idxIt!=idx.end(); idxIt++) {
264     int i = *idxIt;
265     copy(&v[i*elementsPerIndex], &v[(i+1)*elementsPerIndex], arrIt);
266     arrIt += elementsPerIndex;
267     }
268     }
269     v.swap(newArray);
270     }
271    
272     //
273     //
274     //
275     void ElementData::buildReducedElements(const FinleyElementInfo& f)
276     {
277     reducedNodes.clear();
278     reducedNodes.insert(reducedNodes.end(),
279     f.reducedElementSize*count, 0);
280     IntVec::iterator reducedIt = reducedNodes.begin();
281     IntVec::const_iterator origIt;
282     for (origIt=nodes.begin(); origIt!=nodes.end();
283     origIt+=nodesPerElement)
284     {
285     std::copy(origIt, origIt+f.reducedElementSize, reducedIt);
286     reducedIt += f.reducedElementSize;
287     }
288    
289     // Remove comment to save disk space - we don't really need the full
290     // elements except for the main mesh
291     if (f.elementFactor > 1 /*&& name == "Elements"*/) {
292     // replace each element by multiple smaller ones
293     IntVec fullNodes(f.elementSize*f.elementFactor*count);
294     IntVec::iterator cellIt = fullNodes.begin();
295    
296     // copy owner data
297     owner.swap(reducedOwner);
298     owner.clear();
299     for (int i=0; i < count; i++) {
300     owner.insert(owner.end(), f.elementFactor, reducedOwner[i]);
301     for (int j=0; j < f.elementFactor*f.elementSize; j++)
302     *cellIt++ = nodes[
303     nodesPerElement*i+f.multiCellIndices[j]];
304     }
305     nodes.swap(fullNodes);
306     reducedCount = count;
307     reducedNodesPerElement = f.reducedElementSize;
308     nodesPerElement = f.elementSize;
309     count *= f.elementFactor;
310     } else {
311     // we only keep the reduced elements but treat them as regular
312     // ones, so replace the data accordingly
313     nodes.swap(reducedNodes);
314     reducedNodes.clear();
315     nodesPerElement = f.reducedElementSize;
316     type = f.reducedElementType;
317     }
318     }
319    
320     //
321     //
322     //
323     void ElementData::prepareGhostIndices(int ownIndex)
324     {
325     indexArray.clear();
326     reducedIndexArray.clear();
327     numGhostElements = 0;
328     numReducedGhostElements = 0;
329    
330     // move indices of "ghost zones" to the end to be able to reorder
331     // data accordingly
332     for (int i=0; i<count; i++)
333     if (owner[i] == ownIndex)
334     indexArray.push_back(i);
335    
336     for (int i=0; i<count; i++)
337     if (owner[i] != ownIndex) {
338     numGhostElements++;
339     indexArray.push_back(i);
340     }
341    
342     for (int i=0; i<reducedCount; i++)
343     if (reducedOwner[i] == ownIndex)
344     reducedIndexArray.push_back(i);
345    
346     for (int i=0; i<reducedCount; i++)
347     if (reducedOwner[i] != ownIndex) {
348     numReducedGhostElements++;
349     reducedIndexArray.push_back(i);
350     }
351     }
352    
353     //
354     //
355     //
356     void ElementData::handleGhostZones(int ownIndex)
357     {
358     prepareGhostIndices(ownIndex);
359    
360     // move "ghost data" to the end of the arrays
361     if (numGhostElements > 0) {
362     reorderArray(nodes, indexArray, nodesPerElement);
363     reorderArray(owner, indexArray, 1);
364     if (reducedCount == 0) {
365     reorderArray(color, indexArray, 1);
366     reorderArray(ID, indexArray, 1);
367     reorderArray(tag, indexArray, 1);
368     }
369     }
370    
371     if (numReducedGhostElements > 0) {
372     reorderArray(reducedNodes, reducedIndexArray, reducedNodesPerElement);
373     reorderArray(reducedOwner, reducedIndexArray, 1);
374     reorderArray(color, reducedIndexArray, 1);
375     reorderArray(ID, reducedIndexArray, 1);
376     reorderArray(tag, reducedIndexArray, 1);
377     }
378     }
379    
380     //
381     //
382     //
383     void ElementData::removeGhostZones()
384     {
385     if (numGhostElements > 0) {
386     count -= numGhostElements;
387     nodes.resize(count*nodesPerElement);
388     owner.resize(count);
389     if (reducedCount == 0) {
390     color.resize(count);
391     ID.resize(count);
392     tag.resize(count);
393     }
394     numGhostElements = 0;
395     }
396    
397     if (numReducedGhostElements > 0) {
398     reducedCount -= numReducedGhostElements;
399     reducedNodes.resize(reducedCount*reducedNodesPerElement);
400     reducedOwner.resize(reducedCount);
401     color.resize(reducedCount);
402     ID.resize(reducedCount);
403     tag.resize(reducedCount);
404     numReducedGhostElements = 0;
405     }
406     }
407    
408     //
409     //
410     //
411     void ElementData::buildMeshes()
412     {
413     if (count == 0)
414     return;
415    
416     // use existing original mesh for Elements but build a new mesh
417     // containing only the required nodes for other types
418     if (name == "Elements") {
419     fullMesh = const_cast<Mesh*>(originalMesh);
420     fullMeshIsOriginalMesh = true;
421     } else {
422     // first: build a map of required IDs while translating
423     // the original node IDs at the same time
424     IntVec::iterator it;
425     IndexMap nodeID2idx;
426     IntVec nodeIDs;
427     size_t newIdx = 0;
428    
429     for (it = nodes.begin(); it != nodes.end(); it++) {
430     IndexMap::iterator res = nodeID2idx.find(*it);
431     if (res == nodeID2idx.end()) {
432     nodeIDs.push_back(*it);
433     nodeID2idx[*it] = newIdx;
434     *it = newIdx++;
435     } else {
436     *it = res->second;
437     }
438     }
439    
440     // second: use map to fill coordinates for the nodes
441     // newIdx contains the number of nodes
442     CoordArray coords;
443     const CoordArray& origCoords = originalMesh->getCoords();
444     for (int dim=0; dim < originalMesh->getNumDims(); dim++) {
445     float* c = new float[newIdx];
446     coords.push_back(c);
447     IndexMap::const_iterator mIt;
448     for (mIt = nodeID2idx.begin(); mIt != nodeID2idx.end(); mIt++)
449     c[mIt->second] = origCoords[dim][mIt->first];
450     }
451    
452     if (fullMesh)
453     delete fullMesh;
454    
455     fullMesh = new Mesh(coords, originalMesh->getNumDims(), newIdx);
456     fullMesh->setName(name);
457     fullMesh->setIndexMap(nodeID2idx);
458     fullMesh->setNodeIDs(nodeIDs);
459     }
460    
461     #ifdef _DEBUG
462     cout << fullMesh->getName() << " has " << fullMesh->getNumNodes()
463     << " nodes, " << count << " elements" << endl;
464     #endif
465    
466     // build a reduced mesh if necessary
467     if (reducedCount > 0) {
468     IntVec::iterator it;
469     IndexMap nodeID2idx;
470     IntVec reducedNodeIDs;
471     size_t newIdx = 0;
472    
473     for (it = reducedNodes.begin(); it != reducedNodes.end(); it++) {
474     int id = originalMesh->getNodeIDs()[*it];
475     IndexMap::iterator res = nodeID2idx.find(id);
476     if (res == nodeID2idx.end()) {
477     reducedNodeIDs.push_back(id);
478     nodeID2idx[id] = newIdx;
479     *it = newIdx++;
480     } else {
481     *it = res->second;
482     }
483     }
484    
485     CoordArray coords;
486     const CoordArray& origCoords = originalMesh->getCoords();
487     for (int dim=0; dim < originalMesh->getNumDims(); dim++) {
488     float* c = new float[newIdx];
489     coords.push_back(c);
490     IndexMap::const_iterator mIt;
491     for (mIt = nodeID2idx.begin(); mIt != nodeID2idx.end(); mIt++) {
492     IndexMap::const_iterator idx;
493     idx = originalMesh->getIndexMap().find(mIt->first);
494     c[mIt->second] = origCoords[dim][idx->second];
495     }
496     }
497     if (reducedMesh)
498     delete reducedMesh;
499    
500     reducedMesh = new Mesh(coords, originalMesh->getNumDims(), newIdx);
501     reducedMesh->setName(string("Reduced") + name);
502     reducedMesh->setIndexMap(nodeID2idx);
503     reducedMesh->setNodeIDs(reducedNodeIDs);
504    
505     #ifdef _DEBUG
506     cout << reducedMesh->getName() << " has " << newIdx
507     << " nodes, " << reducedCount << " elements" << endl;
508     #endif
509     }
510     }
511    
512     #if HAVE_SILO
513     //
514     //
515     //
516     inline int toSiloElementType(int type)
517     {
518     switch (type) {
519     case ZONETYPE_BEAM: return DB_ZONETYPE_BEAM;
520     case ZONETYPE_HEX: return DB_ZONETYPE_HEX;
521     case ZONETYPE_POLYGON: return DB_ZONETYPE_POLYGON;
522     case ZONETYPE_QUAD: return DB_ZONETYPE_QUAD;
523     case ZONETYPE_TET: return DB_ZONETYPE_TET;
524     case ZONETYPE_TRIANGLE: return DB_ZONETYPE_TRIANGLE;
525     }
526     return 0;
527     }
528     #endif
529    
530     //
531     //
532     //
533     bool ElementData::writeToSilo(DBfile* dbfile, const string& siloPath)
534     {
535     #if HAVE_SILO
536     if (count == 0)
537     return true;
538    
539     const char* meshName;
540     int numCells;
541     string varName, siloMeshName;
542     int ret;
543    
544     if (siloPath != "") {
545     ret = DBSetDir(dbfile, siloPath.c_str());
546     if (ret != 0)
547     return false;
548     }
549    
550     // write out the full mesh in any case
551     if (siloPath == "/")
552     siloMeshName = siloPath + fullMesh->getName();
553     else
554     siloMeshName = siloPath + string("/") + fullMesh->getName();
555    
556     int arraylen = count * nodesPerElement;
557     int eltype = toSiloElementType(type);
558     varName = name + string("_zones");
559     ret = DBPutZonelist2(dbfile, varName.c_str(), count,
560     originalMesh->getNumDims(), &nodes[0], arraylen, 0, 0,
561     numGhostElements, &eltype, &nodesPerElement, &count, 1, NULL);
562 caltinay 2196 if (ret == 0) {
563     CoordArray& coordbase = const_cast<CoordArray&>(fullMesh->getCoords());
564 caltinay 2183 ret = DBPutUcdmesh(dbfile, siloMeshName.c_str(),
565 caltinay 2196 originalMesh->getNumDims(), NULL, &coordbase[0],
566 caltinay 2183 fullMesh->getNumNodes(), count, varName.c_str(),
567     /*"facelist"*/NULL, DB_FLOAT, NULL);
568 caltinay 2196 }
569 caltinay 2183
570     // Point mesh is useful for debugging
571 caltinay 2196 if (0) {
572     CoordArray& coordbase = const_cast<CoordArray&>(fullMesh->getCoords());
573     DBPutPointmesh(dbfile, "/pointmesh",
574     originalMesh->getNumDims(), &coordbase[0],
575     fullMesh->getNumNodes(), DB_FLOAT, NULL);
576     }
577 caltinay 2183
578     if (ret != 0)
579     return false;
580    
581     // decide whether to additionally write out the reduced mesh
582     if (reducedCount > 0) {
583     if (siloPath == "/")
584     siloMeshName = siloPath + reducedMesh->getName();
585     else
586     siloMeshName = siloPath + string("/") + reducedMesh->getName();
587     arraylen = reducedCount * reducedNodesPerElement;
588     eltype = toSiloElementType(reducedType);
589     varName = string("Reduced") + name + string("_zones");
590     ret = DBPutZonelist2(dbfile, varName.c_str(), reducedCount,
591     originalMesh->getNumDims(), &reducedNodes[0], arraylen, 0, 0,
592     numReducedGhostElements, &eltype, &reducedNodesPerElement,
593     &reducedCount, 1, NULL);
594 caltinay 2196 if (ret == 0) {
595     CoordArray& coordbase = const_cast<CoordArray&>(reducedMesh->getCoords());
596 caltinay 2183 ret = DBPutUcdmesh(dbfile, siloMeshName.c_str(),
597 caltinay 2196 originalMesh->getNumDims(), NULL, &coordbase[0],
598 caltinay 2183 reducedMesh->getNumNodes(), reducedCount, varName.c_str(),
599     NULL, DB_FLOAT, NULL);
600 caltinay 2196 }
601 caltinay 2183 if (ret != 0)
602     return false;
603     numCells = reducedCount;
604     } else {
605     numCells = count;
606     }
607     meshName = siloMeshName.c_str();
608    
609     // finally, write out the element-centered variables on the correct mesh
610     varName = name + string("_Color");
611     ret = DBPutUcdvar1(dbfile, varName.c_str(), meshName,
612     (float*)&color[0], numCells, NULL, 0, DB_INT, DB_ZONECENT, NULL);
613     if (ret == 0) {
614     varName = name + string("_Id");
615     ret = DBPutUcdvar1(dbfile, varName.c_str(), meshName,
616     (float*)&ID[0], numCells, NULL, 0, DB_INT, DB_ZONECENT, NULL);
617     }
618     if (ret == 0) {
619     varName = name + string("_Owner");
620     ret = DBPutUcdvar1(dbfile, varName.c_str(), meshName,
621     (float*)&owner[0], numCells, NULL, 0, DB_INT, DB_ZONECENT, NULL);
622     }
623     if (ret == 0) {
624     varName = name + string("_Tag");
625     ret = DBPutUcdvar1(dbfile, varName.c_str(), meshName,
626     (float*)&tag[0], numCells, NULL, 0, DB_INT, DB_ZONECENT, NULL);
627     }
628    
629     DBSetDir(dbfile, "/");
630     return (ret == 0);
631    
632     #else // !HAVE_SILO
633     return false;
634     #endif
635     }
636    
637     //
638     //
639     //
640     FinleyElementInfo ElementData::getFinleyTypeInfo(ElementTypeId typeId)
641     {
642     FinleyElementInfo ret;
643     ret.multiCellIndices = NULL;
644     ret.elementFactor = 1;
645    
646     switch (typeId) {
647     case Point1_Contact://untested
648     case Line2Face_Contact://untested
649     case Line3Face_Contact://untested
650     case Line2Face://untested
651     case Line3Face://untested
652     case Point1://untested
653     cerr << "WARNING: Finley type " <<typeId<< " is untested!" << endl;
654     ret.elementSize = 1;
655     ret.elementType = ZONETYPE_POLYGON;
656     break;
657    
658     case Tri3Face_Contact://untested
659     case Tri3Face://untested
660     cerr << "WARNING: Finley type " <<typeId<< " is untested!" << endl;
661     case Line2_Contact:
662     case Rec4Face_Contact:
663     case Rec4Face:
664     case Line2:
665     ret.elementSize = ret.reducedElementSize = 2;
666     ret.elementType = ret.reducedElementType = ZONETYPE_BEAM;
667     break;
668    
669     case Line3:
670     ret.multiCellIndices = line3indices;
671     ret.elementFactor = 2;
672     // fall through
673     case Line3_Contact:
674     case Tri6Face_Contact://untested
675     case Rec8Face_Contact:
676     case Tri6Face://untested
677     case Rec8Face:
678     //VTK_QUADRATIC_EDGE
679     ret.elementSize = ret.reducedElementSize = 2;
680     ret.elementType = ret.reducedElementType = ZONETYPE_BEAM;
681     break;
682    
683     case Tet4Face_Contact://untested
684     case Tet4Face://untested
685     cerr << "WARNING: Finley type " <<typeId<< " is untested!" << endl;
686     case Tri3_Contact:
687     case Tri3:
688     ret.elementSize = ret.reducedElementSize = 3;
689     ret.elementType = ret.reducedElementType = ZONETYPE_TRIANGLE;
690     break;
691    
692     case Rec4_Contact:
693     case Hex8Face_Contact:
694     case Hex8Face:
695     case Rec4:
696     ret.elementSize = ret.reducedElementSize = 4;
697     ret.elementType = ret.reducedElementType = ZONETYPE_QUAD;
698     break;
699    
700     case Rec9:
701     ret.multiCellIndices = rec9indices;
702     ret.elementFactor = 4;
703     // fall through
704     case Rec9_Contact:
705     ret.elementSize = ret.reducedElementSize = 4;
706     ret.elementType = ret.reducedElementType = ZONETYPE_QUAD;
707     break;
708    
709     case Tet4:
710     ret.elementSize = ret.reducedElementSize = 4;
711     ret.elementType = ret.reducedElementType = ZONETYPE_TET;
712     break;
713    
714     case Tri6:
715     ret.multiCellIndices = tri6indices;
716     ret.elementFactor = 4;
717     // fall through
718     case Tri6_Contact:
719     case Tet10Face_Contact://untested
720     case Tet10Face://untested
721     //VTK_QUADRATIC_TRIANGLE
722     ret.elementSize = ret.reducedElementSize = 3;
723     ret.elementType = ret.reducedElementType = ZONETYPE_TRIANGLE;
724     break;
725    
726     case Rec8:
727     ret.multiCellIndices = rec8indices;
728     ret.elementFactor = 6;
729     // fall through
730     case Hex20Face:
731     case Rec8_Contact:
732     case Hex20Face_Contact:
733     //VTK_QUADRATIC_QUAD
734     ret.elementSize = 3;
735     ret.elementType = ZONETYPE_TRIANGLE;
736     ret.reducedElementSize = 4;
737     ret.reducedElementType = ZONETYPE_QUAD;
738     break;
739    
740     case Tet10:
741     //VTK_QUADRATIC_TETRA
742     ret.multiCellIndices = tet10indices;
743     ret.elementFactor = 8;
744     ret.elementSize = ret.reducedElementSize = 4;
745     ret.elementType = ret.reducedElementType = ZONETYPE_TET;
746     break;
747    
748     case Hex20:
749     //VTK_QUADRATIC_HEXAHEDRON
750     ret.multiCellIndices = hex20indices;
751     ret.elementFactor = 36;
752     ret.elementSize = 3;
753     ret.elementType = ZONETYPE_TRIANGLE;
754     ret.reducedElementSize = 8;
755     ret.reducedElementType = ZONETYPE_HEX;
756     break;
757    
758     case Hex27:
759     ret.multiCellIndices = hex27indices;
760     ret.elementFactor = 8;
761     // fall through
762     case Hex8:
763     ret.elementSize = ret.reducedElementSize = 8;
764     ret.elementType = ret.reducedElementType = ZONETYPE_HEX;
765     break;
766    
767     default:
768     cerr << "WARNING: Unknown Finley Type " << typeId << endl;
769     break;
770     }
771     return ret;
772     }
773    
774 caltinay 2187 } // namespace EscriptReader
775    

  ViewVC Help
Powered by ViewVC 1.1.26