OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-10 02:25:33 -06:00
Code Issues Packages Projects Releases Wiki Activity

overhaul of GridDataOutput class to accept a vector type so that raw pointers are not overrun, this has lead to finding a good way to get the size of the data area that is alloacted by the DamarisVar wrapper class for damaris_alloc()

Browse Source
This commit is contained in:
Josh Bowden 2023-12-11 21:30:18 +01:00
parent 54d6db6f35
commit 9309f5a1bd
5 changed files with 491 additions and 129 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

59
ebos/damariswriter.hh
View File

@ -51,6 +51,8 @@
#include <stdexcept>
#include <string>
#include <omp.h>
#include <fmt/format.h>
#include <Damaris.h>
@ -136,8 +138,9 @@ class DamarisWriter : public EclGenericWriter<GetPropType<TypeTag, Properties::G
using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
using BaseType = EclGenericWriter<Grid,EquilGrid,GridView,ElementMapper,Scalar>;
typedef Opm::DamarisOutput::DamarisVar<int> DamarisVarInt ;
typedef Opm::DamarisOutput::DamarisVar<double> DamarisVarDbl ;
using DamarisVarInt = Opm::DamarisOutput::DamarisVar<int> ;
using DamarisVarChar = Opm::DamarisOutput::DamarisVar<char> ;
using DamarisVarDbl = Opm::DamarisOutput::DamarisVar<double> ;
public:
static void registerParameters()
@ -329,12 +332,11 @@ private:
{std::string("n_elements_local")},
std::string("MPI_RANK"), rank_) ) ;
// N.B. we have not set any offset values, so HDF5 collective and Dask arrays cannot be used.
mpi_rank_var->SetDamarisParameterAndShmem( {this->numElements_ } ) ;
int* shmem_mpi_ptr = mpi_rank_var->data_ptr() ;
mpi_rank_var->setDamarisParameterAndShmem( {this->numElements_ } ) ;
// Fill the created memory area
for (int i = 0 ; i < this->numElements_; i++ )
{
shmem_mpi_ptr[i] = rank_ ; // write the rank vaue to the shared memory area.
mpi_rank_var->data()[i] = rank_ ; // write the rank vaue to the shared memory area.
}
}
@ -412,7 +414,7 @@ private:
std::unique_ptr<DamarisVarInt> mpi_rank_var( new DamarisVarInt(1,
{std::string("n_elements_local")},
std::string("MPI_RANK"), rank_) ) ;
mpi_rank_var->SetDamarisPosition({*temp_int64_t}) ;
mpi_rank_var->setDamarisPosition({*temp_int64_t}) ;
}
@ -444,20 +446,31 @@ private:
// N.B. We have not set any position/offset values (using DamarisVar::SetDamarisPosition).
// They are not needed for mesh data as each process has a local geometric model.
// However, HDF5 collective and Dask arrays cannot be used for this data.
var_x->SetDamarisParameterAndShmem( { geomData.getNVertices() } ) ;
var_x->setDamarisParameterAndShmem( { geomData.getNVertices() } ) ;
std::unique_ptr<Opm::DamarisOutput::DamarisVar<double>> var_y(new Opm::DamarisOutput::DamarisVar<double>(1, {std::string("n_coords_local")}, std::string("coordset/coords/values/y"), rank_)) ;
var_y->ParameterIsSet() ;
var_y->SetPointersToDamarisShmem() ;
var_y->parameterIsSet() ;
var_y->setPointersToDamarisShmem() ;
std::unique_ptr<Opm::DamarisOutput::DamarisVar<double>> var_z(new Opm::DamarisOutput::DamarisVar<double>(1, {std::string("n_coords_local")}, std::string("coordset/coords/values/z"), rank_)) ;
var_z->ParameterIsSet() ;
var_z->SetPointersToDamarisShmem() ;
var_z->parameterIsSet() ;
var_z->setPointersToDamarisShmem() ;
// Now we can return the memory that Damaris has allocated in shmem
if ( geomData.writeGridPoints(var_x->data_ptr(),var_y->data_ptr(),var_z->data_ptr()) < 0)
// Now we can return the memory that Damaris has allocated in shmem and use it to write the X,y,z coordinates
double itime, ftime, exec_time;
itime = omp_get_wtime();
if ( geomData.writeGridPoints(*var_x,*var_y,*var_z) < 0)
DUNE_THROW(Dune::IOError, geomData.getError() );
//if ( geomData.writeGridPoints(var_x->data(),var_y->data(),var_z->data()) < 0)
// DUNE_THROW(Dune::IOError, geomData.getError() );
ftime = omp_get_wtime();
exec_time = ftime - itime;
// OpmLog::info("\n\nTime taken geomData.writeGridPoints(): is " + std::to_string(exec_time) ) ;
std::cout << "\n\n rank_: " << rank_ << " Time taken geomData.writeGridPoints(): is " + std::to_string(exec_time) << std::endl ;
// This is the template XML model for connectivity, offsets and types, as defined in initDamarisXmlFile.cpp which is used to
// build the internally generated Damaris XML configuration file.
// <parameter name="n_connectivity_ph" type="int" value="1" />
@ -471,27 +484,27 @@ private:
// <variable name="types" layout="n_types_layout_ph" type="scalar" visualizable="false" unit="" script="PythonConduitTest" time-varying="false" />
// </group>
std::unique_ptr<Opm::DamarisOutput::DamarisVar<int>> var_connectivity(new Opm::DamarisOutput::DamarisVar<int>(1, {std::string("n_connectivity_ph")}, std::string("topologies/topo/elements/connectivity"), rank_)) ;
var_connectivity->SetDamarisParameterAndShmem({ geomData.getNCorners()}) ;
std::unique_ptr<Opm::DamarisOutput::DamarisVar<int>> var_offsets(new Opm::DamarisOutput::DamarisVar<int>(1, {std::string("n_offsets_types_ph")}, std::string("topologies/topo/elements/offsets"), rank_)) ;
var_offsets->SetDamarisParameterAndShmem({ geomData.getNCells()}) ;
std::unique_ptr<Opm::DamarisOutput::DamarisVar<char>> var_types(new Opm::DamarisOutput::DamarisVar<char>(1, {std::string("n_offsets_types_ph")}, std::string("topologies/topo/elements/types"), rank_)) ;
var_types->ParameterIsSet() ;
var_types->SetPointersToDamarisShmem() ;
std::unique_ptr<DamarisVarInt> var_connectivity(new DamarisVarInt(1, {std::string("n_connectivity_ph")}, std::string("topologies/topo/elements/connectivity"), rank_)) ;
var_connectivity->setDamarisParameterAndShmem({ geomData.getNCorners()}) ;
std::unique_ptr<DamarisVarInt> var_offsets(new DamarisVarInt(1, {std::string("n_offsets_types_ph")}, std::string("topologies/topo/elements/offsets"), rank_)) ;
var_offsets->setDamarisParameterAndShmem({ geomData.getNCells()}) ;
std::unique_ptr<DamarisVarChar> var_types(new DamarisVarChar(1, {std::string("n_offsets_types_ph")}, std::string("topologies/topo/elements/types"), rank_)) ;
var_types->parameterIsSet() ;
var_types->setPointersToDamarisShmem() ;
// Copy the mesh data from the Durne grid
long i = 0 ;
Opm::GridDataOutput::ConnectivityVertexOrder vtkorder = Opm::GridDataOutput::VTK ;
i = geomData.writeConnectivity(var_connectivity->data_ptr(), vtkorder) ;
i = geomData.writeConnectivity(*var_connectivity, vtkorder) ;
if ( i != geomData.getNCorners())
DUNE_THROW(Dune::IOError, geomData.getError());
i = geomData.writeOffsetsCells(var_offsets->data_ptr());
i = geomData.writeOffsetsCells(*var_offsets);
if ( i != geomData.getNCells()+1)
DUNE_THROW(Dune::IOError,geomData.getError());
i = geomData.writeCellTypes(var_types->data_ptr()) ;
i = geomData.writeCellTypes(*var_types) ;
if ( i != geomData.getNCells())
DUNE_THROW(Dune::IOError,geomData.getError());
}

2
opm/simulators/utils/DamarisKeywords.cpp
View File

@ -140,7 +140,7 @@ DamarisSettings::getKeywords([[maybe_unused]] const Parallel::Communication& com
"(--damaris-python-paraview-script command line argument) scripts are valid, however only one "
"type of analysis is supported in a single simulation (due to Paraview installing mpi4py library "
"locally and without header files). "
"Please choose one or the other method of analysis for now. Exiting." )
"Please choose one or the other method of analysis for now. Exiting." );
}
std::string damarisOutputCollective_str;

207
opm/simulators/utils/DamarisVar.hpp
View File

@ -17,8 +17,8 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef DAMARIS_MODEL_DATA_HPP
#define DAMARIS_MODEL_DATA_HPP
#ifndef DAMARISVAR_HPP
#define DAMARISVAR_HPP
#include <string>
#include <iostream>
@ -26,15 +26,19 @@
#include <cassert>
#include <typeinfo>
#include <Damaris.h>
#define BOOST_BIND_GLOBAL_PLACEHOLDERS 1
#include <Damaris.h>
#include <damaris/data/VariableManager.hpp>
#include <damaris/data/Variable.hpp>
#include <damaris/data/Block.hpp>
/*
File: DamarisVar.hpp
Author: Joshua Bowden, Inria
Date: 06/02/2023
The DamarisVar class can be used to allocate memory in the Damaris shared memory region and a user can supply
the data required for the variable to pass to Damaris.
the data required for the variable. This data will then be directly available to the Damaris server side (server core) plugin code.
*/
@ -101,24 +105,17 @@ namespace Opm
class DamarisVarBase {
public:
// DamarisVarBase(int dims, std::vector<std::string>& param_names, std::string& variable_name, int rank=0) = 0;
virtual ~DamarisVarBase( void ) {};
virtual void PrintError ( void ) = 0;
virtual bool HasError( void ) = 0;
// virtual void SetDamarisParameterAndShmem( std::vector<int>& paramSizeVal ) = 0;
virtual void SetDamarisParameterAndShmem( std::vector<int> paramSizeVal ) = 0;
virtual void SetDamarisParameter( std::vector<int>& paramSizeVal ) = 0;
virtual void SetDamarisPosition( std::vector<int64_t> positionsVals ) = 0;
virtual void SetPointersToDamarisShmem( void ) = 0;
virtual void CommitVariableDamarisShmem( void ) = 0;
virtual void ClearVariableDamarisShmem( void ) = 0;
// virtual void * data_ptr( void ) = 0;
virtual void printError ( void ) = 0;
virtual bool hasError( void ) = 0;
virtual void setDamarisParameterAndShmem( std::vector<int> paramSizeVal ) = 0;
virtual void setDamarisParameter( std::vector<int>& paramSizeVal ) = 0;
virtual void setDamarisPosition( std::vector<int64_t> positionsVals ) = 0;
virtual void setPointersToDamarisShmem( void ) = 0;
virtual void commitVariableDamarisShmem( void ) = 0;
virtual void clearVariableDamarisShmem( void ) = 0;
virtual std::string & variable_name( void ) = 0;
}; // class DamarisVar
}; // class DamarisVarBase
/**
* class to store a Damaris variable representation for the XML file (can be used with /ref class DamarisKeywords).
@ -159,6 +156,7 @@ namespace Opm
std::ostringstream dam_err_sstr_; //!< Use dam_err_sstr.str() to return an error string describing detected error
DamarisVarXMLAttributes xml_attributes_; //!< The extra elements that need to be part of a Damaris <variable> type. They are simple string values that may reference other XML elements (and could be checked for existence etc.)
T * data_ptr_; //!< This pointer will be mapped to the Damaris shared memory area for the variable in the SetPointersToDamarisShmem() method. The type T will match the Layout type
size_t current_size_ ; //!< The total number of elements that may be held by this part of the variable - returned by the size() method.
public:
@ -213,7 +211,8 @@ namespace Opm
if ( !TestType(variable_name) ) {
std::exit(-1);
}
current_size_ = 0 ;
num_params_ = param_names_.size();
param_sizes_ = new int(num_params_);
positions_ = new int64_t(dims);
@ -269,7 +268,7 @@ namespace Opm
paramaters_set_ = false;
has_error_ = false;
SetDamarisParameterAndShmem( param_values ); // Initialise the memory size in the constructor.
setDamarisParameterAndShmem( param_values ); // Initialise the memory size in the constructor.
}
~DamarisVar( void ) {
@ -277,11 +276,11 @@ namespace Opm
delete [] positions_;
if (data_ptr_ != nullptr)
{
CommitVariableDamarisShmem();
ClearVariableDamarisShmem();
commitVariableDamarisShmem();
clearVariableDamarisShmem();
}
if (this->HasError())
PrintError(); // flush out any error messages
if (this->hasError())
printError(); // flush out any error messages
}
/**
@ -307,7 +306,7 @@ namespace Opm
double td = 0.0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -316,7 +315,7 @@ namespace Opm
float td = 0.0f;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -325,7 +324,7 @@ namespace Opm
char td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -334,7 +333,7 @@ namespace Opm
unsigned char td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -343,7 +342,7 @@ namespace Opm
short td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -352,7 +351,7 @@ namespace Opm
unsigned short td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -361,7 +360,7 @@ namespace Opm
int td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -370,7 +369,7 @@ namespace Opm
unsigned int td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -379,7 +378,7 @@ namespace Opm
long td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -388,7 +387,7 @@ namespace Opm
unsigned long td = 0;
const std::type_info& t2 = typeid(td);
if (t1 != t2) {
OutputErrorAndAssert(variable_name, t1.name(), t2.name());
outputErrorAndAssert(variable_name, t1.name(), t2.name());
resbool = false;
}
}
@ -408,21 +407,60 @@ namespace Opm
/**
* Allow a user to indicate that the Damaris variable has allocated a size -
* This method is usefull as a single paramater can control one or more layouts and
* This method is usefull as a single parameter can control one or more layouts and
* a single layout can describe the size of multiple <variable> elements.
* i.e. The current variable has had it's paramater(s) set through via another variable.
* i.e. Use when the current variable has had it's paramater(s) set through via another variable.
*/
void ParameterIsSet()
void parameterIsSet()
{
paramaters_set_ = true;
getDataStoreBlockSize() ;
}
long getDataStoreBlockSize( void )
{
long total_size = 1 ;
if (paramaters_set_ == true)
{
std::shared_ptr<damaris::Variable> v = damaris::VariableManager::Search(this->variable_name_.c_str());
damaris::BlocksByIteration::iterator begin;
damaris::BlocksByIteration::iterator end;
int iteration = 0 ;
v->GetBlocksByIteration(iteration, begin, end);
for (damaris::BlocksByIteration::iterator bid = begin; bid != end; bid++) {
std::shared_ptr<damaris::Block> b = *bid;
void PrintError ( void )
// possibly multi-dimensional, so we need a value for each dimension
int blockDimension = b->GetDimensions();
// Compute the size in each dimension
for (int i = 0 ; i < blockDimension ; i++)
{
total_size *= ( b->GetEndIndex(i) - b->GetStartIndex(i) + 1 ); // This includes a variables Ghost zones.
}
}
}
if (total_size > 0) {
current_size_ = total_size ;
} else {
dam_err_sstr_ << " ERROR rank =" << rank_ << " : class DamarisVar::getDataStoreBlockSize() " <<
"The total size of the variable is 0 - please check input paramSizeVal array." << std::endl;
has_error_ = true;
}
return total_size ;
}
void printError ( void )
{
std::cerr << dam_err_sstr_.str();
}
bool HasError( void )
bool hasError( void )
{
return (has_error_);
}
@ -430,9 +468,14 @@ namespace Opm
/**
* Returns the data pointer to shared memory, or nullptr if it has not been allocated
*/
T * data_ptr( void )
T * data( void )
{
return (data_ptr_);
if (paramaters_set_ == true )
{
return (data_ptr_); // This still could be nullptr
} else {
return (nullptr);
}
}
std::string & variable_name( void )
@ -443,7 +486,7 @@ namespace Opm
/**
* Creates the XML representation of the variable from the available strings
*/
std::string ReturnXMLForVariable ( void )
std::string returnXMLForVariable ( void )
{
std::ostringstream var_sstr;
@ -461,28 +504,45 @@ namespace Opm
*
*
*/
void SetDamarisParameterAndShmem( std::vector<int> paramSizeVal )
void setDamarisParameterAndShmem( std::vector<int> paramSizeVal )
{
this->SetDamarisParameter( paramSizeVal );
this->SetPointersToDamarisShmem();
this->setDamarisParameter( paramSizeVal );
this->setPointersToDamarisShmem();
}
/**
* Method to set the Damaris paramater values.
* Returns the number of elements in the memory area.
* Used as a method for compatibility with std::vector
*/
size_t size()
{
if (paramaters_set_ == true )
{
return current_size_ ;
} else {
return 0 ;
}
}
/**
* Method to set the Damaris paramater values. Also calculates the total number
* of elements in the variable (current_size_) that is returned bt size() method.
*
* /param [IN] paramSizeVal : An pointer to a value or array of values to set. One element per param_names_ string
*
* /implicit : Implicitly uses the array of paramater names: \ref param_names_
*/
void SetDamarisParameter( std::vector<int>& paramSizeVal )
void setDamarisParameter( std::vector<int>& paramSizeVal )
{
assert(paramSizeVal.size() == num_params_);
bool resbool = true;
// size_t total_size = 1 ;
for (int varnum = 0; varnum < num_params_; varnum++)
{
param_sizes_[varnum] = paramSizeVal[varnum];
// total_size *= param_sizes_[varnum] ;
dam_err_ = damaris_parameter_set(param_names_[varnum].c_str(), &paramSizeVal[varnum], sizeof(int));
if (dam_err_ != DAMARIS_OK) {
dam_err_sstr_ << " ERROR rank =" << rank_ << " : class DamarisVar : damaris_parameter_set(\"" << param_names_[varnum]
@ -491,8 +551,15 @@ namespace Opm
has_error_ = true;
}
}
if (resbool == true)
paramaters_set_ = true;
{
parameterIsSet() ; // sets paramaters_set_ and gets the size of the variables block storage (as number of elemnts)
}
if (hasError()) {
printError()
}
}
/**
@ -503,7 +570,7 @@ namespace Opm
*
* /implicit : Implicitly uses the variable name: \ref variable_name_
*/
void SetDamarisPosition( std::vector<int64_t> positionsVals )
void setDamarisPosition( std::vector<int64_t> positionsVals )
{
assert(positionsVals.size() == dims_);
@ -517,15 +584,19 @@ namespace Opm
<< variable_name_ << "\", positionsVals); Damaris error = " << damaris_error_string(dam_err_) << std::endl;
has_error_ = true;
}
if (hasError()) {
printError()
}
}
/**
* Method to set the Damaris paramater values.
* Method to set the internal pointer (data_ptr_) to the Damaris shared memory area.
*
* /implicit : Implicitly uses the Damaris variable name string \ref variable_name_
* /implicit : Implicitly uses the class data element : \ref data_ptr_
*/
void SetPointersToDamarisShmem( void )
void setPointersToDamarisShmem( void )
{
if (paramaters_set_ == true )
{
@ -538,12 +609,22 @@ namespace Opm
}
} else {
dam_err_ = -1;
dam_err_sstr_ << "ERROR rank =" << rank_ << " : class DamarisVar : SetDamarisParameter() should be called first so as to define the size of the memory block required for variable : " << variable_name_ << std::endl;
dam_err_sstr_ << "ERROR rank =" << rank_ << " : class DamarisVar : setDamarisParameter() should be called first so as to define the size of the memory block required for variable : " << variable_name_ << std::endl;
has_error_ = true;
}
if (hasError()) {
printError()
}
}
void SetPointersToDamarisShmem( T ** ptr_in )
/**
* Method to set the an externaly sourced pointer to point to the Damaris shared memory.
*
* /implicit : Implicitly uses the Damaris variable name string \ref variable_name_
* /implicit : Implicitly uses the class data element : \ref data_ptr_
*/
void setPointersToDamarisShmem( T ** ptr_in )
{
if (paramaters_set_ == true )
{
@ -560,9 +641,13 @@ namespace Opm
data_ptr_ = temp_ptr;
} else {
dam_err_ = -1;
dam_err_sstr_ << " ERROR rank =" << rank_ << " : class DamarisVar : SetDamarisParameter() should be called first so as to define the size of the memory block required for variable : " << variable_name_ << std::endl;
dam_err_sstr_ << " ERROR rank =" << rank_ << " : class DamarisVar : setDamarisParameter() should be called first so as to define the size of the memory block required for variable : " << variable_name_ << std::endl;
has_error_ = true;
}
if (hasError()) {
printError()
}
}
/**
@ -571,7 +656,7 @@ namespace Opm
*
* /implicit : Implicitly uses the variable name string \ref variable_name_
*/
void CommitVariableDamarisShmem( void )
void commitVariableDamarisShmem( void )
{
// Signal to Damaris we are done writing data for this iteration
dam_err_ = damaris_commit (variable_name_.c_str());
@ -588,7 +673,7 @@ namespace Opm
*
* /implicit : Implicitly uses the variable name string \ref variable_name_
*/
void ClearVariableDamarisShmem( void )
void clearVariableDamarisShmem( void )
{
// Signal to Damaris it has complete charge of the memory area
dam_err_ = damaris_clear(variable_name_.c_str());
@ -600,11 +685,11 @@ namespace Opm
data_ptr_ = nullptr;
}
private:
void OutputErrorAndAssert(std::string& var_name, std::string type_name1, std::string type_name2)
void outputErrorAndAssert(std::string& var_name, std::string type_name1, std::string type_name2)
{
dam_err_sstr_ << "ERROR rank =" << rank_ << " : DamarisVar::DamarisVar () variable_name_: \"" << var_name
<< "\" The template type of Type of DamarisVar<T> in the code: " << type_name1 << " does not match type in XML (float)" << std::endl;
PrintError();
printError();
assert( type_name1 == type_name2 );
}
}; // class DamarisVar

338
opm/simulators/utils/GridDataOutput.hpp
View File

@ -48,6 +48,7 @@
geomData.printGridDetails();
int nvert = geomData.getNVertices();
// example using seperate x, y and z arrays
int nvert = geomData.getNVertices();
double * x_vert = new double[nvert];
@ -163,14 +164,16 @@ namespace Opm::GridDataOutput
}
}
/**
Write the positions of vertices - directly to the pointers given in parameters
Returns the number of vertices written
*/
template <typename T>
long writeGridPoints( T* x_inout, T* y_inout, T* z_inout )
long writeGridPoints( T* x_inout, T* y_inout, T* z_inout, long max_size=0 )
{
if (max_size < nvertices_) {
assert(max_size >= nvertices_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints( T& x_inout, T& y_inout, T& z_inout ) "
+ " Input objects max_size (" + std::to_string(max_size) + ") is not sufficient to fit the nvertices_ values ("
+ std::to_string(nvertices_) + ")" );
}
long i = 0;
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all) )
@ -194,15 +197,70 @@ namespace Opm::GridDataOutput
assert(i == nvertices_); // As we are templated on the Dune::PartitionSet<partitions>, this cannot change
return i;
}
/**
Write the positions of vertices - directly to the pointers given in parameters
Returns the number of vertices written
*/
template <typename T>
long writeGridPoints( T& x_inout, T& y_inout, T& z_inout )
{
size_t check_size = x_inout.size() ;
using VT = decltype(x_inout.data()[0]);
if (check_size < nvertices_) {
// assert(check_size >= nvertices_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints( T& x_inout, T& y_inout, T& z_inout ) "
+ " Input objects size " + std::to_string(check_size) + " is not sufficient to fit the nvertices_ values( "
+ std::to_string(nvertices_) + " )" );
}
long i = 0;
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all) )
{
auto xyz_local = vit.geometry().corner(0); // vertices only have one corner
x_inout.data()[i] = static_cast<VT>(xyz_local[0]);
y_inout.data()[i] = static_cast<VT>(xyz_local[1]);
z_inout.data()[i] = static_cast<VT>(xyz_local[2]);
i++;
}
} else if (dimw_ == 2) {
double td = 0.0 ;
for (const auto& vit : vertices(gridView_, Dune::Partitions::all) )
{
auto xyz_local = vit.geometry().corner(0); // vertices only have one corner
x_inout.data()[i] = static_cast<VT>(xyz_local[0]);
y_inout.data()[i] = static_cast<VT>(xyz_local[1]);
z_inout.data()[i] = static_cast<VT>(td);
i++;
}
}
assert(i == nvertices_); // As we are templated on the Dune::PartitionSet<partitions>, this cannot change
return i;
}
/**
Write positions of vertices as array of structures : x,y,z,x,y,z,x,y,z,...
Returns the number of vertices written
*/
template <typename T>
long writeGridPoints_AOS( T* xyz_inout )
long writeGridPoints_AOS( T* xyz_inout, long max_size=0 )
{
if (max_size < nvertices_ * 3) {
assert(max_size >= nvertices_ * 3);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints_AOS( T* xyz_inout ) "
+ " Input objects max_size (" + std::to_string(max_size) + ") is not sufficient to fit the nvertices_ * 3 values ("
+ std::to_string(nvertices_ * 3) + ")" );
}
long i = 0;
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
@ -224,19 +282,69 @@ namespace Opm::GridDataOutput
}
return ( (i) / 3 );
}
/**
Write positions of vertices as array of structures : x,y,z,x,y,z,x,y,z,...
Returns the number of vertices written
*/
template <typename VectType>
long writeGridPoints_AOS( VectType& xyz_inout )
{
size_t check_size = xyz_inout.size() ;
using VT = decltype(xyz_inout.data()[0]);
if (check_size < nvertices_ * 3) {
assert(check_size >= nvertices_ * 3);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints_AOS( VectType& xyz_inout ) "
+ " Input objects check_size (" + std::to_string(check_size) + ") is not sufficient to fit the nvertices_ * 3 values ("
+ std::to_string(nvertices_ * 3) + ")" );
}
long i = 0;
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
{
auto xyz_local = vit.geometry().corner(0);
xyz_inout.data()[i++] = static_cast<VT>(xyz_local[0]);
xyz_inout[i++] = static_cast<VT>(xyz_local[1]);
xyz_inout[i++] = static_cast<VT>(xyz_local[2]);
}
} else if (dimw_ == 2) {
double td = 0.0 ;
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
{
auto xyz_local = vit.geometry().corner(0);
xyz_inout[i++] = static_cast<VT>(xyz_local[0]);
xyz_inout[i++] = static_cast<VT>(xyz_local[1]);
xyz_inout[i++] = static_cast<VT>(td);
}
}
return ( (i) / 3 );
}
/**
Write positions of vertices as structure of arrays : x,x,x,...,y,y,y,...,z,z,z,...
Returns the number of vertices written
*/
template <typename T>
long writeGridPoints_SOA( T* xyz_inout )
long writeGridPoints_SOA( T* xyz_inout, long max_size=0 )
{
if (max_size < nvertices_ * 3) {
assert(max_size >= nvertices_ * 3);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints_SOA( T& xyz_inout ) "
+ " Input objects max_size (" + std::to_string(max_size) + ") is not sufficient to fit the nvertices_ * 3 values ("
+ std::to_string(nvertices_ * 3) + ")" );
}
long i = 0;
// Get offsets into structure
T* xyz_inout_y = xyz_inout + nvertices_;
T* xyz_inout_z = xyz_inout + (2*nvertices_);
T* xyz_inout_z = xyz_inout + (2 * nvertices_);
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
@ -259,6 +367,53 @@ namespace Opm::GridDataOutput
}
return (i);
}
/**
Write positions of vertices as structure of arrays : x,x,x,...,y,y,y,...,z,z,z,...
Returns the number of vertices written
*/
template <typename VectType>
long writeGridPoints_SOA( VectType& xyz_inout )
{
size_t check_size = xyz_inout.size() ;
if (check_size < nvertices_ * 3) {
assert(check_size >= nvertices_ * 3);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeGridPoints_SOA( VectType& xyz_inout ) "
+ " Input objects check_size (" + std::to_string(check_size) + ") is not sufficient to fit the nvertices_ * 3 values ("
+ std::to_string(nvertices_ * 3) + ")" );
}
using VT = decltype(xyz_inout.data()[0]);
long i = 0;
// Get offsets into structure
VT* xyz_inout_y = xyz_inout.data() + nvertices_;
VT* xyz_inout_z = xyz_inout.data() + (2 * nvertices_);
if (dimw_ == 3) {
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
{
auto xyz_local = vit.geometry().corner(0);
xyz_inout.data()[i] = static_cast<VT>(xyz_local[0]);
xyz_inout_y[i]= static_cast<VT>(xyz_local[1]);
xyz_inout_z[i] = static_cast<VT>(xyz_local[2]);
i++;
}
} else if (dimw_ == 2) {
double td = 0.0 ;
for (const auto& vit : vertices(gridView_, Dune::Partitions::all))
{
auto xyz_local = vit.geometry().corner(0);
xyz_inout.data()[i] = static_cast<VT>(xyz_local[0]);
xyz_inout_y[i]= static_cast<VT>(xyz_local[1]);
xyz_inout_z[i] = static_cast<VT>(td);
i++;
}
}
return (i);
}
/**
* Write the connectivity array - directly to the pointer given in parameter 1
@ -266,9 +421,16 @@ namespace Opm::GridDataOutput
Returns the number of corner indices written.
*/
template <typename I>
long writeConnectivity(I * connectivity_inout, ConnectivityVertexOrder whichOrder)
template <typename Integer>
long writeConnectivity(Integer * connectivity_inout, ConnectivityVertexOrder whichOrder, long max_size=0 )
{
if (max_size < ncorners_ ) {
assert(max_size >= ncorners_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeConnectivity( T* connectivity_inout ) "
+ " Input objects size (" + std::to_string(max_size) + ") is not sufficient to fit the ncorners_ values ("
+ std::to_string(ncorners_) + ")" );
}
long i = 0;
if ( whichOrder == DUNE ) {
// DUNE order
@ -299,15 +461,72 @@ namespace Opm::GridDataOutput
}
return (i);
}
/**
* Write the connectivity array - directly to the pointer given in parameter 1
Reorders the indecies as selected either in DUNE order or VTK order.
Returns the number of corner indices written.
*/
template <typename VectType>
long writeConnectivity(VectType & connectivity_inout, ConnectivityVertexOrder whichOrder)
{
size_t check_size = connectivity_inout.size() ;
if (check_size < ncorners_ ) {
assert(check_size >= ncorners_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeConnectivity( VectType& connectivity_inout ) "
+ " Input objects size (" + std::to_string(check_size) + ") is not sufficient to fit the ncorners_ values ("
+ std::to_string(ncorners_) + ")" );
}
using VT = decltype(connectivity_inout.data()[0]);
long i = 0;
if ( whichOrder == DUNE ) {
// DUNE order
for (const auto& cit : elements(gridView_, dunePartition_))
{
auto cell_corners = cit.geometry().corners();
for( auto vx = 0; vx < cell_corners; ++ vx )
{
const int vxIdx = gridView_.indexSet().subIndex( cit, vx, 3 );
connectivity_inout.data()[i + vx] = vxIdx;
}
i += cell_corners;
}
} else {
// VTK order
for (const auto& cit : elements(gridView_, dunePartition_))
{
auto cell_corners = cit.geometry().corners();
for( auto vx = 0; vx < cell_corners; ++ vx )
{
const int vxIdx = gridView_.indexSet().subIndex( cit, vx, 3 );
int vtkOrder;
vtkOrder = Dune::VTK::renumber(cit.type(), vx);
connectivity_inout.data()[i + vtkOrder] = vxIdx;
}
i += cell_corners;
}
}
return (i);
}
/**
* Write the offsets values - directly to the pointer given in parameter 1
Returns the number of offset values written, which should be 1 greater than ncells_
or -1 if an error was detected
Returns the (number of offset values written + 1)
*/
template <typename I>
long writeOffsetsCells( I* offsets_inout )
template <typename Integer>
long writeOffsetsCells( Integer* offsets_inout, long max_size=0 )
{
if (max_size < ncells_ ) {
assert(max_size >= ncells_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeOffsetsCells( T* offsets_inout ) "
+ " Input objects max_size (" + std::to_string(max_size) + ") is not sufficient to fit the ncells_ values ("
+ std::to_string(ncells_) + ")" );
}
long i = 1;
offsets_inout[0] = 0;
for (const auto& cit : elements(gridView_, dunePartition_))
@ -318,21 +537,80 @@ namespace Opm::GridDataOutput
}
return (i); // This should be 1 greater than ncells_
}
/**
* Write the offsets values - directly to the pointer given in parameter 1
Returns the (number of offset values written + 1)
*/
template <typename VectType>
long writeOffsetsCells( VectType& offsets_inout )
{
size_t check_size = offsets_inout.size() ;
if (check_size < ncells_ ) {
assert(check_size >= ncells_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeOffsetsCells( VectType& offsets_inout ) "
+ " Input objects check_size (" + std::to_string(check_size) + ") is not sufficient to fit the ncells_ values ("
+ std::to_string(ncells_) + ")" );
}
// using VT = decltype(offsets_inout.data()[0]);
long i = 1;
offsets_inout.data()[0] = 0;
for (const auto& cit : elements(gridView_, dunePartition_))
{
auto cell_corners = cit.geometry().corners();
offsets_inout.data()[i] = offsets_inout.data()[i-1] + cell_corners;
i++;
}
return (i); // This should be 1 greater than ncells_
}
/**
* Write the Cell types array - directly to the pointer given in parameter 1
*/
template <typename I>
long writeCellTypes( I* types_inout)
template <typename Integer>
long writeCellTypes( Integer* types_inout, long max_size=0)
{
if (max_size < ncells_ ) {
assert(max_size >= ncells_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeCellTypes( T* types_inout ) "
+ " Input objects max_size (" + std::to_string(max_size) + ") is not sufficient to fit the ncells_ values ("
+ std::to_string(ncells_) + ")" );
}
int i = 0;
for (const auto& cit : elements(gridView_, dunePartition_))
{
I vtktype = static_cast<I>(Dune::VTK::geometryType(cit.type()));
Integer vtktype = static_cast<Integer>(Dune::VTK::geometryType(cit.type()));
types_inout[i++] = vtktype;
}
return (i);
}
/**
* Write the Cell types array - directly to the pointer given in parameter 1
*/
template <typename VectType>
long writeCellTypes( VectType& types_inout)
{
size_t check_size = types_inout.size() ;
if (check_size < ncells_ ) {
assert(check_size >= ncells_);
OPM_THROW(std::runtime_error, "Opm::GridDataOutput::writeCellTypes( VectType& types_inout ) "
+ " Input objects check_size (" + std::to_string(check_size) + ") is not sufficient to fit the ncells_ values ("
+ std::to_string(ncells_) + ")" );
}
using VT = decltype(types_inout.data()[0]);
int i = 0;
for (const auto& cit : elements(gridView_, dunePartition_))
{
int vtktype = static_cast<int>(Dune::VTK::geometryType(cit.type()));
types_inout.data()[i++] = vtktype;
}
return (i);
}
std::string getPartitionTypeString ( )
{
@ -361,29 +639,15 @@ namespace Opm::GridDataOutput
return ( this->dunePartition_ );
}
void printGridDetails()
void printGridDetails(std::ostream& outstr )
{
std::cout << "Dune Partition = " << partition_value_ << ", " << getPartitionTypeString() << std::endl;
printNCells();
printNVertices();
printNCorners();
outstr << "Dune Partition = " << partition_value_ << ", " << getPartitionTypeString() << std::endl;
outstr << "ncells_: " << getNCells() << std::endl;
outstr << "nvertices_: " << getNVertices() << std::endl;
outstr << "ncorners_: " << getNCorners() << std::endl;
}
void printNCells()
{
std::cout << "ncells = " << ncells_ << std::endl;
}
void printNVertices()
{
std::cout << "nvertices = " << nvertices_ << std::endl;
}
void printNCorners()
{
std::cout << "ncorners = " << ncorners_ << std::endl;
}
int getNCells()
{
return(ncells_);

14
opm/simulators/utils/initDamarisXmlFile.cpp
View File

@ -108,20 +108,20 @@ std::string initDamarisXmlFile()
<parameter name="n_coords_local" type="int" value="1" />
<layout name="n_coords_layout" type="double" dimensions="n_coords_local" comment="For the individual x, y and z coordinates of the mesh vertices, these values are referenced in the topologies/topo/subelements/connectivity_pg data" />
<group name="coordset/coords/values">
<variable name="x" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="y" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="z" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="x" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
<variable name="y" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
<variable name="z" layout="n_coords_layout" type="scalar" visualizable="false" unit="m" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
</group>
<parameter name="n_connectivity_ph" type="int" value="1" />
<layout name="n_connections_layout_ph" type="int" dimensions="n_connectivity_ph" comment="Layout for connectivities " />
<parameter name="n_offsets_types_ph" type="int" value="1" />
<layout name="n_offsets_layout_ph" type="int" dimensions="n_offsets_types_ph + 1" comment="Layout for the offsets_ph" />
<layout name="n_offsets_layout_ph" type="int" dimensions="n_offsets_types_ph" comment="Layout for the offsets_ph" />
<layout name="n_types_layout_ph" type="char" dimensions="n_offsets_types_ph" comment="Layout for the types_ph " />
<group name="topologies/topo/elements">
<variable name="connectivity" layout="n_connections_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="offsets" layout="n_offsets_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="types" layout="n_types_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" />
<variable name="connectivity" layout="n_connections_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
<variable name="offsets" layout="n_offsets_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
<variable name="types" layout="n_types_layout_ph" type="scalar" visualizable="false" script="_MAKE_AVAILABLE_IN_PYTHON_" time-varying="false" store="MyStore" />
</group>
<mesh name="us_mesh" type="unstructured" topology="3" time-varying="false"