OilfieldToolsNet/LBPM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Fixing issues with enums and the summary files

Browse Source
This commit is contained in:
Mark Berrill 2021-03-18 09:41:39 -04:00
parent dfa97a0132
commit 59b7b9a0fe
10 changed files with 480 additions and 351 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

170
IO/Mesh.cpp
View File

@ -1,4 +1,5 @@
#include "Mesh.h"
#include "IO/IOHelpers.h"
#include "common/Utilities.h"
#include <limits>
@ -28,16 +29,23 @@ Mesh::~Mesh() {}
/****************************************************
* MeshDataStruct *
****************************************************/
bool MeshDataStruct::check() const
#define checkResult( pass, msg ) \
do { \
if ( !( pass ) ) { \
if ( abort ) \
ERROR( msg ); \
return false; \
} \
} while ( 0 )
bool MeshDataStruct::check( bool abort ) const
{
bool pass = mesh != nullptr;
for ( const auto &var : vars ) {
pass = pass && static_cast<int>( var->type ) >= 1 && static_cast<int>( var->type ) <= 3;
pass = pass && !var->data.empty();
}
if ( !pass ) {
std::cerr << "Invalid variable detected\n";
return false;
checkResult( var->type == VariableType::NodeVariable ||
var->type == VariableType::EdgeVariable ||
var->type == VariableType::SurfaceVariable ||
var->type == VariableType::VolumeVariable,
"Invalid data type" );
checkResult( !var->data.empty(), "Variable data is empty" );
}
const std::string &meshClass = mesh->className();
if ( meshClass == "PointList" ) {
@ -45,7 +53,9 @@ bool MeshDataStruct::check() const
ASSERT( mesh2 );
for ( const auto &var : vars ) {
if ( var->type == IO::VariableType::NodeVariable ) {
pass = pass && var->data.size() == ArraySize( mesh2->points.size(), var->dim );
size_t N_points = mesh2->points.size();
checkResult( var->data.size( 0 ) == N_points, "sizeof NodeVariable" );
checkResult( var->data.size( 1 ) == var->dim, "sizeof NodeVariable" );
} else if ( var->type == IO::VariableType::EdgeVariable ) {
ERROR( "Invalid type for PointList" );
} else if ( var->type == IO::VariableType::SurfaceVariable ) {
@ -61,15 +71,16 @@ bool MeshDataStruct::check() const
ASSERT( mesh2 );
for ( const auto &var : vars ) {
if ( var->type == IO::VariableType::NodeVariable ) {
pass = pass &&
var->data.size() == ArraySize( mesh2->vertices->points.size(), var->dim );
size_t N_points = mesh2->vertices->points.size();
checkResult( var->data.size( 0 ) == N_points, "sizeof NodeVariable" );
checkResult( var->data.size( 1 ) == var->dim, "sizeof NodeVariable" );
} else if ( var->type == IO::VariableType::EdgeVariable ) {
ERROR( "Not finished" );
} else if ( var->type == IO::VariableType::SurfaceVariable ) {
ERROR( "Not finished" );
} else if ( var->type == IO::VariableType::VolumeVariable ) {
pass = pass && var->data.size( 0 ) == mesh2->A.size() &&
var->data.size( 1 ) == var->dim;
checkResult( var->data.size( 0 ) == mesh2->A.size(), "sizeof VolumeVariable" );
checkResult( var->data.size( 1 ) == var->dim, "sizeof VolumeVariable" );
} else {
ERROR( "Invalid variable type" );
}
@ -90,14 +101,16 @@ bool MeshDataStruct::check() const
} else {
ERROR( "Invalid variable type" );
}
if ( var->data.size() == ArraySize( varSize[0] * varSize[1] * varSize[2], varSize[3] ) )
if ( var->data.size( 0 ) == varSize[0] * varSize[1] * varSize[2] &&
var->data.size( 1 ) == varSize[3] )
var->data.resize( varSize );
pass = pass && var->data.size() == varSize;
for ( int d = 0; d < 4; d++ )
checkResult( var->data.size( d ) == varSize[d], "DomainMesh Variable" );
}
} else {
ERROR( "Unknown mesh class: " + mesh->className() );
}
return pass;
return true;
}
@ -479,4 +492,129 @@ std::shared_ptr<const TriList> getTriList( std::shared_ptr<const Mesh> mesh )
}
/****************************************************
* Convert enum values *
****************************************************/
std::string getString( VariableType type )
{
if ( type == VariableType::NodeVariable )
return "node";
else if ( type == VariableType::EdgeVariable )
return "edge";
else if ( type == VariableType::SurfaceVariable )
return "face";
else if ( type == VariableType::VolumeVariable )
return "cell";
else if ( type == VariableType::NullVariable )
return "null";
else
ERROR( "Invalid type" );
return "";
}
VariableType getVariableType( const std::string &type_in )
{
auto type = deblank( type_in );
if ( type == "node" )
return VariableType::NodeVariable;
else if ( type == "edge" || type == "1" )
return VariableType::EdgeVariable;
else if ( type == "face" )
return VariableType::SurfaceVariable;
else if ( type == "cell" || type == "3" )
return VariableType::VolumeVariable;
else if ( type == "null" )
return VariableType::NullVariable;
else
ERROR( "Invalid type: " + type );
return VariableType::NullVariable;
}
std::string getString( DataType type )
{
if ( type == DataType::Double )
return "double";
else if ( type == DataType::Float )
return "float";
else if ( type == DataType::Int )
return "int";
else if ( type == DataType::Null )
return "null";
else
ERROR( "Invalid type" );
return "";
}
DataType getDataType( const std::string &type_in )
{
auto type = deblank( type_in );
if ( type == "double" )
return DataType::Double;
else if ( type == "float" )
return DataType::Float;
else if ( type == "int" )
return DataType::Int;
else if ( type == "null" )
return DataType::Null;
else
ERROR( "Invalid type: " + type );
return DataType::Null;
}
std::string getString( MeshType type )
{
if ( type == MeshType::PointMesh )
return "PointMesh";
else if ( type == MeshType::SurfaceMesh )
return "SurfaceMesh";
else if ( type == MeshType::VolumeMesh )
return "VolumeMesh";
else if ( type == MeshType::Unknown )
return "unknown";
else
ERROR( "Invalid type" );
return "";
}
MeshType getMeshType( const std::string &type_in )
{
auto type = deblank( type_in );
if ( type == "PointMesh" || type == "1" )
return MeshType::PointMesh;
else if ( type == "SurfaceMesh" || type == "2" )
return MeshType::SurfaceMesh;
else if ( type == "VolumeMesh" || type == "3" )
return MeshType::VolumeMesh;
else if ( type == "unknown" || type == "-1" )
return MeshType::Unknown;
else
ERROR( "Invalid type: " + type );
return MeshType::Unknown;
}
std::string getString( FileFormat type )
{
if ( type == FileFormat::OLD )
return "old";
else if ( type == FileFormat::NEW )
return "new";
else if ( type == FileFormat::NEW_SINGLE )
return "new(single)";
else if ( type == FileFormat::SILO )
return "silo";
else
ERROR( "Invalid type" );
return "";
}
FileFormat getFileFormat( const std::string &type_in )
{
auto type = deblank( type_in );
if ( type == "old" || type == "1" )
return FileFormat::OLD;
else if ( type == "new" || type == "2" )
return FileFormat::NEW;
else if ( type == "new(single)" || type == "3" )
return FileFormat::NEW_SINGLE;
else if ( type == "silo" || type == "4" )
return FileFormat::SILO;
else
ERROR( "Invalid type: " + type );
return FileFormat::SILO;
}
} // namespace IO

31
IO/Mesh.h
View File

