Dear deal.II Community,
Within the context of an MPI parallel code, I need to broadcast the
ConstitutiveParameters struct which contains std::unique pointer to the
base class StrainEnergyFuncParameters; cf. MWE file.
Using the ideas discussed in the post
https://groups.google.com/g/dealii/c/bmqTbs4ofHQ/m/aCuDpUekBAAJ, the code
works fine if the members are concrete derived classes since there the
serialize function is enough. However, I get segmentation fault during the
unpack() call since the deserialize() function is not called.
Running the attached MWE would explain the issue.
Could someone please help in understanding how I could solve this issue?
Thanks and best regards,
Paras Kumar
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/5f719546-d28f-4670-981b-00e04283c546n%40googlegroups.com.
#include <deal.II/base/utilities.h>
#include <deal.II/base/mpi.h>
#include <deal.II/base/conditional_ostream.h>
using UnsignedIntType = unsigned int;
using FloatType= double;
FloatType Epsilon_Tol = 1e-8;
template <typename NumberType = double>
bool
Check_ifEqual(const NumberType &a,
const NumberType &b,
const NumberType &tolFac = 1e-5,
const NumberType &baseTol = Epsilon_Tol)
{
const NumberType tol =
(0.5 * (std::fabs(a) + std::fabs(b)) * tolFac) + baseTol;
bool result = false;
if (std::fabs(a - b) < tol)
result = true;
return (result);
}
enum StrainEnergyFuncType
{
NH = 1,
Og = 2,
};
template <typename NumberType = double>
struct StrainEnergyFuncParameters
{
StrainEnergyFuncType strainEnergyFuncType_;
StrainEnergyFuncParameters(){};
~StrainEnergyFuncParameters() = default;
virtual void
Print_details(std::ostream &outStream) const = 0;
template <class Archive>
void
serialize(Archive &ar, const unsigned int version)
{
ar &strainEnergyFuncType_;
}
};
template <typename NumberType = double>
struct NeoHookeParameters : public StrainEnergyFuncParameters<NumberType>
{
NumberType shearModulus_;
NumberType poissonRatio_;
NeoHookeParameters();
NeoHookeParameters(const StrainEnergyFuncType strainEnergyFuncType);
NeoHookeParameters(const NeoHookeParameters<NumberType> &other);
NeoHookeParameters<NumberType>& operator=(const NeoHookeParameters<NumberType> &other);
bool operator==(const NeoHookeParameters<NumberType> &other) const;
virtual void Print_details(std::ostream &outStream) const override;
template <class Archive>
void
serialize(Archive &ar, const unsigned int version)
{
std::cout << "my ranK: " << dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
<< "Entered NeooHookeParameters::serialize()"
<< std::endl;
ar &boost::serialization::base_object<StrainEnergyFuncParameters<NumberType>>(*this);
ar &shearModulus_;
ar &poissonRatio_;
}
template <class Archive>
void
deserialize(Archive &ar, const unsigned int version)
{
std::cout << "my ranK: " << dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
<< "Entered NeooHookeParameters::deserialize()"
<< std::endl;
ar &boost::serialization::base_object<StrainEnergyFuncParameters<NumberType>>(*this);
ar &shearModulus_;
ar &poissonRatio_;
}
};
template <typename NumberType>
NeoHookeParameters<NumberType>::NeoHookeParameters()
{}
template <typename NumberType>
NeoHookeParameters<NumberType>::NeoHookeParameters(const StrainEnergyFuncType strainEnergyFuncType)
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = strainEnergyFuncType;
}
template <typename NumberType>
NeoHookeParameters<NumberType>::NeoHookeParameters(const NeoHookeParameters<NumberType> &other)
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = other.strainEnergyFuncType_;
this->shearModulus_ = other.shearModulus_;
this->poissonRatio_ = other.poissonRatio_;
}
template <typename NumberType>
NeoHookeParameters<NumberType> &
NeoHookeParameters<NumberType>::operator=(const NeoHookeParameters<NumberType> &other)
{
if (this != &other) // confirm there is no self-assignment
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = other.strainEnergyFuncType_;
this->shearModulus_ = other.shearModulus_;
this->poissonRatio_ = other.poissonRatio_;
}
return (*this);
}
template <typename NumberType>
bool NeoHookeParameters<NumberType>::operator==(const NeoHookeParameters<NumberType> &other) const
{
return (
(StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ ==
other.strainEnergyFuncType_) &&
(fabs(shearModulus_ - other.shearModulus_) < Epsilon_Tol) &&
(fabs(poissonRatio_ - other.poissonRatio_) < Epsilon_Tol));
}
template <typename NumberType>
void NeoHookeParameters<NumberType>::Print_details(std::ostream &outStream) const
{
outStream << "Strain energy function type: "
<< StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_
<< std::endl;
outStream << "Shear modulus: " << this->shearModulus_ << std::endl;
outStream << "Poisson's ratio: " << this->poissonRatio_ << std::endl;
}
template <typename NumberType = double>
struct OgdenParameters : public StrainEnergyFuncParameters<NumberType>
{
NumberType volumetricConstant_;
NumberType bulkModulus_;
UnsignedIntType nOgdenTerms_;
std::vector<NumberType> shearModuli_;
std::vector<NumberType> shearConstants_;
OgdenParameters();
OgdenParameters(const OgdenParameters<NumberType> &other);
OgdenParameters<NumberType>& operator=(const OgdenParameters<NumberType> &other);
bool operator==(const OgdenParameters<NumberType> &other) const;
virtual void Print_details(std::ostream &outStream) const override;
template <class Archive>
void
serialize(Archive &ar, const unsigned int version)
{
ar &boost::serialization::base_object<StrainEnergyFuncParameters<NumberType>>(*this);
ar &volumetricConstant_;
ar &bulkModulus_;
ar &nOgdenTerms_;
ar &shearModuli_;
ar &shearConstants_;
}
template <class Archive>
void
deserialize(Archive &ar, const unsigned int version)
{
ar &boost::serialization::base_object<StrainEnergyFuncParameters<NumberType>>(*this);
ar &volumetricConstant_;
ar &bulkModulus_;
ar &nOgdenTerms_;
ar &shearModuli_;
ar &shearConstants_;
}
};
template <typename NumberType>
OgdenParameters<NumberType>::OgdenParameters()
: volumetricConstant_(-2.0)
, bulkModulus_(2.0e3)
, nOgdenTerms_(3)
, shearModuli_({0.63, -0.01, 0.0012})
, shearConstants_({1.3, -2.0, 5.0})
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = StrainEnergyFuncType::Og;
}
template <typename NumberType>
OgdenParameters<NumberType>::OgdenParameters(const OgdenParameters<NumberType> &other)
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = other.strainEnergyFuncType_;
this->volumetricConstant_ = other.volumetricConstant_;
this->bulkModulus_ = other.bulkModulus_;
this->nOgdenTerms_ = other.nOgdenTerms_;
this->shearModuli_.resize(this->nOgdenTerms_);
this->shearConstants_.resize(this->nOgdenTerms_);
for (UnsignedIntType isoTermCtr = 0; isoTermCtr < nOgdenTerms_;++isoTermCtr)
{
this->shearModuli_[isoTermCtr] = other.shearModuli_[isoTermCtr];
this->shearConstants_[isoTermCtr] = other.shearConstants_[isoTermCtr];
}
}
template <typename NumberType>
OgdenParameters<NumberType> &
OgdenParameters<NumberType>::operator=(const OgdenParameters<NumberType> &other)
{
if (this != &other) // confirm there is no self-assignment
{
StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ = other.strainEnergyFuncType_;
this->volumetricConstant_ = other.volumetricConstant_;
this->bulkModulus_ = other.bulkModulus_;
this->nOgdenTerms_ = other.nOgdenTerms_;
this->shearModuli_.resize(this->nOgdenTerms_);
this->shearConstants_.resize(this->nOgdenTerms_);
for (UnsignedIntType isoTermCtr = 0; isoTermCtr < nOgdenTerms_;++isoTermCtr)
{
this->shearModuli_[isoTermCtr] = other.shearModuli_[isoTermCtr];
this->shearConstants_[isoTermCtr] = other.shearConstants_[isoTermCtr];
}
}
return (*this);
}
template <typename NumberType>
bool OgdenParameters<NumberType>::operator==(const OgdenParameters<NumberType> &other) const
{
bool result =
(StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_ ==
other.strainEnergyFuncType_) &&
Check_ifEqual(volumetricConstant_, other.volumetricConstant_) &&
Check_ifEqual(bulkModulus_, other.bulkModulus_) &&
(nOgdenTerms_ == other.nOgdenTerms_);
if (result)
for (UnsignedIntType isoTermCtr = 0; isoTermCtr < nOgdenTerms_;
++isoTermCtr)
{
result = Check_ifEqual(shearModuli_[isoTermCtr],
other.shearModuli_[isoTermCtr]) &&
Check_ifEqual(shearConstants_[isoTermCtr],
other.shearConstants_[isoTermCtr]);
if (!result)
break;
}
return (result);
}
template <typename NumberType>
void OgdenParameters<NumberType>::Print_details(std::ostream &outStream) const
{
outStream << "Strain energy function type: "
<< StrainEnergyFuncParameters<NumberType>::strainEnergyFuncType_
<< std::endl;
outStream << "Volumetric constant: " << this->volumetricConstant_
<< std::endl;
outStream << "Bulk Modulus: " << this->bulkModulus_ << std::endl;
outStream << "Number of Ogden terms: " << this->nOgdenTerms_
<< std::endl;
outStream << "Shear moduli: ";
for (const auto &shearModulus : this->shearModuli_)
outStream << shearModulus << "|";
outStream << std::endl;
outStream << "Shear constants: ";
for (const auto &shearConstant : this->shearConstants_)
outStream << shearConstant << "|";
outStream << std::endl;
}
template <UnsignedIntType dim, typename NumberType = double>
struct ConstitutiveParameters
{
std::string materialFileName_;
std::string materialName_;
FloatType materialDensity_;
std::unique_ptr<StrainEnergyFuncParameters<NumberType>> strainEnergyFuncParams_;
ConstitutiveParameters(){};
~ConstitutiveParameters(){};
ConstitutiveParameters(const ConstitutiveParameters<dim, NumberType> &other);
ConstitutiveParameters<dim, NumberType>& operator=(const ConstitutiveParameters<dim, NumberType> &other);
bool operator==(const ConstitutiveParameters &other) const;
template <class Archive>
void
serialize(Archive &ar, const unsigned int version)
{
ar &materialFileName_;
ar &materialName_;
ar &materialDensity_;
if (strainEnergyFuncParams_)
{
int strainEnergyFuncTypeVal = static_cast<int>(
strainEnergyFuncParams_->strainEnergyFuncType_);
ar &strainEnergyFuncTypeVal;
std::cout << "my ranK: " << dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
<< " Inside ConstitutiveParameters::serialize(); strainEnergyFuncTypeVal: "
<< strainEnergyFuncTypeVal << " strainEnergyFuncType: "
<< strainEnergyFuncParams_->strainEnergyFuncType_
<< std::endl;
switch (strainEnergyFuncParams_->strainEnergyFuncType_)
{
case StrainEnergyFuncType::NH:
dynamic_cast<NeoHookeParameters<NumberType> &>(*strainEnergyFuncParams_).serialize(ar, version);
break;
case StrainEnergyFuncType::Og:
dynamic_cast<OgdenParameters<NumberType> &>(*strainEnergyFuncParams_).serialize(ar, version);
break;
default:
Assert(false,dealii::ExcMessage("Unsupported strain energy function type!!"));
}
}
}
template <class Archive>
void
deserialize(Archive &ar, const unsigned int version)
{
std::cout << "my ranK: " << dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)
<< "Entered ConstitutiveParameters::deserialize()"
<< std::endl;
ar &materialFileName_;
ar &materialName_;
ar &materialDensity_;
int strainEnergyFuncTypeVal;
ar & strainEnergyFuncTypeVal;
StrainEnergyFuncType strainEnergyFuncType =
static_cast<StrainEnergyFuncType>(
strainEnergyFuncTypeVal);
std::cout
<< "Inside ConstitutiveParameters::deserialize(); strainEnergyFuncTypeVal: "
<< strainEnergyFuncTypeVal << " strainEnergyFuncType: "
<< strainEnergyFuncType
<< std::endl;
switch (strainEnergyFuncType)
{
case StrainEnergyFuncType::NH:
{
NeoHookeParameters<NumberType> neoHookeParameters;
neoHookeParameters.deserialize(ar,version);
strainEnergyFuncParams_ = std::make_unique<NeoHookeParameters<NumberType>>(neoHookeParameters);
break;
}
case StrainEnergyFuncType::Og:
{
OgdenParameters<NumberType> ogdenParameters;
ogdenParameters.deserialize(ar,version);
strainEnergyFuncParams_ = std::make_unique<OgdenParameters<NumberType>>(ogdenParameters);
break;
}
default:
Assert(false,dealii::ExcMessage("Unsupported strain energy function type!!"));
}
}
};
template <UnsignedIntType dim, typename NumberType>
ConstitutiveParameters<dim, NumberType>::ConstitutiveParameters(const ConstitutiveParameters<dim, NumberType> &other)
{
this->materialFileName_ = other.materialFileName_;
this->materialName_ = other.materialName_;
this->materialDensity_ = other.materialDensity_;
switch (other.strainEnergyFuncParams_->strainEnergyFuncType_)
{
case StrainEnergyFuncType::NH:
this->strainEnergyFuncParams_ = std::make_unique<NeoHookeParameters<NumberType>>(*dynamic_cast<NeoHookeParameters<NumberType> *>(other.strainEnergyFuncParams_.get()));
break;
case StrainEnergyFuncType::Og:
this->strainEnergyFuncParams_ = std::make_unique<OgdenParameters<NumberType>>(*dynamic_cast<OgdenParameters<NumberType> *>(other.strainEnergyFuncParams_.get()));
break;
default:
Assert(false,dealii::ExcMessage("Unsupported strain energy function type!!"));
}
}
template <UnsignedIntType dim, typename NumberType>
ConstitutiveParameters<dim, NumberType> &
ConstitutiveParameters<dim, NumberType>::operator=(const ConstitutiveParameters<dim, NumberType> &other)
{
if (this != &other) // confirm there is no self-assignment
{
this->materialFileName_ = other.materialFileName_;
this->materialName_ = other.materialName_;
this->materialDensity_ = other.materialDensity_;
switch (other.strainEnergyFuncParams_->strainEnergyFuncType_)
{
case StrainEnergyFuncType::NH:
this->strainEnergyFuncParams_ = std::make_unique<NeoHookeParameters<NumberType>>(*dynamic_cast<NeoHookeParameters<NumberType> *>(other.strainEnergyFuncParams_.get()));
break;
case StrainEnergyFuncType::Og:
this->strainEnergyFuncParams_ = std::make_unique<OgdenParameters<NumberType>>(*dynamic_cast<OgdenParameters<NumberType> *>(other.strainEnergyFuncParams_.get()));
break;
default:
Assert(false,dealii::ExcMessage("Unsupported strain energy function type!!"));
}
}
return (*this);
}
template <UnsignedIntType dim, typename NumberType>
bool
ConstitutiveParameters<dim, NumberType>::
operator==(const ConstitutiveParameters<dim, NumberType> &other) const
{
bool result = (materialFileName_ == other.materialFileName_) &&
(materialName_ == other.materialName_) &&
(fabs(materialDensity_ - other.materialDensity_) < Epsilon_Tol);
Assert(strainEnergyFuncParams_,
dealii::ExcMessage(
"strainEnergyFuncParams_ pointer for self is null"));
Assert(other.strainEnergyFuncParams_,
dealii::ExcMessage(
"strainEnergyFuncParams_ pointer for other is null"));
result =
result && (strainEnergyFuncParams_->strainEnergyFuncType_ ==
other.strainEnergyFuncParams_->strainEnergyFuncType_);
if (result)
{
if (strainEnergyFuncParams_->strainEnergyFuncType_ == StrainEnergyFuncType::NH)
{
const NeoHookeParameters<NumberType>
*strainEnergyParamsSelf = dynamic_cast<NeoHookeParameters<NumberType> *>(this->strainEnergyFuncParams_.get()),
*strainEnergyParamsOther = dynamic_cast<NeoHookeParameters<NumberType> *>(other.strainEnergyFuncParams_.get());
result = *strainEnergyParamsSelf == *strainEnergyParamsOther;
}
else if (strainEnergyFuncParams_->strainEnergyFuncType_ == StrainEnergyFuncType::Og)
{
const OgdenParameters<NumberType>
*strainEnergyParamsSelf = dynamic_cast<OgdenParameters<NumberType> *>(this->strainEnergyFuncParams_.get()),
*strainEnergyParamsOther = dynamic_cast<OgdenParameters<NumberType> *>(other.strainEnergyFuncParams_.get());
result = *strainEnergyParamsSelf == *strainEnergyParamsOther;
}
}
return (result);
}
int
main(int argc, char *argv[])
{
const UnsignedIntType masterProcessRank = 0;
// create MPI object
dealii::Utilities::MPI::MPI_InitFinalize mpiInitialization(argc, argv, 1);
MPI_Comm const & mpiCommunicator(MPI_COMM_WORLD);
UnsignedIntType processorRank = dealii::Utilities::MPI::this_mpi_process(mpiCommunicator);
dealii::ConditionalOStream pCout(std::cout, processorRank == masterProcessRank);
const UnsignedIntType nProcs = dealii::Utilities::MPI::n_mpi_processes(mpiCommunicator);
pCout << "Running test with " << nProcs << " processes" << std::endl;
// prepare reference results
constexpr UnsignedIntType dim = 2;
std::string materialFileNameRef = "material.prm";
std::string materialNameRef = "material";
FloatType materialDensityRef = 1e-6;
StrainEnergyFuncType strainEnergyFuncTypeRef = NH;
NeoHookeParameters<FloatType> neoHookeParamsRef(StrainEnergyFuncType::NH);
{
neoHookeParamsRef.shearModulus_ = 225e3;
neoHookeParamsRef.poissonRatio_ = 0.3;
}
OgdenParameters<FloatType> ogdenParamsRef;
{
ogdenParamsRef.volumetricConstant_ = 12.0;
ogdenParamsRef.bulkModulus_ = 0.2;
ogdenParamsRef.nOgdenTerms_ = 2;
ogdenParamsRef.shearModuli_ = {25.5, 98.0};
ogdenParamsRef.shearConstants_ = {0.02, -0.00025};
}
// generate a default ConstitutiveParameters object and fill in its contents
ConstitutiveParameters<dim, FloatType> constiParamsRef;
{
constiParamsRef.materialFileName_ = materialFileNameRef;
constiParamsRef.materialName_ = materialNameRef;
constiParamsRef.materialDensity_ = materialDensityRef;
constiParamsRef.strainEnergyFuncParams_ = std::make_unique<OgdenParameters<FloatType>>(ogdenParamsRef);
}
if (processorRank == masterProcessRank)
{
std::cout << "Now displaying strainEnergyFuncParams_ form master proc: " << std::endl;
constiParamsRef.strainEnergyFuncParams_->Print_details(std::cout);
}
MPI_Barrier(mpiCommunicator);
// pack data on master process
std::vector<char> buffer;
UnsignedIntType bufferSize = 0;
if (processorRank == masterProcessRank)
{
buffer = dealii::Utilities::pack<ConstitutiveParameters<dim, FloatType>>(constiParamsRef);
bufferSize = buffer.size();
}
std::cout << "about to b_cast buffer size; my rank: " << processorRank << std::endl;
MPI_Barrier(mpiCommunicator);
// broadcast buffer size and resize buffer on receiving procs
MPI_Bcast(&bufferSize, 1, MPI_INT, masterProcessRank, mpiCommunicator);
if (processorRank != masterProcessRank)
buffer.resize(bufferSize);
std::cout << "about to b_cast the full buffer; my rank: " << processorRank << std::endl;
MPI_Barrier(mpiCommunicator);
// broadcast packed data to slave process and unpack there
MPI_Bcast(&buffer[0],
static_cast<int>(buffer.size()),
MPI_CHAR,
masterProcessRank,
mpiCommunicator);
std::cout << "buffer received; my rank: " << processorRank << std::endl;
MPI_Barrier(mpiCommunicator);
ConstitutiveParameters<dim, FloatType> constiParamsTest;
if (processorRank != masterProcessRank)
constiParamsTest = dealii::Utilities::unpack<ConstitutiveParameters<dim, FloatType>>(buffer);
std::cout << "unpacking done; my rank: " << processorRank << std::endl;
MPI_Barrier(mpiCommunicator);
if (processorRank == masterProcessRank + 1)
{
std::cout << "Now displaying strainEnergyFuncParams_ form a slave proc: " << std::endl;
constiParamsTest.strainEnergyFuncParams_->Print_details(std::cout);
}
MPI_Barrier(mpiCommunicator);
// check for equality
bool result = true;
if (processorRank != masterProcessRank)
if (!(constiParamsRef == constiParamsTest))
result = false;
// reduce and broadcast bool
MPI_Allreduce(&result, &result, 1, MPI_CXX_BOOL, MPI_LAND, mpiCommunicator);
if (result)
pCout << "Test for serialization of ConstitutiveParameters passed:)" << std::endl;
else
pCout << "Test for serialization of ConstitutiveParameters failed:(" << std::endl;
return (0);
}
