Skip to content
Snippets Groups Projects
Commit 1c482538 authored by Oliver Sander's avatar Oliver Sander
Browse files

Fused 2d- and 3d-writer into one method. And the completely ridiculous

#ifdef _3 stuff disappeared.

Added code the writes the leaf level (and not a fixed level) into an
AmiraMesh file.  Not working yet, though.

[[Imported from SVN: r2082]]
parent 9a543d60
Branches
Tags
No related merge requests found
......@@ -6,307 +6,367 @@
#include <dune/istl/bvector.hh>
#if _3
/** \todo Make sure that the grid is three-dimensional */
// ///////////////////////////////////////////////////////////
// Write the leaf level of the grid as an AmiraMesh file
// ///////////////////////////////////////////////////////////
template<class GridType>
void Dune::AmiraMeshWriter<GridType>::writeGrid(const GridType& grid,
const std::string& filename)
{
// Temporary: we write this level
int level = grid.maxlevel();
// Find out whether the grid contains only tetrahedra. If yes, then
// it is written in TetraGrid format. If not, it is written in
// hexagrid format.
bool containsOnlyTetrahedra = true;
const int maxlevel = grid.maxlevel();
typedef typename GridType::template codim<0>::LevelIterator ElementIterator;
#ifndef UGGRID_WITH_INDEX_SETS
writeGrid(grid, filename, maxlevel);
#else
// Find out whether the grid contains only triangles. If yes, then
// it is written as a HxTriangularGrid. If not, it is written
// but it cannot (currently) be read by Amira
bool containsOnlySimplices = true;
ElementIterator eIt = grid.template lbegin<0>(level);
ElementIterator eEndIt = grid.template lend<0>(level);
LeafIterator element = grid.leafbegin(maxlevel);
LeafIterator end = grid.leafend(maxlevel);
for (; eIt!=eEndIt; ++eIt) {
if (eIt->geometry().type() != tetrahedron) {
containsOnlyTetrahedra = false;
for (; element!=end; ++element) {
if (element->geometry().type() != triangle) {
containsOnlySimplices = false;
break;
}
}
int maxVerticesPerElement = (containsOnlyTetrahedra) ? 4 : 8;
int maxVerticesPerElement = (containsOnlySimplices) ? 3 : 4;
const int DIM = 3;
// /////////////////////////////////////////////////////
// Extract the leaf vertices
// /////////////////////////////////////////////////////
int noOfNodes = grid.size(level, 3);
int noOfElem = grid.size(level, 0);
LeafIterator lIt = grid.leafbegin(maxlevel);
LeafIterator lEndIt = grid.leafend(maxlevel);
int i;
int noOfElem = 0;
Array<int> leafVertices(grid.hierarchicIndexSet().size(0,dim));
Array<FieldVector<double,dim> > leafVertexCoords(grid.hierarchicIndexSet().size(0,dim));
// create amiramesh object
AmiraMesh am;
leafVertices.set(-1);
// write grid vertex coordinates
AmiraMesh::Location* geo_nodes = new AmiraMesh::Location("Nodes", noOfNodes);
am.insert(geo_nodes);
for (; lIt!=lEndIt; ++lIt) {
AmiraMesh::Data* geo_node_data = new AmiraMesh::Data("Coordinates", geo_nodes,
McPrimType::mc_float, DIM);
am.insert(geo_node_data);
typedef typename GridType::template codim<3>::LevelIterator VertexIterator;
VertexIterator vertex = grid.template lbegin<DIM>(level);
VertexIterator endvertex = grid.template lend<DIM>(level);
for (; vertex!=endvertex; ++vertex) {
#ifndef UGGRID_WITH_INDEX_SETS
int index = vertex->index();
#else
int index = grid.levelIndexSet().index(*vertex);
#endif
const FieldVector<double, 3>& coords = vertex->geometry()[0];
((float*)geo_node_data->dataPtr())[3*index+0] = coords[0];
((float*)geo_node_data->dataPtr())[3*index+1] = coords[1];
((float*)geo_node_data->dataPtr())[3*index+2] = coords[2];
for (int i=0; i<lIt->template count<dim>(); i++) {
leafVertices[grid.hierarchicIndexSet().template subIndex<dim>(*lIt,i)] = 1;
leafVertexCoords[grid.hierarchicIndexSet().template subIndex<dim>(*lIt,i)] = lIt->geometry()[i];
}
noOfElem++;
}
/* write element section to file */
AmiraMesh::Location* element_loc = NULL;
int noOfNodes = 0;
for (int i=0; i<leafVertices.size(); i++)
if (leafVertices[i]==1)
leafVertices[i] = noOfNodes++;
if (containsOnlyTetrahedra)
element_loc = new AmiraMesh::Location("Tetrahedra", noOfElem);
else
element_loc = new AmiraMesh::Location("Hexahedra", noOfElem);
// create amiramesh object
AmiraMesh am;
am.insert(element_loc);
// Set the appropriate content type
am.parameters.set("ContentType", "HxTriangularGrid");
AmiraMesh::Data* element_data = new AmiraMesh::Data("Nodes", element_loc,
McPrimType::mc_int32, maxVerticesPerElement);
am.insert(element_data);
// //////////////////////////////////
// Write grid vertex coordinates
// //////////////////////////////////
AmiraMesh::Location* geoNodes = new AmiraMesh::Location("Nodes", noOfNodes);
am.insert(geoNodes);
int *dPtr = (int*)element_data->dataPtr();
AmiraMesh::Data* geo_node_data = new AmiraMesh::Data("Coordinates", geoNodes,
McPrimType::mc_float, dim);
am.insert(geo_node_data);
eIt = grid.template lbegin<0>(level);
int vertexIdx = 0;
for (int i=0; i<leafVertices.size(); i++)
if (leafVertices[i] != -1) {
if (containsOnlyTetrahedra) {
for (int j=0; j<dim; j++)
((float*)geo_node_data->dataPtr())[dim*vertexIdx+j] = leafVertexCoords[i][j];
for (i=0; eIt!=eEndIt; ++eIt, i++) {
vertexIdx++;
for (int j=0; j<4; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[i*4+j] = eIt->template subIndex<3>(j)+1;
#else
dPtr[i*4+j] = grid.levelIndexSet().template subIndex<3>(*eIt,j)+1;
#endif
}
} else {
for (i=0; eIt!=eEndIt; ++eIt, i++) {
switch (eIt->geometry().type()) {
// ///////////////////////////////////////////
// Write element section
// ///////////////////////////////////////////
AmiraMesh::Location* elementLoc = (containsOnlySimplices)
? new AmiraMesh::Location("Triangles", noOfElem)
: new AmiraMesh::Location("Quadrangles", noOfElem);
case hexahedron : {
am.insert(elementLoc);
const int hexaReordering[8] = {0, 1, 3, 2, 4, 5, 7, 6};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(hexaReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, hexaReordering[j])+1;
#endif
break;
}
case prism : {
const int prismReordering[8] = {0, 1, 1, 2, 3, 4, 4, 5};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(prismReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, prismReordering[j])+1;
#endif
break;
}
case pyramid : {
const int pyramidReordering[8] = {0, 1, 2, 3, 4, 4, 4, 4};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(pyramidReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, pyramidReordering[j])+1;
#endif
AmiraMesh::Data* elementData = new AmiraMesh::Data("Nodes", elementLoc,
McPrimType::mc_int32, maxVerticesPerElement);
am.insert(elementData);
break;
}
int *dPtr = (int*)elementData->dataPtr();
case tetrahedron : {
LeafIterator lIt2 = grid.leafbegin(maxlevel);
//LeafIterator lEndIt = grid.leafend(maxlevel);
const int tetraReordering[8] = {0, 1, 2, 2, 3, 3, 3, 3};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(tetraReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, tetraReordering[j])+1;
#endif
for (int i=0; lIt2!=lEndIt; ++lIt2, i++) {
break;
}
switch (lIt2->geometry().type()) {
default :
DUNE_THROW(NotImplemented, "Unknown element type encountered");
default :
}
for (int j=0; j<lIt2->template count<dim>(); j++)
dPtr[i*maxVerticesPerElement+j] = leafVertices[grid.hierarchicIndexSet().template subIndex<dim>(*lIt2,j)]+1;
// If the element has less than four vertices use the last value
// to fill up the remaining slot
for (int j=lIt2->template count<dim>(); j<maxVerticesPerElement; j++)
dPtr[i*maxVerticesPerElement+j] = dPtr[i*maxVerticesPerElement+lIt2->template count<dim>()-1];
}
}
// write material section to grid file
AmiraMesh::Data* element_materials = new AmiraMesh::Data("Materials", element_loc, McPrimType::mc_uint8, 1);
am.insert(element_materials);
// for(i=0; i<noOfElem; i++)
// ((unsigned char*)element_materials->dataPtr())[i] = SUBDOMAIN(elemList[i]);
// write material section to AmiraMesh object
AmiraMesh::Data* elementMaterials = new AmiraMesh::Data("Materials", elementLoc, McPrimType::mc_uint8, 1);
am.insert(elementMaterials);
for(i=0; i<noOfElem; i++)
((unsigned char*)element_materials->dataPtr())[i] = 0;
for(int i=0; i<noOfElem; i++)
((unsigned char*)elementMaterials->dataPtr())[i] = 0;
// Actually write the file
if(!am.write(filename.c_str(), 1))
DUNE_THROW(IOError, "Writing geometry file failed!");
DUNE_THROW(IOError, "Writing geometry file failed");
std::cout << "Grid written successfully to: " << filename << std::endl;
#endif
}
#endif // #if _3
/*************************************************************************/
/* Write one level of a 3d-grid into an AmiraMesh file */
/*************************************************************************/
#ifdef _2
template<class GridType>
void Dune::AmiraMeshWriter<GridType>::writeGrid(const GridType& grid,
const std::string& filename)
const std::string& filename, int level)
{
if ((dim!=2 && dim!=3) || int(dim) != int(GridType::dimensionworld))
DUNE_THROW(IOError, "You can only write grids as AmiraMesh if dim==dimworld==2"
<< " or dim==dimworld==3.");
// Temporary: we write this level
int level = grid.maxlevel();
// Find out whether the grid contains only triangles. If yes, then
// it is written as a HxTriangularGrid. If not, it cannot be
// written (so far).
bool containsOnlyTetrahedra = true;
// Find out whether the grid contains only tetrahedra. If yes, then
// it is written in TetraGrid format. If not, it is written in
// hexagrid format.
bool containsOnlySimplices = true;
typename GridType::template codim<0>::LevelIterator element = grid.template lbegin<0>(level);
typename GridType::template codim<0>::LevelIterator end = grid.template lend<0>(level);
ElementIterator eIt = grid.template lbegin<0>(level);
ElementIterator eEndIt = grid.template lend<0>(level);
for (; element!=end; ++element) {
if (element->geometry().type() != triangle) {
containsOnlyTetrahedra = false;
for (; eIt!=eEndIt; ++eIt) {
if (eIt->geometry().type() != tetrahedron &&
eIt->geometry().type() != triangle &&
eIt->geometry().type() != simplex) {
containsOnlySimplices = false;
break;
}
}
int maxVerticesPerElement = (containsOnlyTetrahedra) ? 3 : 4;
int maxVerticesPerElement = (dim==3)
? ((containsOnlySimplices) ? 4 : 8)
: ((containsOnlySimplices) ? 3 : 4);
const int DIM = 2;
int noOfNodes = grid.size(level, DIM);
int noOfNodes = grid.size(level, dim);
int noOfElem = grid.size(level, 0);
int i;
// Construct the name for the file
std::string geoFilename(filename);
geoFilename += ".am";
// create amiramesh object
AmiraMesh am_geometry;
AmiraMesh am;
// Set the appropriate content type
am_geometry.parameters.set("ContentType", "HxTriangularGrid");
if (dim==2)
am.parameters.set("ContentType", "HxTriangularGrid");
// write grid vertex coordinates
AmiraMesh::Location* geo_nodes = new AmiraMesh::Location("Nodes", noOfNodes);
am_geometry.insert(geo_nodes);
am.insert(geo_nodes);
AmiraMesh::Data* geo_node_data = new AmiraMesh::Data("Coordinates", geo_nodes,
McPrimType::mc_float, DIM);
am_geometry.insert(geo_node_data);
McPrimType::mc_float, dim);
am.insert(geo_node_data);
typedef typename GridType::template codim<dim>::LevelIterator VertexIterator;
VertexIterator vertex = grid.template lbegin<dim>(level);
VertexIterator endvertex = grid.template lend<dim>(level);
typename GridType::template codim<2>::LevelIterator vertex = grid.template lbegin<DIM>(level);
typename GridType::template codim<2>::LevelIterator endvertex = grid.template lend<DIM>(level);
for (; vertex!=endvertex; ++vertex) {
for (; vertex!=endvertex; ++vertex)
{
#ifndef UGGRID_WITH_INDEX_SETS
int index = vertex->index();
#else
int index = grid.levelIndexSet().index(*vertex);
#endif
const FieldVector<double, DIM>& coords = vertex->geometry()[0];
const FieldVector<double, dim>& coords = vertex->geometry()[0];
((float*)geo_node_data->dataPtr())[2*index+0] = coords[0];
((float*)geo_node_data->dataPtr())[2*index+1] = coords[1];
// Copy coordinates
for (int i=0; i<dim; i++)
((float*)geo_node_data->dataPtr())[dim*index+i] = coords[i];
}
/* write element section to geo - file */
/* write element section to file */
AmiraMesh::Location* element_loc = NULL;
if (containsOnlyTetrahedra)
element_loc = new AmiraMesh::Location("Triangles", noOfElem);
else
element_loc = new AmiraMesh::Location("Quadrangles", noOfElem);
if (dim==3) {
am_geometry.insert(element_loc);
if (containsOnlySimplices)
element_loc = new AmiraMesh::Location("Tetrahedra", noOfElem);
else
element_loc = new AmiraMesh::Location("Hexahedra", noOfElem);
} else {
if (containsOnlySimplices)
element_loc = new AmiraMesh::Location("Triangles", noOfElem);
else
element_loc = new AmiraMesh::Location("Quadrangles", noOfElem);
}
am.insert(element_loc);
AmiraMesh::Data* element_data = new AmiraMesh::Data("Nodes", element_loc,
McPrimType::mc_int32, maxVerticesPerElement);
am_geometry.insert(element_data);
am.insert(element_data);
int *dPtr = (int*)element_data->dataPtr();
typename GridType::template codim<0>::LevelIterator element2 = grid.template lbegin<0>(level);
typename GridType::template codim<0>::LevelIterator endelement = grid.template lend<0>(level);
eIt = grid.template lbegin<0>(level);
for (i=0; element2!=endelement; ++element2, i++) {
switch (element2->geometry().type()) {
if (dim==3) {
default :
// //////////////////////////////////////////////////
// Write elements of a 3D-grid
// //////////////////////////////////////////////////
if (containsOnlySimplices) {
for (i=0; eIt!=eEndIt; ++eIt, i++) {
for (int j=0; j<4; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[i*4+j] = eIt->template subIndex<3>(j)+1;
#else
dPtr[i*4+j] = grid.levelIndexSet().template subIndex<3>(*eIt,j)+1;
#endif
}
} else {
for (i=0; eIt!=eEndIt; ++eIt, i++) {
switch (eIt->geometry().type()) {
case hexahedron : {
const int hexaReordering[8] = {0, 1, 3, 2, 4, 5, 7, 6};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(hexaReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, hexaReordering[j])+1;
#endif
break;
}
case prism : {
const int prismReordering[8] = {0, 1, 1, 2, 3, 4, 4, 5};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(prismReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, prismReordering[j])+1;
#endif
break;
}
case pyramid : {
const int pyramidReordering[8] = {0, 1, 2, 3, 4, 4, 4, 4};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(pyramidReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, pyramidReordering[j])+1;
#endif
break;
}
case tetrahedron : {
const int tetraReordering[8] = {0, 1, 2, 2, 3, 3, 3, 3};
for (int j=0; j<8; j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[8*i + j] = eIt->template subIndex<3>(tetraReordering[j])+1;
#else
dPtr[8*i + j] = grid.levelIndexSet().template subIndex<3>(*eIt, tetraReordering[j])+1;
#endif
break;
}
default :
DUNE_THROW(NotImplemented, "Unknown element type encountered");
}
for (int j=0; j<element2->geometry().corners(); j++)
}
}
} else {
typename GridType::template codim<0>::LevelIterator element2 = grid.template lbegin<0>(level);
typename GridType::template codim<0>::LevelIterator endelement = grid.template lend<0>(level);
for (i=0; element2!=endelement; ++element2, i++) {
switch (element2->geometry().type()) {
default :
for (int j=0; j<element2->geometry().corners(); j++)
#ifndef UGGRID_WITH_INDEX_SETS
dPtr[i*maxVerticesPerElement+j] = element2->template subIndex<DIM>(j)+1;
dPtr[i*maxVerticesPerElement+j] = element2->template subIndex<dim>(j)+1;
#else
dPtr[i*maxVerticesPerElement+j] = grid.levelIndexSet().template subIndex<DIM>(*element2, j)+1;
dPtr[i*maxVerticesPerElement+j] = grid.levelIndexSet().template subIndex<dim>(*element2, j)+1;
#endif
// If the element has less than 8 vertices use the last value
// to fill up the remaining slots
for (int j=element2->geometry().corners(); j<maxVerticesPerElement; j++)
dPtr[i*maxVerticesPerElement+j] = dPtr[i*maxVerticesPerElement+element2->geometry().corners()-1];
// If the element has less than 8 vertices use the last value
// to fill up the remaining slots
for (int j=element2->geometry().corners(); j<maxVerticesPerElement; j++)
dPtr[i*maxVerticesPerElement+j] = dPtr[i*maxVerticesPerElement+element2->geometry().corners()-1];
}
}
}
// write material section to geo-file
// write material section to grid file
AmiraMesh::Data* element_materials = new AmiraMesh::Data("Materials", element_loc, McPrimType::mc_uint8, 1);
am_geometry.insert(element_materials);
am.insert(element_materials);
for(i=0; i<noOfElem; i++)
((unsigned char*)element_materials->dataPtr())[i] = 0;
// Actually write the file
if(!am.write(filename.c_str(), 1))
DUNE_THROW(IOError, "Writing geometry file failed!");
std::cout << "Grid written successfully to: " << filename << std::endl;
}
if(!am_geometry.write(geoFilename.c_str(), 1))
DUNE_THROW(IOError, "Writing geometry file failed");
std::cout << "Grid written successfully to: " << geoFilename << std::endl;
}
#endif // #ifdef _2
template<class GridType>
template<class DiscFuncType>
......@@ -371,8 +431,8 @@ void Dune::AmiraMeshWriter<GridType>::writeBlockVector(const GridType& grid,
// hexagrid format.
bool containsOnlyTetrahedra = true;
typename GridType::template codim<0>::LevelIterator element = grid.template lbegin<0>(level);
typename GridType::template codim<0>::LevelIterator end = grid.template lend<0>(level);
ElementIterator element = grid.template lbegin<0>(level);
ElementIterator end = grid.template lend<0>(level);
for (; element!=end; ++element) {
if (element->geometry().type() != tetrahedron) {
......
......@@ -15,9 +15,15 @@ namespace Dune {
template<class GridType>
class AmiraMeshWriter {
enum {dim = GridType::dimension};
typedef typename GridType::template codim<0>::LevelIterator ElementIterator;
typedef typename GridType::template codim<dim>::LevelIterator VertexIterator;
typedef typename GridType::LeafIterator LeafIterator;
public:
/** \brief Writes a grid in AmiraMesh format
/** \brief Write the leaf level of a grid in AmiraMesh format
*
* @param grid The grid objects that is to be written
* @param filename The filename
......@@ -25,6 +31,16 @@ namespace Dune {
static void writeGrid(const GridType& grid,
const std::string& filename);
/** \brief Write one level of a grid in AmiraMesh format
*
* @param grid The grid objects that is to be written
* @param filename The filename
* @param level The level to be written
*/
static void writeGrid(const GridType& grid,
const std::string& filename,
int level);
/** \brief Writes a discrete function in AmiraMesh format
*
* @param f Function that should be written
......@@ -49,6 +65,6 @@ namespace Dune {
#include "amiramesh/amirameshwriter.cc"
// the amiramesh writer for SGrid
//#include "amiramesh/amsgridwriter.cc"
#include "amiramesh/amsgridwriter.cc"
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment