Hello fellow palabos users,
dear developers,
For my current project, I want to monitor the development of wall shear stress during a simulation in the form of a histogram.
At first I thought, I just pass the adress of a single std:vector to my processing functional which would grow as all processes ‘push_back’ to it:
[code=“cpp”]
/* ******** MaskedBoxStressNormVectorFunctional2D.h ******* /
template<typename T1, typename T2>
class MaskedBoxStressNormVectorFunctional2D : public BoxProcessingFunctional2D_SS<T1, T2>
{
public:
MaskedBoxStressNormVectorFunctional2D(std::vector wallStressVec_, int flag_);
virtual void process( Box2D domain,
ScalarField2D& mask,
ScalarField2D& stressNormField );
virtual MaskedBoxStressNormVectorFunctional2D<T1,T2>* clone() const;
virtual void getTypeOfModification(std::vectormodif::ModifT& modified) const {
modified[0] = modif::nothing;
modified[1] = modif::nothing;
}
private:
T1 flag;
std::vector* wallStressVec;
};
/* ******** MaskedBoxStressNormVectorFunctional2D.hh ******* /
template<typename T1, typename T2>
MaskedBoxStressNormVectorFunctional2D<T1,T2>::MaskedBoxStressNormVectorFunctional2D(std::vector wallStressVec_, int flag_)
: flag(flag_), wallStressVec(wallStressVec_)
{ }
template<typename T1, typename T2>
void MaskedBoxStressNormVectorFunctional2D<T1,T2>::process (
Box2D domain,
ScalarField2D& mask,
ScalarField2D& stressNormField )
{
Dot2D offset = computeRelativeDisplacement(mask, stressNormField);
for (plint iX=domain.x0; iX<=domain.x1; ++iX) {
for (plint iY=domain.y0; iY<=domain.y1; ++iY) {
if (mask.get(iX,iY)==flag) {
// pcout << "get wall stress, which is = " << stressNormField.get(iX,iY) << std::endl;
wallStressVec->push_back(stressNormField.get(iX+offset.x, iY+offset.y));
// pcout << " process: size of wall stress vec = " << wallStressVec->size() << std::endl;
}
}
}
}
template<typename T1, typename T2>
MaskedBoxStressNormVectorFunctional2D<T1,T2>* MaskedBoxStressNormVectorFunctional2D<T1,T2>::clone() const
{
return new MaskedBoxStressNormVectorFunctional2D<T1,T2>(*this);
}
For 1 or 2 processes this seem to work fine.
However, as I increase the number of processes, strange things happen and I end up with a segfault.
I think this is due to the vector adress not being accessible for all processes.
Palabos uses its 'reductive processing functionals' to calculate/compare values from across all processes- this seems to be the preferred approach.
However, only the statistics 'average', 'sum', 'min' & 'max' are implemented in the palabos source code...
So I thought I might just as well add a new functionality:
I tried to adapt the implemented [i]MaskedBoxScalarAverageFunctional2D[/i] in a way so it adds all values to a private vector which should be accessible using a public function similar to [i]getAverageScalar()[/i].
[code="cpp"]
/* ******** MaskedBoxScalarHistogramFunctional2D ******* */
template<typename T>
class MaskedBoxScalarHistogramFunctional2D : public ReductiveBoxProcessingFunctional2D_SS<T,int> {
public:
MaskedBoxScalarHistogramFunctional2D(int flag_)
: histogramScalarId(this->getStatistics().subscribeHistogram()),
flag(flag_)
{ }
void process ( Box2D domain,
ScalarField2D<T>& scalarField,
ScalarField2D<int>& mask ) {
Dot2D offset = computeRelativeDisplacement(scalarField, mask);
BlockStatistics& statistics = this->getStatistics();
for (plint iX=domain.x0; iX<=domain.x1; ++iX) {
for (plint iY=domain.y0; iY<=domain.y1; ++iY) {
if (mask.get(iX+offset.x, iY+offset.y)==flag) {
statistics.gatherHistogram(histogramScalarId, (T)scalarField.get(iX,iY));
statistics.incrementStats();
}
}
}
}
MaskedBoxScalarHistogramFunctional2D<T>* clone() const {
return new MaskedBoxScalarHistogramFunctional2D<T>(*this);
}
virtual void getTypeOfModification(std::vector<modif::ModifT>& modified) const {
modified[0] = modif::nothing;
modified[1] = modif::nothing;
}
std::vector<T> getHistogramScalar() const {
std::vector<T> histogram = this->getStatistics().getHistogram(histogramScalarId);
// The histogram is internally computed on floating-point values. If T is
// integer, the value must be rounded at the end.
if (std::numeric_limits<T>::is_integer) {
return std::vector<T>(histogram.begin(), histogram.end());
}
return histogram;
}
private:
plint histogramScalarId;
int flag;
};
In the course of this I found myself changing the structure of ReductiveBoxProcessingFunctionals in general, and now I’m stuck somewhere in the deepest, darkest parts of palabos: the generic implementation of how MPI is used to communicate between processes “src/parallelism/parallelStatistics.cpp”
So instead of learning how to use MPI and altering the palabos source code (another time maybe), I hope that one of you might suggest a simpler, more elegant approach.
How would you gather values from a ScalarField (e.g. norm of shear stress) from all processes in a single std::vector ?
Cheers,
david