@ -14,15 +14,28 @@
namespace IO {
//! Possible variable types
enum class VariableType : unsigned char {
NodeVariable = 1,
EdgeVariable = 2,
SurfaceVariable = 3,
VolumeVariable = 4,
NullVariable = 0
//! Enums to define types
enum class VariableType {
NodeVariable,
EdgeVariable,
SurfaceVariable,
VolumeVariable,
NullVariable
};
enum class DataType : unsigned char { Double = 1, Float = 2, Int = 2, Null = 0 };
enum class DataType { Double, Float, Int, Null };
enum class MeshType { PointMesh, SurfaceMesh, VolumeMesh, Unknown };
enum class FileFormat { OLD, NEW, NEW_SINGLE, SILO };
//! Convert enums to/from strings (more future-proof than static_cast<int>)
std::string getString( VariableType );
std::string getString( DataType );
std::string getString( MeshType );
std::string getString( FileFormat );
VariableType getVariableType( const std::string & );
DataType getDataType( const std::string & );
MeshType getMeshType( const std::string & );
FileFormat getFileFormat( const std::string & );
/*! \class Mesh
@ -216,7 +229,7 @@ struct MeshDataStruct {
//! Empty constructor
MeshDataStruct() : precision( DataType::Double ) {}
//! Check the data
bool check() const;
bool check( bool abort = true ) const;
};

76
IO/MeshDatabase.cpp
View File

@ -13,38 +13,31 @@
#include <ProfilerApp.h>
// MeshType
template<>
size_t packsize<IO::MeshType>( const IO::MeshType &rhs )
{
return sizeof( IO::MeshType );
}
template<>
void pack<IO::MeshType>( const IO::MeshType &rhs, char *buffer )
{
memcpy( buffer, &rhs, sizeof( IO::MeshType ) );
}
template<>
void unpack<IO::MeshType>( IO::MeshType &data, const char *buffer )
{
memcpy( &data, buffer, sizeof( IO::MeshType ) );
}
// Variable::VariableType
template<>
size_t packsize<IO::VariableType>( const IO::VariableType &rhs )
{
return sizeof( IO::VariableType );
}
template<>
void pack<IO::VariableType>( const IO::VariableType &rhs, char *buffer )
{
memcpy( buffer, &rhs, sizeof( IO::VariableType ) );
}
template<>
void unpack<IO::VariableType>( IO::VariableType &data, const char *buffer )
{
memcpy( &data, buffer, sizeof( IO::VariableType ) );
}
// Default pack/unpack
// clang-format off
#define INSTANTIATE_PACK( TYPE ) \
template<> \
size_t packsize<TYPE>( const TYPE &rhs ) \
{ \
return sizeof( TYPE ); \
} \
template<> \
void pack<TYPE>( const TYPE &rhs, char *buffer ) \
{ \
memcpy( buffer, &rhs, sizeof( IO::MeshType ) ); \
} \
template<> \
void unpack<TYPE>( TYPE &data, const char *buffer ) \
{ \
memcpy( &data, buffer, sizeof( IO::MeshType ) ); \
}
INSTANTIATE_PACK( IO::VariableType )
INSTANTIATE_PACK( IO::DataType )
INSTANTIATE_PACK( IO::MeshType )
INSTANTIATE_PACK( IO::FileFormat )
// clang-format on
// DatabaseEntry
template<>
size_t packsize<IO::DatabaseEntry>( const IO::DatabaseEntry &rhs )
@ -327,8 +320,7 @@ std::vector<MeshDatabase> gatherAll(
// Return the results
std::vector<MeshDatabase> data2( data.size() );
size_t i = 0;
for ( std::map<std::string, MeshDatabase>::iterator it = data.begin(); it != data.end();
++it, ++i )
for ( auto it = data.begin(); it != data.end(); ++it, ++i )
data2[i] = it->second;
PROFILE_STOP( "gatherAll-unpack", 2 );
PROFILE_STOP( "gatherAll" );
@ -343,19 +335,19 @@ void write( const std::vector<MeshDatabase> &meshes, const std::string &filename
FILE *fid = fopen( filename.c_str(), "wb" );
for ( size_t i = 0; i < meshes.size(); i++ ) {
fprintf( fid, "%s\n", meshes[i].name.c_str() );
fprintf( fid, " type: %i\n", static_cast<int>( meshes[i].type ) );
fprintf( fid, " type: %s\n", getString( meshes[i].type ).data() );
fprintf( fid, " meshClass: %s\n", meshes[i].meshClass.c_str() );
fprintf( fid, " format: %i\n", static_cast<int>( meshes[i].format ) );
fprintf( fid, " format: %s\n", getString( meshes[i].format ).data() );
for ( size_t j = 0; j < meshes[i].domains.size(); j++ )
fprintf( fid, " domain: %s\n", meshes[i].domains[j].write().c_str() );
fprintf( fid, " variables: " );
for ( size_t j = 0; j < meshes[i].variables.size(); j++ ) {
const VariableDatabase &var = meshes[i].variables[j];
fprintf( fid, "%s|%i|%i; ", var.name.c_str(), static_cast<int>( var.type ), var.dim );
fprintf( fid, "%s|%s|%i; ", var.name.data(), getString( var.type ).data(), var.dim );
}
fprintf( fid, "\n" );
std::map<std::pair<std::string, std::string>, DatabaseEntry>::const_iterator it;
for ( it = meshes[i].variable_data.begin(); it != meshes[i].variable_data.end(); ++it ) {
for ( auto it = meshes[i].variable_data.begin(); it != meshes[i].variable_data.end();
++it ) {
const char *domain = it->first.first.c_str();
const char *variable = it->first.second.c_str();
fprintf(
@ -386,9 +378,9 @@ std::vector<MeshDatabase> read( const std::string &filename )
name.resize( name.size() - 1 );
meshes.back().name = name;
} else if ( strncmp( line, " format:", 10 ) == 0 ) {
meshes.back().format = static_cast<unsigned char>( atoi( &line[10] ) );
meshes.back().format = getFileFormat( &line[10] );
} else if ( strncmp( line, " type:", 8 ) == 0 ) {
meshes.back().type = static_cast<MeshType>( atoi( &line[8] ) );
meshes.back().type = getMeshType( &line[8] );
} else if ( strncmp( line, " meshClass:", 13 ) == 0 ) {
meshes.back().meshClass = deblank( std::string( &line[13] ) );
} else if ( strncmp( line, " domain:", 10 ) == 0 ) {
@ -402,7 +394,7 @@ std::vector<MeshDatabase> read( const std::string &filename )
std::vector<std::string> tmp = splitList( variables[i].c_str(), '|' );
ASSERT( tmp.size() == 3 );
mesh.variables[i].name = tmp[0];
mesh.variables[i].type = static_cast<VariableType>( atoi( tmp[1].c_str() ) );
mesh.variables[i].type = getVariableType( tmp[1] );
mesh.variables[i].dim = atoi( tmp[2].c_str() );
}
} else if ( strncmp( line, " variable(", 12 ) == 0 ) {

9
IO/MeshDatabase.h
View File

@ -13,13 +13,6 @@
namespace IO {
class Mesh;
//! Enum to identify mesh type
// enum class MeshType : char { PointMesh=1, SurfaceMesh=2, VolumeMesh=3, Unknown=-1 };
enum class MeshType { PointMesh = 1, SurfaceMesh = 2, VolumeMesh = 3, Unknown = -1 };
//! Helper struct for containing offsets for the mesh info
struct DatabaseEntry {
@ -56,7 +49,7 @@ struct MeshDatabase {
std::string name; //!< Name of the mesh
MeshType type; //!< Mesh type
std::string meshClass; //!< Mesh class
unsigned char format; //!< Data format (1: old, 2: new, 3: new (single), 4: silo)
FileFormat format; //!< Data format (1: old, 2: new, 3: new (single), 4: silo)
std::vector<DatabaseEntry> domains; //!< List of the domains
std::vector<VariableDatabase> variables; //!< List of the variables
std::map<variable_id, DatabaseEntry> variable_data; //!< Data for the variables

18
IO/Reader.cpp
View File

@ -155,7 +155,7 @@ std::shared_ptr<IO::Mesh> IO::getMesh( const std::string &path, const std::strin
{
PROFILE_START( "getMesh" );
std::shared_ptr<IO::Mesh> mesh;
if ( meshDatabase.format == 1 ) {
if ( meshDatabase.format == FileFormat::OLD ) {
// Old format (binary doubles)
std::string filename = path + "/" + timestep + "/" + meshDatabase.domains[domain].file;
FILE *fid = fopen( filename.c_str(), "rb" );
@ -206,7 +206,8 @@ std::shared_ptr<IO::Mesh> IO::getMesh( const std::string &path, const std::strin
ERROR( "Unknown mesh type" );
}
delete[] data;
} else if ( meshDatabase.format == 2 ) {
} else if ( meshDatabase.format == FileFormat::NEW ||
meshDatabase.format == FileFormat::NEW_SINGLE ) {
const DatabaseEntry &database = meshDatabase.domains[domain];
std::string filename = path + "/" + timestep + "/" + database.file;
FILE *fid = fopen( filename.c_str(), "rb" );
@ -233,7 +234,7 @@ std::shared_ptr<IO::Mesh> IO::getMesh( const std::string &path, const std::strin
}
mesh->unpack( std::pair<size_t, void *>( bytes, data ) );
delete[] data;
} else if ( meshDatabase.format == 4 ) {
} else if ( meshDatabase.format == FileFormat::SILO ) {
// Reading a silo file
#ifdef USE_SILO
const DatabaseEntry &database = meshDatabase.domains[domain];
@ -301,12 +302,11 @@ std::shared_ptr<IO::Variable> IO::getVariable( const std::string &path, const st
const MeshDatabase &meshDatabase, int domain, const std::string &variable )
{
std::pair<std::string, std::string> key( meshDatabase.domains[domain].name, variable );
std::map<std::pair<std::string, std::string>, DatabaseEntry>::const_iterator it;
it = meshDatabase.variable_data.find( key );
auto it = meshDatabase.variable_data.find( key );
if ( it == meshDatabase.variable_data.end() )
return std::shared_ptr<IO::Variable>();
std::shared_ptr<IO::Variable> var;
if ( meshDatabase.format == 2 ) {
if ( meshDatabase.format == FileFormat::NEW || meshDatabase.format == FileFormat::NEW_SINGLE ) {
const DatabaseEntry &database = it->second;
std::string filename = path + "/" + timestep + "/" + database.file;
FILE *fid = fopen( filename.c_str(), "rb" );
@ -318,13 +318,13 @@ std::shared_ptr<IO::Variable> IO::getVariable( const std::string &path, const st
std::vector<std::string> values = splitList( &line[i2 + 1], ',' );
ASSERT( values.size() == 5 );
int dim = atoi( values[0].c_str() );
int type = atoi( values[1].c_str() );
auto type = values[1];
size_t N = atol( values[2].c_str() );
size_t bytes = atol( values[3].c_str() );
std::string precision = values[4];
var = std::shared_ptr<IO::Variable>( new IO::Variable() );
var->dim = dim;
var->type = static_cast<IO::VariableType>( type );
var->type = getVariableType( type );
var->name = variable;
var->data.resize( N, dim );
if ( precision == "double" ) {
@ -334,7 +334,7 @@ std::shared_ptr<IO::Variable> IO::getVariable( const std::string &path, const st
ERROR( "Format not implimented" );
}
fclose( fid );
} else if ( meshDatabase.format == 4 ) {
} else if ( meshDatabase.format == FileFormat::SILO ) {
// Reading a silo file
#ifdef USE_SILO
const auto &database = meshDatabase.domains[domain];

90
IO/Writer.cpp
View File

@ -110,7 +110,7 @@ static std::vector<IO::MeshDatabase> writeMeshesOrigFormat(
mesh_entry.name = meshData[i].meshName;
mesh_entry.type = meshType( *mesh );
mesh_entry.meshClass = meshData[i].mesh->className();
mesh_entry.format = 1;
mesh_entry.format = IO::FileFormat::OLD;
IO::DatabaseEntry domain;
domain.name = domainname;
domain.file = filename;
@ -171,7 +171,7 @@ static std::vector<IO::MeshDatabase> writeMeshesOrigFormat(
// Create the database entry for the mesh data
static IO::MeshDatabase getDatabase(
const std::string &filename, const IO::MeshDataStruct &mesh, int format, int rank )
const std::string &filename, const IO::MeshDataStruct &mesh, IO::FileFormat format, int rank )
{
char domainname[100];
sprintf( domainname, "%s_%05i", mesh.meshName.c_str(), rank );
@ -209,8 +209,8 @@ static IO::MeshDatabase getDatabase(
// Write a mesh (and variables) to a file
static IO::MeshDatabase write_domain(
FILE *fid, const std::string &filename, const IO::MeshDataStruct &mesh, int format, int rank )
static IO::MeshDatabase write_domain( FILE *fid, const std::string &filename,
const IO::MeshDataStruct &mesh, IO::FileFormat format, int rank )
{
const int level = 0;
// Create the MeshDatabase
@ -225,19 +225,17 @@ static IO::MeshDatabase write_domain(
delete[]( char * ) data.second;
// Write the variables
for ( size_t i = 0; i < mesh.vars.size(); i++ ) {
ASSERT( mesh.vars[i]->type != IO::VariableType::NullVariable );
std::pair<std::string, std::string> key( domain.name, mesh.vars[i]->name );
IO::DatabaseEntry &variable = database.variable_data[key];
variable.offset = ftell( fid );
int dim = mesh.vars[i]->dim;
int type = static_cast<int>( mesh.vars[i]->type );
size_t N = mesh.vars[i]->data.length();
if ( type == static_cast<int>( IO::VariableType::NullVariable ) ) {
ERROR( "Variable type not set" );
}
size_t N_mesh = mesh.mesh->numberPointsVar( mesh.vars[i]->type );
auto &variable = database.variable_data[key];
variable.offset = ftell( fid );
int dim = mesh.vars[i]->dim;
auto type = getString( mesh.vars[i]->type );
size_t N = mesh.vars[i]->data.length();
size_t N_mesh = mesh.mesh->numberPointsVar( mesh.vars[i]->type );
ASSERT( N == dim * N_mesh );
fprintf( fid, "Var: %s-%05i-%s: %i, %i, %lu, %lu, double\n", database.name.c_str(), rank,
variable.name.c_str(), dim, type, N_mesh, N * sizeof( double ) );
fprintf( fid, "Var: %s-%05i-%s: %i, %s, %lu, %lu, double\n", database.name.c_str(), rank,
variable.name.c_str(), dim, type.data(), N_mesh, N * sizeof( double ) );
fwrite( mesh.vars[i]->data.data(), sizeof( double ), N, fid );
fprintf( fid, "\n" );
}
@ -259,7 +257,7 @@ static void writeSiloPointMesh(
z[i] = points[i].z;
}
const TYPE *coords[] = { x.data(), y.data(), z.data() };
silo::writePointMesh<TYPE>( fid, meshname, 3, points.size(), coords );
IO::silo::writePointMesh<TYPE>( fid, meshname, 3, points.size(), coords );
}
static void writeSiloPointList(
DBfile *fid, const IO::MeshDataStruct &meshData, IO::MeshDatabase database )
@ -281,19 +279,19 @@ static void writeSiloPointList(
z[i] = points[i].z;
}
const double *coords[] = { x.data(), y.data(), z.data() };
silo::writePointMesh( fid, meshname, 3, points.size(), coords );
IO::silo::writePointMesh( fid, meshname, 3, points.size(), coords );
for ( size_t i = 0; i < meshData.vars.size(); i++ ) {
const IO::Variable &var = *meshData.vars[i];
if ( var.precision == IO::DataType::Double ) {
silo::writePointMeshVariable( fid, meshname, var.name, var.data );
IO::silo::writePointMeshVariable( fid, meshname, var.name, var.data );
} else if ( var.precision == IO::DataType::Float ) {
Array<float> data2( var.data.size() );
data2.copy( var.data );
silo::writePointMeshVariable( fid, meshname, var.name, data2 );
IO::silo::writePointMeshVariable( fid, meshname, var.name, data2 );
} else if ( var.precision == IO::DataType::Int ) {
Array<int> data2( var.data.size() );
data2.copy( var.data );
silo::writePointMeshVariable( fid, meshname, var.name, data2 );
IO::silo::writePointMeshVariable( fid, meshname, var.name, data2 );
} else {
ERROR( "Unsupported format" );
}
@ -312,7 +310,7 @@ static void writeSiloTriMesh( DBfile *fid, const IO::TriMesh &mesh, const std::s
}
const TYPE *coords[] = { x.data(), y.data(), z.data() };
const int *tri[] = { mesh.A.data(), mesh.B.data(), mesh.C.data() };
silo::writeTriMesh<TYPE>( fid, meshname, 3, 2, points.size(), coords, mesh.A.size(), tri );
IO::silo::writeTriMesh<TYPE>( fid, meshname, 3, 2, points.size(), coords, mesh.A.size(), tri );
}
static void writeSiloTriMesh2( DBfile *fid, const IO::MeshDataStruct &meshData,
const IO::TriMesh &mesh, IO::MeshDatabase database )
@ -327,17 +325,16 @@ static void writeSiloTriMesh2( DBfile *fid, const IO::MeshDataStruct &meshData,
}
for ( size_t i = 0; i < meshData.vars.size(); i++ ) {
const IO::Variable &var = *meshData.vars[i];
auto type = static_cast<silo::VariableType>( var.type );
if ( var.precision == IO::DataType::Double ) {
silo::writeTriMeshVariable( fid, 3, meshname, var.name, var.data, type );
IO::silo::writeTriMeshVariable( fid, 3, meshname, var.name, var.data, var.type );
} else if ( var.precision == IO::DataType::Float ) {
Array<float> data2( var.data.size() );
data2.copy( var.data );
silo::writeTriMeshVariable( fid, 3, meshname, var.name, data2, type );
IO::silo::writeTriMeshVariable( fid, 3, meshname, var.name, data2, var.type );
} else if ( var.precision == IO::DataType::Int ) {
Array<int> data2( var.data.size() );
data2.copy( var.data );
silo::writeTriMeshVariable( fid, 3, meshname, var.name, data2, type );
IO::silo::writeTriMeshVariable( fid, 3, meshname, var.name, data2, var.type );
} else {
ERROR( "Unsupported format" );
}
@ -367,30 +364,29 @@ static void writeSiloDomainMesh(
( info.kz + 1 ) * mesh.Lz / info.nz };
std::array<int, 3> N = { mesh.nx, mesh.ny, mesh.nz };
auto meshname = database.domains[0].name;
silo::writeUniformMesh<3>( fid, meshname, range, N );
silo::write<int>(
IO::silo::writeUniformMesh<3>( fid, meshname, range, N );
IO::silo::write<int>(
fid, meshname + "_rankinfo", { mesh.rank, mesh.nprocx, mesh.nprocy, mesh.nprocz } );
for ( size_t i = 0; i < meshData.vars.size(); i++ ) {
const auto &var = *meshData.vars[i];
auto type = static_cast<silo::VariableType>( var.type );
if ( var.precision == IO::DataType::Double ) {
silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, var.data, type );
IO::silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, var.data, var.type );
} else if ( var.precision == IO::DataType::Float ) {
Array<float> data2( var.data.size() );
data2.copy( var.data );
silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, data2, type );
IO::silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, data2, var.type );
} else if ( var.precision == IO::DataType::Int ) {
Array<int> data2( var.data.size() );
data2.copy( var.data );
silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, data2, type );
IO::silo::writeUniformMeshVariable<3>( fid, meshname, N, var.name, data2, var.type );
} else {
ERROR( "Unsupported format" );
}
}
}
// Write a mesh (and variables) to a file
static IO::MeshDatabase write_domain_silo(
DBfile *fid, const std::string &filename, const IO::MeshDataStruct &mesh, int format, int rank )
static IO::MeshDatabase write_domain_silo( DBfile *fid, const std::string &filename,
const IO::MeshDataStruct &mesh, IO::FileFormat format, int rank )
{
// Create the MeshDatabase
auto database = getDatabase( filename, mesh, format, rank );
@ -432,7 +428,7 @@ std::pair<int, int> getSiloMeshType( const std::string &meshClass )
void writeSiloSummary(
const std::vector<IO::MeshDatabase> &meshes_written, const std::string &filename )
{
auto fid = silo::open( filename, silo::CREATE );
auto fid = IO::silo::open( filename, IO::silo::CREATE );
for ( const auto &data : meshes_written ) {
auto type = getSiloMeshType( data.meshClass );
std::vector<int> meshTypes( data.domains.size(), type.first );
@ -440,22 +436,23 @@ void writeSiloSummary(
std::vector<std::string> meshNames;
for ( const auto &tmp : data.domains )
meshNames.push_back( tmp.file + ":" + tmp.name );
silo::writeMultiMesh( fid, data.name, meshNames, meshTypes );
IO::silo::writeMultiMesh( fid, data.name, meshNames, meshTypes );
for ( const auto &variable : data.variables ) {
std::vector<std::string> varnames;
for ( const auto &tmp : data.domains )
varnames.push_back( tmp.file + ":" + variable.name );
silo::writeMultiVar( fid, variable.name, varnames, varTypes );
IO::silo::writeMultiVar( fid, variable.name, varnames, varTypes );
}
}
silo::close( fid );
IO::silo::close( fid );
}
#endif
// Write the mesh data in the new format
static std::vector<IO::MeshDatabase> writeMeshesNewFormat(
const std::vector<IO::MeshDataStruct> &meshData, const std::string &path, int format, int rank )
const std::vector<IO::MeshDataStruct> &meshData, const std::string &path, IO::FileFormat format,
int rank )
{
std::vector<IO::MeshDatabase> meshes_written;
char filename[100], fullpath[200];
@ -473,19 +470,20 @@ static std::vector<IO::MeshDatabase> writeMeshesNewFormat(
// Write the mesh data to silo
static std::vector<IO::MeshDatabase> writeMeshesSilo(
const std::vector<IO::MeshDataStruct> &meshData, const std::string &path, int format, int rank )
const std::vector<IO::MeshDataStruct> &meshData, const std::string &path, IO::FileFormat format,
int rank )
{
#ifdef USE_SILO
std::vector<IO::MeshDatabase> meshes_written;
char filename[100], fullpath[200];
sprintf( filename, "%05i.silo", rank );
sprintf( fullpath, "%s/%s", path.c_str(), filename );
auto fid = silo::open( fullpath, silo::CREATE );
auto fid = IO::silo::open( fullpath, IO::silo::CREATE );
for ( size_t i = 0; i < meshData.size(); i++ ) {
auto mesh = meshData[i].mesh;
meshes_written.push_back( write_domain_silo( fid, filename, meshData[i], format, rank ) );
}
silo::close( fid );
IO::silo::close( fid );
return meshes_written;
#else
NULL_USE( meshData );
@ -509,10 +507,8 @@ void IO::writeData( const std::string &subdir, const std::vector<IO::MeshDataStr
PROFILE_START( "writeData" );
int rank = Utilities::MPI( MPI_COMM_WORLD ).getRank();
// Check the meshData before writing
for ( const auto &data : meshData ) {
if ( !data.check() )
ERROR( "Error in meshData" );
}
for ( const auto &data : meshData )
ASSERT( data.check() );
// Create the output directory
std::string path = global_IO_path + "/" + subdir;
recursiveMkdir( path, S_IRWXU | S_IRGRP );
@ -523,10 +519,10 @@ void IO::writeData( const std::string &subdir, const std::vector<IO::MeshDataStr
meshes_written = writeMeshesOrigFormat( meshData, path, rank );
} else if ( global_IO_format == Format::NEW ) {
// Write the new format (double precision)
meshes_written = writeMeshesNewFormat( meshData, path, 2, rank );
meshes_written = writeMeshesNewFormat( meshData, path, IO::FileFormat::NEW, rank );
} else if ( global_IO_format == Format::SILO ) {
// Write silo
meshes_written = writeMeshesSilo( meshData, path, 4, rank );
meshes_written = writeMeshesSilo( meshData, path, IO::FileFormat::SILO, rank );
} else {
ERROR( "Unknown format" );
}

16
IO/silo.cpp
View File

@ -10,7 +10,7 @@
#include <silo.h>
namespace silo {
namespace IO::silo {
/****************************************************
@ -34,16 +34,16 @@ void close( DBfile *fid ) { DBClose( fid ); }
/****************************************************
* Helper functions *
****************************************************/
VariableDataType varDataType( DBfile *fid, const std::string &name )
DataType varDataType( DBfile *fid, const std::string &name )
{
auto type = DBGetVarType( fid, name.c_str() );
VariableDataType type2 = VariableDataType::UNKNOWN;
auto type = DBGetVarType( fid, name.c_str() );
DataType type2 = DataType::Null;
if ( type == DB_DOUBLE )
type2 = VariableDataType::DOUBLE;
type2 = DataType::Double;
else if ( type == DB_FLOAT )
type2 = VariableDataType::FLOAT;
type2 = DataType::Float;
else if ( type == DB_INT )
type2 = VariableDataType::INT;
type2 = DataType::Int;
return type2;
}
@ -99,7 +99,7 @@ void writeMultiVar( DBfile *fid, const std::string &varname,
}
}; // namespace silo
}; // namespace IO::silo
#else

