6 #include "parcer/Eigen.h"
12 DistributedCollection::DistributedCollection(std::shared_ptr<SampleCollection> collection, std::shared_ptr<parcer::Communicator> comm) :
SampleCollection(), collection(collection), comm(comm) {}
23 assert(i<collection->
size());
28 assert(i<collection->
size());
35 std::shared_ptr<SamplingState> state =
nullptr;
38 for(
unsigned int j=0; j<
comm->GetSize(); ++j ) {
40 comm->Bcast(localSize, j);
42 if( i<
size+localSize ) {
44 comm->Bcast(state, j);
57 std::shared_ptr<SamplingState> state =
nullptr;
60 for(
unsigned int j=0; j<
comm->GetSize(); ++j ) {
62 comm->Bcast(localSize, j);
64 if( i<
size+localSize ) {
66 comm->Bcast(state, j);
80 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
82 comm->Bcast(localSize, i);
119 const Eigen::VectorXd& local =
LocalESS(blockDim);
120 Eigen::VectorXd global = Eigen::VectorXd::Zero(local.size());
122 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
123 Eigen::VectorXd l(local.size());
124 if(
comm->GetRank()==i ) { l = local; }
139 Eigen::MatrixXd global(local.rows(),
GlobalSize());
142 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
144 comm->Bcast(localSize, i);
146 Eigen::MatrixXd l(local.rows(), localSize);
147 if(
comm->GetRank()==i ) { l = local; }
150 global.block(0, numSamps, local.rows(), localSize) = l;
152 numSamps += localSize;
164 Eigen::VectorXd global = Eigen::VectorXd::Constant(
GlobalSize(), std::numeric_limits<double>::quiet_NaN());
167 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
169 comm->Bcast(localSize, i);
171 Eigen::VectorXd l(localSize);
172 if(
comm->GetRank()==i ) { l = local; }
175 global.segment(numSamps, localSize) = l;
177 numSamps += localSize;
193 Eigen::VectorXd global = Eigen::VectorXd::Zero(f->outputSizes(0));
196 for(
unsigned int i=0; i<
comm->GetSize(); ++i ) {
197 Eigen::VectorXd l(f->outputSizes(0));
198 if(
comm->GetRank()==i ) { l = local; }
204 return global/(double)
comm->GetSize();
std::shared_ptr< parcer::Communicator > comm
The communicator for this collection.
Eigen::VectorXd LocalExpectedValue(std::shared_ptr< muq::Modeling::ModPiece > const &f, std::vector< std::string > const &metains=std::vector< std::string >()) const
std::shared_ptr< SampleCollection > collection
The local sample collection (stored on this processor)
virtual Eigen::VectorXd ESS(int blockDim=-1) const override
Compute the ESS using all of the samples.
Eigen::Matrix< scalar, rows, cols, options, maxRows, maxCols > GlobalEigenMean(Eigen::Matrix< scalar, rows, cols, options, maxRows, maxCols > const &local) const
virtual Eigen::VectorXd Weights() const override
Return all of the weights as a vector.
Eigen::MatrixXd GlobalCovariance(int blockDim=-1) const
Compute the covariance using all of the samples.
virtual void Add(std::shared_ptr< SamplingState > newSamp) override
Add a sample that is stored on this processor.
Eigen::VectorXd GlobalMean(int blockDim=-1) const
Compute the mean using all of the samples.
Eigen::VectorXd LocalESS(int blockDim=-1) const
Compute the ESS using only local samples.
virtual void WriteToFile(std::string const &filename, std::string const &dataset="/") const override
virtual std::shared_ptr< SamplingState > at(unsigned i) override
Get the local state at the index.
virtual Eigen::VectorXd ExpectedValue(std::shared_ptr< muq::Modeling::ModPiece > const &f, std::vector< std::string > const &metains=std::vector< std::string >()) const override
Eigen::MatrixXd AsGlobalMatrix(int blockDim=-1) const
Return all of the samples as a matrix.
Eigen::VectorXd LocalMean(int blockDim=-1) const
Compute the mean using only local samples.
virtual Eigen::VectorXd CentralMoment(unsigned order, int blockDim=-1) const override
Computes the componentwise central moments (e.g., variance, skewness, kurtosis, etc....
std::shared_ptr< SamplingState > GlobalAt(unsigned i)
Get the global state at the index.
virtual Eigen::VectorXd Mean(int blockDim=-1) const override
Compute the mean using all of the samples.
Eigen::VectorXd GlobalCentralMoment(unsigned order, int blockDim=-1) const
Computes the componentwise central moments (e.g., variance, skewness, kurtosis, etc....
unsigned int GlobalSize() const
The total number of samples.
virtual unsigned int size() const override
The total number of samples.
Eigen::VectorXd LocalCentralMoment(unsigned order, int blockDim=-1) const
Computes the componentwise central moments (e.g., variance, skewness, kurtosis, etc....
std::shared_ptr< SamplingState > LocalAt(unsigned i)
Get the local state at the index.
Eigen::VectorXd GlobalWeights() const
Return all of the weights as a vector.
virtual Eigen::MatrixXd AsMatrix(int blockDim=-1) const override
Return all of the samples as a matrix.
unsigned int LocalSize() const
The number of samples stored locally.
virtual Eigen::MatrixXd Covariance(int blockDim=-1) const override
Compute the covariance using all of the samples.
virtual Eigen::VectorXd Variance(int blockDim=-1) const override
Compute the variance using all of the samples.
Eigen::VectorXd LocalWeights() const
Return all of the weights on this processor as a vector.
Eigen::VectorXd LocalVariance(int blockDim=-1) const
Compute the variance using only local samples.
Eigen::VectorXd GlobalVariance(int blockDim=-1) const
Compute the variance using all of the samples.
Eigen::VectorXd GlobalESS(int blockDim=-1) const
Compute the ESS using all of the samples.
Eigen::MatrixXd AsLocalMatrix(int blockDim=-1) const
Return all of the samples on the local processor as a matrix.
Eigen::VectorXd GlobalExpectedValue(std::shared_ptr< muq::Modeling::ModPiece > const &f, std::vector< std::string > const &metains=std::vector< std::string >()) const
Eigen::MatrixXd LocalCovariance(int blockDim=-1) const
Compute the covariance using only local samples.
A class to hold and analyze a collection of SamplingState objects.