MUQ  0.4.3
DistributedCollection.cpp
Go to the documentation of this file.
2 
3 #if MUQ_HAS_MPI
4 
5 #if MUQ_HAS_PARCER
6 #include "parcer/Eigen.h"
7 #endif
8 
9 using namespace muq::Utilities;
10 using namespace muq::SamplingAlgorithms;
11 
12 DistributedCollection::DistributedCollection(std::shared_ptr<SampleCollection> collection, std::shared_ptr<parcer::Communicator> comm) : SampleCollection(), collection(collection), comm(comm) {}
13 
14 void DistributedCollection::Add(std::shared_ptr<SamplingState> newSamp) {
15  collection->Add(newSamp);
16 }
17 
18 std::shared_ptr<SamplingState> DistributedCollection::at(unsigned i) { return GlobalAt(i); }
19 
20 const std::shared_ptr<SamplingState> DistributedCollection::at(unsigned i) const { return GlobalAt(i); }
21 
22 std::shared_ptr<SamplingState> DistributedCollection::LocalAt(unsigned i) {
23  assert(i<collection->size());
24  return collection->at(i);
25 }
26 
27 const std::shared_ptr<SamplingState> DistributedCollection::LocalAt(unsigned i) const {
28  assert(i<collection->size());
29  return collection->at(i);
30 }
31 
32 std::shared_ptr<SamplingState> DistributedCollection::GlobalAt(unsigned i) {
33  assert(i<GlobalSize());
34 
35  std::shared_ptr<SamplingState> state = nullptr;
36 
37  int size = 0;
38  for( unsigned int j=0; j<comm->GetSize(); ++j ) {
39  int localSize = j==comm->GetRank() ? LocalSize() : 0;
40  comm->Bcast(localSize, j);
41 
42  if( i<size+localSize ) {
43  state = j==comm->GetRank() ? LocalAt(i-size) : nullptr;
44  comm->Bcast(state, j);
45  break;
46  }
47 
48  size += localSize;
49  }
50 
51  return state;
52 }
53 
54 const std::shared_ptr<SamplingState> DistributedCollection::GlobalAt(unsigned i) const {
55  assert(i<GlobalSize());
56 
57  std::shared_ptr<SamplingState> state = nullptr;
58 
59  int size = 0;
60  for( unsigned int j=0; j<comm->GetSize(); ++j ) {
61  int localSize = j==comm->GetRank() ? LocalSize() : 0;
62  comm->Bcast(localSize, j);
63 
64  if( i<size+localSize ) {
65  state = j==comm->GetRank() ? LocalAt(i-size) : nullptr;
66  comm->Bcast(state, j);
67  break;
68  }
69 
70  size += localSize;
71  }
72 
73  return state;
74 }
75 
76 unsigned int DistributedCollection::LocalSize() const { return collection->size(); }
77 
78 unsigned int DistributedCollection::GlobalSize() const {
79  int size = 0;
80  for( unsigned int i=0; i<comm->GetSize(); ++i ) {
81  int localSize = i==comm->GetRank() ? LocalSize() : 0;
82  comm->Bcast(localSize, i);
83 
84  size += localSize;
85  }
86 
87  return size;
88 }
89 
90 unsigned int DistributedCollection::size() const { return GlobalSize(); }
91 
92 Eigen::VectorXd DistributedCollection::LocalCentralMoment(unsigned order, int blockDim) const { return collection->CentralMoment(order, GlobalMean(blockDim), blockDim); }
93 
94 Eigen::VectorXd DistributedCollection::GlobalCentralMoment(unsigned order, int blockDim) const { return GlobalEigenMean(LocalCentralMoment(order, blockDim)); }
95 
96 Eigen::VectorXd DistributedCollection::CentralMoment(unsigned order, int blockDim) const { return GlobalCentralMoment(order, blockDim); }
97 
98 Eigen::VectorXd DistributedCollection::LocalMean(int blockDim) const { return collection->Mean(blockDim); }
99 
100 Eigen::VectorXd DistributedCollection::GlobalMean(int blockDim) const { return GlobalEigenMean(LocalMean(blockDim)); }
101 
102 Eigen::VectorXd DistributedCollection::Mean(int blockDim) const { return GlobalMean(blockDim); }
103 
104 Eigen::VectorXd DistributedCollection::LocalVariance(int blockDim) const { return collection->Variance(blockDim); }
105 
106 Eigen::VectorXd DistributedCollection::GlobalVariance(int blockDim) const { return GlobalEigenMean(LocalVariance(blockDim)); }
107 
108 Eigen::VectorXd DistributedCollection::Variance(int blockDim) const { return GlobalVariance(blockDim); }
109 
110 Eigen::MatrixXd DistributedCollection::LocalCovariance(int blockDim) const { return collection->Covariance(GlobalMean(blockDim), blockDim); }
111 
112 Eigen::MatrixXd DistributedCollection::GlobalCovariance(int blockDim) const { return GlobalEigenMean(LocalCovariance(blockDim)); }
113 
114 Eigen::MatrixXd DistributedCollection::Covariance(int blockDim) const { return GlobalCovariance(blockDim); }
115 
116 Eigen::VectorXd DistributedCollection::LocalESS(int blockDim) const { return collection->ESS(blockDim); }
117 
118 Eigen::VectorXd DistributedCollection::GlobalESS(int blockDim) const {
119  const Eigen::VectorXd& local = LocalESS(blockDim);
120  Eigen::VectorXd global = Eigen::VectorXd::Zero(local.size());
121 
122  for( unsigned int i=0; i<comm->GetSize(); ++i ) {
123  Eigen::VectorXd l(local.size());
124  if( comm->GetRank()==i ) { l = local; }
125  comm->Bcast(l, i);
126 
127  global += l;
128  }
129 
130  return global;
131 }
132 
133 Eigen::VectorXd DistributedCollection::ESS(int blockDim) const { return GlobalESS(blockDim); }
134 
135 Eigen::MatrixXd DistributedCollection::AsLocalMatrix(int blockDim) const { return collection->AsMatrix(blockDim); }
136 
137 Eigen::MatrixXd DistributedCollection::AsGlobalMatrix(int blockDim) const {
138  const Eigen::MatrixXd& local = AsLocalMatrix(blockDim);
139  Eigen::MatrixXd global(local.rows(), GlobalSize());
140 
141  int numSamps = 0;
142  for( unsigned int i=0; i<comm->GetSize(); ++i ) {
143  int localSize = i==comm->GetRank() ? LocalSize() : 0;
144  comm->Bcast(localSize, i);
145 
146  Eigen::MatrixXd l(local.rows(), localSize);
147  if( comm->GetRank()==i ) { l = local; }
148 
149  comm->Bcast(l, i);
150  global.block(0, numSamps, local.rows(), localSize) = l;
151 
152  numSamps += localSize;
153  }
154 
155  return global;
156 }
157 
158 Eigen::MatrixXd DistributedCollection::AsMatrix(int blockDim) const { return AsGlobalMatrix(blockDim); }
159 
160 Eigen::VectorXd DistributedCollection::LocalWeights() const { return collection->Weights(); }
161 
162 Eigen::VectorXd DistributedCollection::GlobalWeights() const {
163  const Eigen::VectorXd& local = LocalWeights();
164  Eigen::VectorXd global = Eigen::VectorXd::Constant(GlobalSize(), std::numeric_limits<double>::quiet_NaN());
165 
166  int numSamps = 0;
167  for( unsigned int i=0; i<comm->GetSize(); ++i ) {
168  int localSize = i==comm->GetRank() ? LocalSize() : 0;
169  comm->Bcast(localSize, i);
170 
171  Eigen::VectorXd l(localSize);
172  if( comm->GetRank()==i ) { l = local; }
173 
174  comm->Bcast(l, i);
175  global.segment(numSamps, localSize) = l;
176 
177  numSamps += localSize;
178  }
179 
180  return global;
181 }
182 
183 Eigen::VectorXd DistributedCollection::Weights() const {
184  return GlobalWeights();
185 }
186 
187 void DistributedCollection::WriteToFile(std::string const& filename, std::string const& dataset) const {
188  collection->WriteToFile(filename, dataset);
189 }
190 
191 Eigen::VectorXd DistributedCollection::GlobalExpectedValue(std::shared_ptr<muq::Modeling::ModPiece> const& f, std::vector<std::string> const& metains) const {
192  const Eigen::VectorXd& local = LocalExpectedValue(f, metains);
193  Eigen::VectorXd global = Eigen::VectorXd::Zero(f->outputSizes(0));
194 
195  int numSamps = 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; }
199 
200  comm->Bcast(l, i);
201  global += l;
202  }
203 
204  return global/(double)comm->GetSize();
205 }
206 
207 Eigen::VectorXd DistributedCollection::LocalExpectedValue(std::shared_ptr<muq::Modeling::ModPiece> const& f, std::vector<std::string> const& metains) const {
208  return collection->ExpectedValue(f, metains);
209 }
210 
211 Eigen::VectorXd DistributedCollection::ExpectedValue(std::shared_ptr<muq::Modeling::ModPiece> const& f, std::vector<std::string> const& metains) const {
212  return GlobalExpectedValue(f, metains);
213 }
214 
215 #endif // end MUQ_HAS_MPI
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.