17
IO/silo.h
View File

@ -5,6 +5,7 @@
#include <string>
#include <vector>
#include "IO/Mesh.h"
#include "common/Array.h"
#include "common/Communication.h"
#include "common/MPI.h"
@ -17,21 +18,11 @@ typedef int DBfile;
#endif
namespace silo {
namespace IO::silo {
enum FileMode { READ, WRITE, CREATE };
enum class VariableType : int {
NodeVariable = 1,
EdgeVariable = 2,
SurfaceVariable = 2,
VolumeVariable = 3,
NullVariable = 0
};
enum class VariableDataType { DOUBLE, FLOAT, INT, UNKNOWN };
/*!
* @brief Open silo file
@ -57,7 +48,7 @@ void close( DBfile *fid );
* @param[in] fid Handle to the open file
* @param[in] name Name of variable
*/
VariableDataType varDataType( DBfile *dbfile, const std::string &name );
DataType varDataType( DBfile *dbfile, const std::string &name );
/*!
@ -265,7 +256,7 @@ void writeMultiVar( DBfile *fid, const std::string &varname,
const std::vector<std::string> &subVarNames, const std::vector<int> &subVarTypes );
}; // namespace silo
}; // namespace IO::silo
#endif
#include "IO/silo.hpp"

4
IO/silo.hpp
View File

@ -13,7 +13,7 @@
#include <silo.h>
namespace silo {
namespace IO::silo {
/****************************************************
@ -413,7 +413,7 @@ Array<TYPE> readTriMeshVariable( DBfile *fid, const std::string &varname )
}
}; // namespace silo
}; // namespace IO::silo
#endif

400
tests/TestWriter.cpp
View File

@ -1,115 +1,117 @@
#include <exception>
#include <fstream>
#include <iostream>
#include <memory>
#include <stdexcept>
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include <exception>
#include <stdexcept>
#include <fstream>
#include <memory>
#include "common/UnitTest.h"
#include "common/Utilities.h"
#include "common/MPI.h"
#include "IO/MeshDatabase.h"
#include "IO/Reader.h"
#include "IO/Writer.h"
#include "ProfilerApp.h"
#include "common/MPI.h"
#include "common/UnitTest.h"
#include "common/Utilities.h"
inline bool approx_equal( const Point& A, const Point& B )
inline bool approx_equal( const Point &A, const Point &B )
{
double tol = 1e-7*sqrt(A.x*A.x+A.y*A.y+A.z*A.z);
return fabs(A.x-B.x)<=tol && fabs(A.y-B.y)<=tol && fabs(A.z-B.z)<=tol;
double tol = 1e-7 * sqrt( A.x * A.x + A.y * A.y + A.z * A.z );
return fabs( A.x - B.x ) <= tol && fabs( A.y - B.y ) <= tol && fabs( A.z - B.z ) <= tol;
}
inline bool approx_equal( const double& A, const double& B )
inline bool approx_equal( const double &A, const double &B )
{
return fabs(A-B) <= std::max<double>(1e-7*fabs(A+B),1e-20);
return fabs( A - B ) <= std::max<double>( 1e-7 * fabs( A + B ), 1e-20 );
}
inline double distance( const Point& p )
{
return sqrt(p.x*p.x+p.y*p.y+p.z*p.z);
}
inline double distance( const Point &p ) { return sqrt( p.x * p.x + p.y * p.y + p.z * p.z ); }
bool checkMesh( const std::vector<IO::MeshDataStruct>& meshData, const std::string& format, std::shared_ptr<IO::Mesh> mesh )
bool checkMesh( const std::vector<IO::MeshDataStruct> &meshData, const std::string &format,
std::shared_ptr<IO::Mesh> mesh )
{
// Get direct access to the meshes used to test the reader
const auto pointmesh = dynamic_cast<IO::PointList*>( meshData[0].mesh.get() );
const auto trimesh = dynamic_cast<IO::TriMesh*>( meshData[1].mesh.get() );
const auto trilist = dynamic_cast<IO::TriList*>( meshData[2].mesh.get() );
const auto domain = dynamic_cast<IO::DomainMesh*>( meshData[3].mesh.get() );
const size_t N_tri = trimesh->A.size();
const auto pointmesh = dynamic_cast<IO::PointList *>( meshData[0].mesh.get() );
const auto trimesh = dynamic_cast<IO::TriMesh *>( meshData[1].mesh.get() );
const auto trilist = dynamic_cast<IO::TriList *>( meshData[2].mesh.get() );
const auto domain = dynamic_cast<IO::DomainMesh *>( meshData[3].mesh.get() );
const size_t N_tri = trimesh->A.size();
if ( mesh->className() == "pointmesh" ) {
// Check the pointmesh
auto pmesh = IO::getPointList(mesh);
if ( pmesh.get()==NULL )
auto pmesh = IO::getPointList( mesh );
if ( pmesh.get() == NULL )
return false;
if ( pmesh->points.size() != pointmesh->points.size() )
return false;
}
if ( mesh->className() == "trimesh" || mesh->className() == "trilist" ) {
// Check the trimesh/trilist
auto mesh1 = IO::getTriMesh(mesh);
auto mesh2 = IO::getTriList(mesh);
if ( mesh1.get()==NULL || mesh2.get()==NULL )
return false;
if ( mesh1->A.size()!=N_tri || mesh1->B.size()!=N_tri || mesh1->C.size()!=N_tri ||
mesh2->A.size()!=N_tri || mesh2->B.size()!=N_tri || mesh2->C.size()!=N_tri )
return false;
const std::vector<Point>& P1 = mesh1->vertices->points;
const std::vector<int>& A1 = mesh1->A;
const std::vector<int>& B1 = mesh1->B;
const std::vector<int>& C1 = mesh1->C;
const std::vector<Point>& A2 = mesh2->A;
const std::vector<Point>& B2 = mesh2->B;
const std::vector<Point>& C2 = mesh2->C;
const std::vector<Point>& A = trilist->A;
const std::vector<Point>& B = trilist->B;
const std::vector<Point>& C = trilist->C;
for (size_t i=0; i<N_tri; i++) {
if ( !approx_equal(P1[A1[i]],A[i]) || !approx_equal(P1[B1[i]],B[i]) || !approx_equal(P1[C1[i]],C[i]) )
auto mesh1 = IO::getTriMesh( mesh );
auto mesh2 = IO::getTriList( mesh );
if ( mesh1.get() == NULL || mesh2.get() == NULL )
return false;
if ( mesh1->A.size() != N_tri || mesh1->B.size() != N_tri || mesh1->C.size() != N_tri ||
mesh2->A.size() != N_tri || mesh2->B.size() != N_tri || mesh2->C.size() != N_tri )
return false;
const std::vector<Point> &P1 = mesh1->vertices->points;
const std::vector<int> &A1 = mesh1->A;
const std::vector<int> &B1 = mesh1->B;
const std::vector<int> &C1 = mesh1->C;
const std::vector<Point> &A2 = mesh2->A;
const std::vector<Point> &B2 = mesh2->B;
const std::vector<Point> &C2 = mesh2->C;
const std::vector<Point> &A = trilist->A;
const std::vector<Point> &B = trilist->B;
const std::vector<Point> &C = trilist->C;
for ( size_t i = 0; i < N_tri; i++ ) {
if ( !approx_equal( P1[A1[i]], A[i] ) || !approx_equal( P1[B1[i]], B[i] ) ||
!approx_equal( P1[C1[i]], C[i] ) )
return false;
if ( !approx_equal(A2[i],A[i]) || !approx_equal(B2[i],B[i]) || !approx_equal(C2[i],C[i]) )
if ( !approx_equal( A2[i], A[i] ) || !approx_equal( B2[i], B[i] ) ||
!approx_equal( C2[i], C[i] ) )
return false;
}
}
if ( mesh->className() == "domain" && format!="old" ) {
if ( mesh->className() == "domain" && format != "old" ) {
// Check the domain mesh
const IO::DomainMesh& mesh1 = *std::dynamic_pointer_cast<IO::DomainMesh>(mesh);
if ( mesh1.nprocx!=domain->nprocx || mesh1.nprocy!=domain->nprocy || mesh1.nprocz!=domain->nprocz )
const IO::DomainMesh &mesh1 = *std::dynamic_pointer_cast<IO::DomainMesh>( mesh );
if ( mesh1.nprocx != domain->nprocx || mesh1.nprocy != domain->nprocy ||
mesh1.nprocz != domain->nprocz )
return false;
if ( mesh1.nx!=domain->nx || mesh1.ny!=domain->ny || mesh1.nz!=domain->nz )
if ( mesh1.nx != domain->nx || mesh1.ny != domain->ny || mesh1.nz != domain->nz )
return false;
if ( mesh1.Lx!=domain->Lx || mesh1.Ly!=domain->Ly || mesh1.Lz!=domain->Lz )
if ( mesh1.Lx != domain->Lx || mesh1.Ly != domain->Ly || mesh1.Lz != domain->Lz )
return false;
}
return true;
}
bool checkVar( const std::string& format, std::shared_ptr<IO::Mesh> mesh,
bool checkVar( const std::string &format, std::shared_ptr<IO::Mesh> mesh,
std::shared_ptr<IO::Variable> variable1, std::shared_ptr<IO::Variable> variable2 )
{
if ( format=="new" )
if ( format == "new" )
IO::reformatVariable( *mesh, *variable2 );
bool pass = true;
const IO::Variable& var1 = *variable1;
const IO::Variable& var2 = *variable2;
pass = var1.name == var2.name;
pass = pass && var1.dim == var2.dim;
pass = pass && var1.type == var2.type;
pass = pass && var1.data.length() == var2.data.length();
bool pass = true;
const IO::Variable &var1 = *variable1;
const IO::Variable &var2 = *variable2;
pass = var1.name == var2.name;
pass = pass && var1.dim == var2.dim;
pass = pass && var1.type == var2.type;
pass = pass && var1.data.length() == var2.data.length();
if ( pass ) {
for (size_t m=0; m<var1.data.length(); m++)
pass = pass && approx_equal(var1.data(m),var2.data(m));
for ( size_t m = 0; m < var1.data.length(); m++ )
pass = pass && approx_equal( var1.data( m ), var2.data( m ) );
}
return pass;
}
// Test writing and reading the given format
void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& meshData, UnitTest& ut )
void testWriter(
const std::string &format, std::vector<IO::MeshDataStruct> &meshData, UnitTest &ut )
{
PROFILE_SCOPED( path, 0, timer );
@ -124,20 +126,20 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
// Get the format
std::string format2 = format;
auto precision = IO::DataType::Double;
auto precision = IO::DataType::Double;
if ( format == "silo-double" ) {
format2 = "silo";
format2 = "silo";
precision = IO::DataType::Double;
} else if ( format == "silo-float" ) {
format2 = "silo";
format2 = "silo";
precision = IO::DataType::Float;
}
// Set the precision for the variables
for ( auto& data : meshData ) {
for ( auto &data : meshData ) {
data.precision = precision;
for ( auto& var : data.vars )
for ( auto &var : data.vars )
var->precision = precision;
}
@ -150,18 +152,18 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
// Get a list of the timesteps
auto timesteps = IO::readTimesteps( path, format2 );
if ( timesteps.size()==2 )
ut.passes(format+": Corrent number of timesteps");
if ( timesteps.size() == 2 )
ut.passes( format + ": Corrent number of timesteps" );
else
ut.failure(format+": Incorrent number of timesteps");
ut.failure( format + ": Incorrent number of timesteps" );
// Test the simple read interface
bool pass = true;
for ( const auto& timestep : timesteps ) {
for ( const auto &timestep : timesteps ) {
auto data = IO::readData( path, timestep, comm.getRank() );
pass = pass && data.size() == meshData.size();
for ( size_t i=0; i<data.size(); i++ ) {
pass = pass && data.size() == meshData.size();
for ( size_t i = 0; i < data.size(); i++ ) {
pass = pass && checkMesh( meshData, format, data[i].mesh );
}
}
@ -172,15 +174,15 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
// Test reading each mesh domain
for ( const auto& timestep : timesteps ) {
for ( const auto &timestep : timesteps ) {
// Load the list of meshes and check its size
auto databaseList = IO::getMeshList(path,timestep);
if ( databaseList.size()==meshData.size() )
ut.passes(format+": Corrent number of meshes found");
auto databaseList = IO::getMeshList( path, timestep );
if ( databaseList.size() == meshData.size() )
ut.passes( format + ": Corrent number of meshes found" );
else
ut.failure(format+": Incorrent number of meshes found");
ut.failure( format + ": Incorrent number of meshes found" );
// Check the number of domains for each mesh
for ( const auto& database : databaseList ) {
for ( const auto &database : databaseList ) {
int N_domains = database.domains.size();
if ( N_domains != nprocs ) {
ut.failure( format + ": Incorrent number of domains for mesh" );
@ -188,8 +190,8 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
}
// For each domain, load the mesh and check its data
bool pass = true;
for (int k=0; k<N_domains; k++) {
auto mesh = IO::getMesh(path,timestep,database,k);
for ( int k = 0; k < N_domains; k++ ) {
auto mesh = IO::getMesh( path, timestep, database, k );
if ( !mesh ) {
ut.failure( "Failed to load " + database.name );
pass = false;
@ -198,31 +200,33 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
}
}
if ( pass ) {
ut.passes(format+": Mesh \"" + database.name + "\" loaded correctly");
ut.passes( format + ": Mesh \"" + database.name + "\" loaded correctly" );
} else {
ut.failure(format+": Mesh \"" + database.name + "\" did not load correctly");
ut.failure( format + ": Mesh \"" + database.name + "\" did not load correctly" );
continue;
}
// Load the variables and check their data
if ( format=="old" )
continue; // Old format does not support variables
const IO::MeshDataStruct* mesh0 = nullptr;
for (size_t k=0; k<meshData.size(); k++) {
if ( format == "old" )
continue; // Old format does not support variables
const IO::MeshDataStruct *mesh0 = nullptr;
for ( size_t k = 0; k < meshData.size(); k++ ) {
if ( meshData[k].meshName == database.name ) {
mesh0 = &meshData[k];
break;
}
}
for (int k=0; k<N_domains; k++) {
auto mesh = IO::getMesh(path,timestep,database,k);
for (size_t v=0; v<mesh0->vars.size(); v++) {
PROFILE_START(format+"-read-getVariable");
auto variable = IO::getVariable(path,timestep,database,k,mesh0->vars[v]->name);
for ( int k = 0; k < N_domains; k++ ) {
auto mesh = IO::getMesh( path, timestep, database, k );
for ( size_t v = 0; v < mesh0->vars.size(); v++ ) {
PROFILE_START( format + "-read-getVariable" );
auto variable =
IO::getVariable( path, timestep, database, k, mesh0->vars[v]->name );
pass = checkVar( format, mesh, mesh0->vars[v], variable );
if ( pass ) {
ut.passes(format+": Variable \"" + variable->name + "\" matched");
ut.passes( format + ": Variable \"" + variable->name + "\" matched" );
} else {
ut.failure(format+": Variable \"" + variable->name + "\" did not match");
ut.failure(
format + ": Variable \"" + variable->name + "\" did not match" );
break;
}
}
@ -233,157 +237,161 @@ void testWriter( const std::string& format, std::vector<IO::MeshDataStruct>& mes
// Main
int main(int argc, char **argv)
int main( int argc, char **argv )
{
Utilities::startup( argc, argv );
Utilities::MPI comm( MPI_COMM_WORLD );
int rank = comm.getRank();
int rank = comm.getRank();
int nprocs = comm.getSize();
Utilities::setAbortBehavior(true,2);
Utilities::setAbortBehavior( true, 2 );
Utilities::setErrorHandlers();
UnitTest ut;
// Create some points
const int N_points = 8;
const int N_tri = 12;
double x[8] = { 0, 1, 0, 1, 0, 1, 0, 1 };
double y[8] = { 0, 0, 1, 1, 0, 0, 1, 1 };
double z[8] = { 0, 0, 0, 0, 1, 1, 1, 1 };
int tri[N_tri][3] = {
{0,1,3}, {0,3,2}, {4,5,7}, {4,7,6}, // z faces
{0,1,4}, {1,4,5}, {2,3,6}, {3,6,7}, // y faces
{0,2,4}, {2,4,6}, {1,3,5}, {3,5,7} // x faces
const int N_tri = 12;
double x[8] = { 0, 1, 0, 1, 0, 1, 0, 1 };
double y[8] = { 0, 0, 1, 1, 0, 0, 1, 1 };
double z[8] = { 0, 0, 0, 0, 1, 1, 1, 1 };
int tri[N_tri][3] = {
{ 0, 1, 3 }, { 0, 3, 2 }, { 4, 5, 7 }, { 4, 7, 6 }, // z faces
{ 0, 1, 4 }, { 1, 4, 5 }, { 2, 3, 6 }, { 3, 6, 7 }, // y faces
{ 0, 2, 4 }, { 2, 4, 6 }, { 1, 3, 5 }, { 3, 5, 7 } // x faces
};
// Create the meshes
auto set1 = std::make_shared<IO::PointList>(N_points);
for (int i=0; i<N_points; i++) {
auto set1 = std::make_shared<IO::PointList>( N_points );
for ( int i = 0; i < N_points; i++ ) {
set1->points[i].x = x[i];
set1->points[i].y = y[i];
set1->points[i].z = z[i];
}
auto trimesh = std::make_shared<IO::TriMesh>(N_tri,set1);
for (int i=0; i<N_tri; i++) {
auto trimesh = std::make_shared<IO::TriMesh>( N_tri, set1 );
for ( int i = 0; i < N_tri; i++ ) {
trimesh->A[i] = tri[i][0];
trimesh->B[i] = tri[i][1];
trimesh->C[i] = tri[i][2];
}
auto trilist = std::make_shared<IO::TriList>(*trimesh);
for (int i=0; i<N_tri; i++) {
Point A(x[tri[i][0]],y[tri[i][0]],z[tri[i][0]]);
Point B(x[tri[i][1]],y[tri[i][1]],z[tri[i][1]]);
Point C(x[tri[i][2]],y[tri[i][2]],z[tri[i][2]]);
if ( !approx_equal(trilist->A[i],A) || !approx_equal(trilist->B[i],B) || !approx_equal(trilist->C[i],C) )
{
printf("Failed to create trilist\n");
auto trilist = std::make_shared<IO::TriList>( *trimesh );
for ( int i = 0; i < N_tri; i++ ) {
Point A( x[tri[i][0]], y[tri[i][0]], z[tri[i][0]] );
Point B( x[tri[i][1]], y[tri[i][1]], z[tri[i][1]] );
Point C( x[tri[i][2]], y[tri[i][2]], z[tri[i][2]] );
if ( !approx_equal( trilist->A[i], A ) || !approx_equal( trilist->B[i], B ) ||
!approx_equal( trilist->C[i], C ) ) {
printf( "Failed to create trilist\n" );
return -1;
}
}
RankInfoStruct rank_data( rank, nprocs, 1, 1 );
auto domain = std::make_shared<IO::DomainMesh>(rank_data,6,7,8,1.0,1.0,1.0);
auto domain = std::make_shared<IO::DomainMesh>( rank_data, 6, 7, 8, 1.0, 1.0, 1.0 );
// Create the variables
const auto NodeVar = IO::VariableType::NodeVariable;
const auto VolVar = IO::VariableType::VolumeVariable;
auto set_node_mag = std::make_shared<IO::Variable>(1,NodeVar,"Node_set_mag");
auto set_node_vec = std::make_shared<IO::Variable>(3,NodeVar,"Node_set_vec");
auto list_node_mag = std::make_shared<IO::Variable>(1,NodeVar,"Node_list_mag");
auto list_node_vec = std::make_shared<IO::Variable>(3,NodeVar,"Node_list_vec");
auto point_node_mag = std::make_shared<IO::Variable>(1,NodeVar,"Node_point_mag");
auto point_node_vec = std::make_shared<IO::Variable>(3,NodeVar,"Node_point_vec");
auto domain_node_mag = std::make_shared<IO::Variable>(1,NodeVar,"Node_domain_mag");
auto domain_node_vec = std::make_shared<IO::Variable>(3,NodeVar,"Node_domain_vec");
auto set_cell_mag = std::make_shared<IO::Variable>(1,VolVar,"Cell_set_mag");
auto set_cell_vec = std::make_shared<IO::Variable>(3,VolVar,"Cell_set_vec");
auto list_cell_mag = std::make_shared<IO::Variable>(1,VolVar,"Cell_list_mag");
auto list_cell_vec = std::make_shared<IO::Variable>(3,VolVar,"Cell_list_vec");
auto domain_cell_mag = std::make_shared<IO::Variable>(1,VolVar,"Cell_domain_mag");
auto domain_cell_vec = std::make_shared<IO::Variable>(3,VolVar,"Cell_domain_vec");
auto set_node_mag = std::make_shared<IO::Variable>( 1, NodeVar, "Node_set_mag" );
auto set_node_vec = std::make_shared<IO::Variable>( 3, NodeVar, "Node_set_vec" );
auto list_node_mag = std::make_shared<IO::Variable>( 1, NodeVar, "Node_list_mag" );
auto list_node_vec = std::make_shared<IO::Variable>( 3, NodeVar, "Node_list_vec" );
auto point_node_mag = std::make_shared<IO::Variable>( 1, NodeVar, "Node_point_mag" );
auto point_node_vec = std::make_shared<IO::Variable>( 3, NodeVar, "Node_point_vec" );
auto domain_node_mag = std::make_shared<IO::Variable>( 1, NodeVar, "Node_domain_mag" );
auto domain_node_vec = std::make_shared<IO::Variable>( 3, NodeVar, "Node_domain_vec" );
auto set_cell_mag = std::make_shared<IO::Variable>( 1, VolVar, "Cell_set_mag" );
auto set_cell_vec = std::make_shared<IO::Variable>( 3, VolVar, "Cell_set_vec" );
auto list_cell_mag = std::make_shared<IO::Variable>( 1, VolVar, "Cell_list_mag" );
auto list_cell_vec = std::make_shared<IO::Variable>( 3, VolVar, "Cell_list_vec" );
auto domain_cell_mag = std::make_shared<IO::Variable>( 1, VolVar, "Cell_domain_mag" );
auto domain_cell_vec = std::make_shared<IO::Variable>( 3, VolVar, "Cell_domain_vec" );
point_node_mag->data.resize( N_points );
point_node_vec->data.resize( N_points, 3 );
for (int i=0; i<N_points; i++) {
point_node_mag->data(i) = distance(set1->points[i]);
point_node_vec->data(i,0) = set1->points[i].x;
point_node_vec->data(i,1) = set1->points[i].y;
point_node_vec->data(i,2) = set1->points[i].z;
for ( int i = 0; i < N_points; i++ ) {
point_node_mag->data( i ) = distance( set1->points[i] );
point_node_vec->data( i, 0 ) = set1->points[i].x;
point_node_vec->data( i, 1 ) = set1->points[i].y;
point_node_vec->data( i, 2 ) = set1->points[i].z;
}
set_node_mag->data = point_node_mag->data;
set_node_vec->data = point_node_vec->data;
list_node_mag->data.resize( 3*N_tri );
list_node_vec->data.resize( 3*N_tri, 3 );
for (int i=0; i<N_points; i++) {
list_node_mag->data(3*i+0) = distance(trilist->A[i]);
list_node_mag->data(3*i+1) = distance(trilist->B[i]);
list_node_mag->data(3*i+2) = distance(trilist->C[i]);
list_node_vec->data(3*i+0,0) = trilist->A[i].x;
list_node_vec->data(3*i+0,1) = trilist->A[i].y;
list_node_vec->data(3*i+0,2) = trilist->A[i].z;
list_node_vec->data(3*i+1,0) = trilist->B[i].x;
list_node_vec->data(3*i+1,1) = trilist->B[i].y;
list_node_vec->data(3*i+1,2) = trilist->B[i].z;
list_node_vec->data(3*i+2,0) = trilist->C[i].x;
list_node_vec->data(3*i+2,1) = trilist->C[i].y;
list_node_vec->data(3*i+2,2) = trilist->C[i].z;
list_node_mag->data.resize( 3 * N_tri );
list_node_vec->data.resize( 3 * N_tri, 3 );
for ( int i = 0; i < N_points; i++ ) {
list_node_mag->data( 3 * i + 0 ) = distance( trilist->A[i] );
list_node_mag->data( 3 * i + 1 ) = distance( trilist->B[i] );
list_node_mag->data( 3 * i + 2 ) = distance( trilist->C[i] );
list_node_vec->data( 3 * i + 0, 0 ) = trilist->A[i].x;
list_node_vec->data( 3 * i + 0, 1 ) = trilist->A[i].y;
list_node_vec->data( 3 * i + 0, 2 ) = trilist->A[i].z;
list_node_vec->data( 3 * i + 1, 0 ) = trilist->B[i].x;
list_node_vec->data( 3 * i + 1, 1 ) = trilist->B[i].y;
list_node_vec->data( 3 * i + 1, 2 ) = trilist->B[i].z;
list_node_vec->data( 3 * i + 2, 0 ) = trilist->C[i].x;
list_node_vec->data( 3 * i + 2, 1 ) = trilist->C[i].y;
list_node_vec->data( 3 * i + 2, 2 ) = trilist->C[i].z;
}
domain_node_mag->data.resize(domain->nx+1,domain->ny+1,domain->nz+1);
domain_node_vec->data.resize({(size_t)domain->nx+1,(size_t)domain->ny+1,(size_t)domain->nz+1,3});
for (int i=0; i<domain->nx+1; i++) {
for (int j=0; j<domain->ny+1; j++) {
for (int k=0; k<domain->nz+1; k++) {
domain_node_mag->data(i,j,k) = distance(Point(i,j,k));
domain_node_vec->data(i,j,k,0) = Point(i,j,k).x;
domain_node_vec->data(i,j,k,1) = Point(i,j,k).y;
domain_node_vec->data(i,j,k,2) = Point(i,j,k).z;
domain_node_mag->data.resize( domain->nx + 1, domain->ny + 1, domain->nz + 1 );
domain_node_vec->data.resize(
{ (size_t) domain->nx + 1, (size_t) domain->ny + 1, (size_t) domain->nz + 1, 3 } );
for ( int i = 0; i < domain->nx + 1; i++ ) {
for ( int j = 0; j < domain->ny + 1; j++ ) {
for ( int k = 0; k < domain->nz + 1; k++ ) {
domain_node_mag->data( i, j, k ) = distance( Point( i, j, k ) );
domain_node_vec->data( i, j, k, 0 ) = Point( i, j, k ).x;
domain_node_vec->data( i, j, k, 1 ) = Point( i, j, k ).y;
domain_node_vec->data( i, j, k, 2 ) = Point( i, j, k ).z;
}
}
}
set_cell_mag->data.resize( N_tri );
set_cell_vec->data.resize( N_tri, 3 );
for (int i=0; i<N_tri; i++) {
set_cell_mag->data(i) = i;
set_cell_vec->data(i,0) = 3*i+0;
set_cell_vec->data(i,1) = 3*i+1;
set_cell_vec->data(i,2) = 3*i+2;
for ( int i = 0; i < N_tri; i++ ) {
set_cell_mag->data( i ) = i;
set_cell_vec->data( i, 0 ) = 3 * i + 0;
set_cell_vec->data( i, 1 ) = 3 * i + 1;
set_cell_vec->data( i, 2 ) = 3 * i + 2;
}
list_cell_mag->data = set_cell_mag->data;
list_cell_vec->data = set_cell_vec->data;
domain_cell_mag->data.resize(domain->nx,domain->ny,domain->nz);
domain_cell_vec->data.resize({(size_t)domain->nx,(size_t)domain->ny,(size_t)domain->nz,3});
for (int i=0; i<domain->nx; i++) {
for (int j=0; j<domain->ny; j++) {
for (int k=0; k<domain->nz; k++) {
domain_cell_mag->data(i,j,k) = distance(Point(i,j,k));
domain_cell_vec->data(i,j,k,0) = Point(i,j,k).x;
domain_cell_vec->data(i,j,k,1) = Point(i,j,k).y;
domain_cell_vec->data(i,j,k,2) = Point(i,j,k).z;
domain_cell_mag->data.resize( domain->nx, domain->ny, domain->nz );
domain_cell_vec->data.resize(
{ (size_t) domain->nx, (size_t) domain->ny, (size_t) domain->nz, 3 } );
for ( int i = 0; i < domain->nx; i++ ) {
for ( int j = 0; j < domain->ny; j++ ) {
for ( int k = 0; k < domain->nz; k++ ) {
domain_cell_mag->data( i, j, k ) = distance( Point( i, j, k ) );
domain_cell_vec->data( i, j, k, 0 ) = Point( i, j, k ).x;
domain_cell_vec->data( i, j, k, 1 ) = Point( i, j, k ).y;
domain_cell_vec->data( i, j, k, 2 ) = Point( i, j, k ).z;
}
}
}
// Create the MeshDataStruct
std::vector<IO::MeshDataStruct> meshData(4);
std::vector<IO::MeshDataStruct> meshData( 4 );
meshData[0].meshName = "pointmesh";
meshData[0].mesh = set1;
meshData[0].vars.push_back(point_node_mag);
meshData[0].vars.push_back(point_node_vec);
meshData[0].mesh = set1;
meshData[0].vars.push_back( point_node_mag );
meshData[0].vars.push_back( point_node_vec );
meshData[1].meshName = "trimesh";
meshData[1].mesh = trimesh;
meshData[1].vars.push_back(set_node_mag);
meshData[1].vars.push_back(set_node_vec);
meshData[1].vars.push_back(set_cell_mag);
meshData[1].vars.push_back(set_cell_vec);
meshData[1].mesh = trimesh;
meshData[1].vars.push_back( set_node_mag );
meshData[1].vars.push_back( set_node_vec );
meshData[1].vars.push_back( set_cell_mag );
meshData[1].vars.push_back( set_cell_vec );
meshData[2].meshName = "trilist";
meshData[2].mesh = trilist;
meshData[2].vars.push_back(list_node_mag);
meshData[2].vars.push_back(list_node_vec);
meshData[2].vars.push_back(list_cell_mag);
meshData[2].vars.push_back(list_cell_vec);
meshData[2].mesh = trilist;
meshData[2].vars.push_back( list_node_mag );
meshData[2].vars.push_back( list_node_vec );
meshData[2].vars.push_back( list_cell_mag );
meshData[2].vars.push_back( list_cell_vec );
meshData[3].meshName = "domain";
meshData[3].mesh = domain;
meshData[3].vars.push_back(domain_node_mag);
meshData[3].vars.push_back(domain_node_vec);
meshData[3].vars.push_back(domain_cell_mag);
meshData[3].vars.push_back(domain_cell_vec);
meshData[3].mesh = domain;
meshData[3].vars.push_back( domain_node_mag );
meshData[3].vars.push_back( domain_node_vec );
meshData[3].vars.push_back( domain_cell_mag );
meshData[3].vars.push_back( domain_cell_vec );
for ( const auto &data : meshData )
ASSERT( data.check( true ) );
// Run the tests
testWriter( "old", meshData, ut );
@ -393,11 +401,9 @@ int main(int argc, char **argv)
// Finished
ut.report();
PROFILE_SAVE("TestWriter",true);
PROFILE_SAVE( "TestWriter", true );
int N_errors = ut.NumFailGlobal();
comm.barrier();
Utilities::shutdown();
return N_errors;
}