herwig
is hosted by
Hepforge
,
IPPP Durham
Herwig
7.3.0
Sampling
CellGrids
SimpleCellGrid.h
1
// -*- C++ -*-
2
//
3
// SimpleCellGrid.hpp is a part of ExSample
4
// Copyright (C) 2012-2019 Simon Platzer, The Herwig Collaboration
5
//
6
// ExSample is licenced under version 3 of the GPL, see COPYING for details.
7
//
8
9
#ifndef EXSAMPLE_SimpleCellGrid_hpp_included
10
#define EXSAMPLE_SimpleCellGrid_hpp_included
11
12
#include "CellGrid.h"
13
#include <cmath>
14
15
namespace
ExSample {
16
17
22
class
SimpleCellGrid
23
:
public
CellGrid
{
24
25
public
:
26
30
SimpleCellGrid
()
31
:
CellGrid
() {}
32
36
SimpleCellGrid
(
const
std::vector<double>& newLowerLeft,
37
const
std::vector<double>& newUpperRight,
38
bool
keepWeightInformation =
true
,
39
double
newWeight = 0.0);
40
44
virtual
CellGrid
*
makeInstance
()
const
;
45
49
virtual
CellGrid
*
makeInstance
(
const
std::vector<double>& newLowerLeft,
50
const
std::vector<double>& newUpperRight,
51
double
newWeight = 0.0)
const
;
52
57
virtual
void
split
(std::size_t newSplitDimension,
double
newSplitCoordinate);
58
59
60
virtual
void
splitter(
size_t
dim,
int
rat);
61
62
public
:
63
67
const
SimpleCellGrid
&
firstChild
()
const
{
68
return
dynamic_cast<
const
SimpleCellGrid
&
>
(
CellGrid::firstChild
());
69
}
70
74
SimpleCellGrid
&
firstChild
() {
75
return
dynamic_cast<
SimpleCellGrid
&
>
(
CellGrid::firstChild
());
76
}
77
81
const
SimpleCellGrid
&
secondChild
()
const
{
82
return
dynamic_cast<
const
SimpleCellGrid
&
>
(
CellGrid::secondChild
());
83
}
84
88
SimpleCellGrid
&
secondChild
() {
89
return
dynamic_cast<
SimpleCellGrid
&
>
(
CellGrid::secondChild
());
90
}
91
92
public
:
93
97
struct
Counter
{
98
102
Counter
()
103
:
nPoints
(0.0),
sumOfWeights
(0.0),
104
sumOfSquaredWeights
(0.0),
105
maxWeight
(0.0) {}
106
110
double
nPoints
;
111
115
double
sumOfWeights
;
116
120
double
sumOfSquaredWeights
;
121
125
double
maxWeight
;
126
130
void
book
(
double
weight
) {
131
nPoints
+= 1.0;
132
sumOfWeights
+= std::abs(
weight
);
133
sumOfSquaredWeights
+=
sqr
(
weight
);
134
maxWeight
= std::max(std::abs(
weight
),
maxWeight
);
135
}
136
140
double
averageWeight
()
const
{
return
nPoints
!= 0.0 ?
sumOfWeights
/
nPoints
: 0.0; }
141
145
double
varianceOfAverage
()
const
{
146
return
147
nPoints
> 1.0 ?
148
fabs(
sumOfSquaredWeights
/
nPoints
-
sqr
(
sumOfWeights
/
nPoints
))/(
nPoints
-1) : 0.0;
149
}
150
151
};
152
156
const
std::vector<std::pair<Counter,Counter> >&
weightInformation
()
const
{
return
theWeightInformation
; }
157
161
std::vector<std::pair<Counter,Counter> >&
weightInformation
() {
return
theWeightInformation
; }
162
166
virtual
void
updateWeightInformation
(
const
std::vector<double>& p,
167
double
w);
168
172
void
adjustReferenceWeight
(
double
w) {
173
theReferenceWeight
= std::max(
theReferenceWeight
,std::abs(w));
174
}
175
179
double
getReferenceWeight
()
const
{
180
return
theReferenceWeight
;
181
}
182
188
virtual
void
adapt
(
double
gain,
double
epsilon,
189
std::set<SimpleCellGrid*>& newCells);
190
195
virtual
void
setWeights
();
196
197
public
:
198
202
template
<
class
RndGenerator>
203
void
sampleFlatPoint
(std::vector<double>& p,
204
RndGenerator& rnd)
const
{
205
assert(p.size() ==
lowerLeft
().
size
());
206
for
(
size_t
k = 0; k < p.size(); ++k ) {
207
p[k] =
lowerLeft
()[k] + rnd.rnd()*(
upperRight
()[k]-
lowerLeft
()[k]);
208
}
209
}
210
214
template
<
class
RndGenerator>
215
void
sampleFlatPoint
(std::vector<double>& p,
216
const
std::vector<bool>& parameterFlags,
217
RndGenerator& rnd)
const
{
218
assert(p.size() ==
lowerLeft
().
size
());
219
for
(
size_t
k = 0; k < p.size(); ++k ) {
220
if
( parameterFlags[k] )
221
continue
;
222
p[k] =
lowerLeft
()[k] + rnd.rnd()*(
upperRight
()[k]-
lowerLeft
()[k]);
223
}
224
}
225
232
template
<
class
RndGenerator,
class
Function>
233
void
explore
(std::size_t nPoints,
234
RndGenerator& rnd,
235
Function& f,
236
std::set<SimpleCellGrid*>& newCells,
237
std::ostream& warn) {
238
unsigned
long
nanPoints = 0;
239
if
( !
isLeaf
() ) {
240
firstChild
().
explore
(nPoints,rnd,f,newCells,warn);
241
secondChild
().
explore
(nPoints,rnd,f,newCells,warn);
242
return
;
243
}
244
if
( !newCells.empty() ) {
245
if
( newCells.find(
this
) == newCells.end() )
246
return
;
247
}
248
std::vector<double> point(
lowerLeft
().
size
());
249
for
( std::size_t k = 0; k < nPoints; ++k ) {
250
sampleFlatPoint
(point,rnd);
251
double
w = f.evaluate(point);
252
if
( ! std::isfinite(w) ) {
253
++nanPoints;
254
continue
;
255
}
256
updateWeightInformation
(point,std::abs(w));
257
}
258
if
( nanPoints ) {
259
warn <<
"Warning: "
<< nanPoints <<
" out of "
260
<< nPoints <<
" points with nan or inf weight encountered while "
261
<<
"exploring a cell.\n"
<< std::flush;
262
}
263
}
264
268
template
<
class
RndGenerator>
269
SimpleCellGrid
*
selectCell
(RndGenerator& rnd) {
270
if
(
isLeaf
() )
271
return
this
;
272
if
(
firstChild
().
active
() &&
273
secondChild
().
active
() ) {
274
double
p =
firstChild
().
integral
()/
integral
();
275
if
( rnd.rnd() <= p )
276
return
firstChild
().
selectCell
(rnd);
277
else
278
return
secondChild
().
selectCell
(rnd);
279
}
280
if
(
firstChild
().
active
() &&
281
!
secondChild
().
active
() )
282
return
firstChild
().
selectCell
(rnd);
283
else
284
return
secondChild
().
selectCell
(rnd);
285
}
286
290
template
<
class
RndGenerator,
class
Function>
291
double
sample
(RndGenerator& rnd,
292
Function& f,
293
std::vector<double>& p,
294
bool
unweight,
295
bool
adjustReference) {
296
SimpleCellGrid
* selected =
selectCell
(rnd);
297
selected->
sampleFlatPoint
(p,rnd);
298
double
w = f.evaluate(p);
299
selected->
updateWeightInformation
(p,w);
300
double
xw =
integral
()*w/selected->
weight
();
301
if
( adjustReference ) {
302
selected->
adjustReferenceWeight
(xw);
303
}
304
if
( unweight ) {
305
double
r = selected->
getReferenceWeight
();
306
if
( r == 0. )
307
return
xw;
308
double
p = std::min(std::abs(xw),r)/r;
309
double
sign = xw >= 0. ? 1. : -1.;
310
if
( p < 1 && rnd.rnd() > p )
311
xw = 0.;
312
else
313
xw = sign*std::max(std::abs(xw),r);
314
}
315
return
xw;
316
}
317
321
template
<
class
RndGenerator,
class
Function>
322
std::pair<double,double>
generate
(RndGenerator& rnd,
323
Function& f,
324
std::vector<double>& p) {
325
SimpleCellGrid
* selected =
selectCell
(rnd);
326
selected->
sampleFlatPoint
(p,rnd);
327
double
w = f.evaluate(p);
328
selected->
updateWeightInformation
(p,w);
329
return
std::make_pair(w,selected->
weight
());
330
}
331
335
template
<
class
RndGenerator,
class
Function>
336
std::pair<double,double>
generate
(RndGenerator& rnd,
337
Function& f,
338
std::vector<double>& p,
339
const
std::vector<bool>& parameterFlags) {
340
SimpleCellGrid
* selected =
selectCell
(rnd);
341
selected->
sampleFlatPoint
(p,parameterFlags,rnd);
342
double
w = f.evaluate(p);
343
selected->
updateWeightInformation
(p,w);
344
return
std::make_pair(w,selected->
weight
());
345
}
346
347
public
:
348
352
virtual
void
fromXML
(
const
XML::Element
&);
353
357
virtual
XML::Element
toXML
()
const
;
358
359
private
:
360
364
std::vector<std::pair<Counter,Counter> >
theWeightInformation
;
365
369
double
theReferenceWeight
;
370
371
};
372
373
}
374
375
#endif
// EXSAMPLE_SimpleCellGrid_hpp_included
376
ExSample::CellGrid
A binary cell grid.
Definition:
CellGrid.h:42
ExSample::CellGrid::lowerLeft
const std::vector< double > & lowerLeft() const
Return the lower left corner of the cell grid.
Definition:
CellGrid.h:87
ExSample::CellGrid::firstChild
const CellGrid & firstChild() const
Return the first child.
ExSample::CellGrid::integral
double integral() const
Return the integral relevant for the last parameter point set.
ExSample::CellGrid::isLeaf
bool isLeaf() const
Return true, if this is a leaf in the tree.
Definition:
CellGrid.h:108
ExSample::CellGrid::active
bool active() const
Return true, if this grid is active with respect to the last parameter point passed.
ExSample::CellGrid::size
std::size_t size() const
Return the number of nodes contained in this grid.
ExSample::CellGrid::weight
double weight() const
Return the weight.
Definition:
CellGrid.h:178
ExSample::CellGrid::upperRight
const std::vector< double > & upperRight() const
Return the upper right corner of the cell grid.
Definition:
CellGrid.h:92
ExSample::CellGrid::secondChild
const CellGrid & secondChild() const
Return the second child.
ExSample::SimpleCellGrid
A simple cell grid providing basic adaption and sampling.
Definition:
SimpleCellGrid.h:23
ExSample::SimpleCellGrid::getReferenceWeight
double getReferenceWeight() const
Return the reference weight.
Definition:
SimpleCellGrid.h:179
ExSample::SimpleCellGrid::secondChild
SimpleCellGrid & secondChild()
Access the second child.
Definition:
SimpleCellGrid.h:88
ExSample::SimpleCellGrid::sampleFlatPoint
void sampleFlatPoint(std::vector< double > &p, const std::vector< bool > ¶meterFlags, RndGenerator &rnd) const
Sample a point flat in this cell, keeping parameters fixed.
Definition:
SimpleCellGrid.h:215
ExSample::SimpleCellGrid::makeInstance
virtual CellGrid * makeInstance(const std::vector< double > &newLowerLeft, const std::vector< double > &newUpperRight, double newWeight=0.0) const
Produce a new instance of a cell grid.
ExSample::SimpleCellGrid::generate
std::pair< double, double > generate(RndGenerator &rnd, Function &f, std::vector< double > &p)
Sample a point and return its weight.
Definition:
SimpleCellGrid.h:322
ExSample::SimpleCellGrid::SimpleCellGrid
SimpleCellGrid()
Default constructor.
Definition:
SimpleCellGrid.h:30
ExSample::SimpleCellGrid::firstChild
SimpleCellGrid & firstChild()
Access the first child.
Definition:
SimpleCellGrid.h:74
ExSample::SimpleCellGrid::weightInformation
const std::vector< std::pair< Counter, Counter > > & weightInformation() const
Return weight information for adaption steps.
Definition:
SimpleCellGrid.h:156
ExSample::SimpleCellGrid::firstChild
const SimpleCellGrid & firstChild() const
Return the first child.
Definition:
SimpleCellGrid.h:67
ExSample::SimpleCellGrid::sampleFlatPoint
void sampleFlatPoint(std::vector< double > &p, RndGenerator &rnd) const
Sample a point flat in this cell.
Definition:
SimpleCellGrid.h:203
ExSample::SimpleCellGrid::explore
void explore(std::size_t nPoints, RndGenerator &rnd, Function &f, std::set< SimpleCellGrid * > &newCells, std::ostream &warn)
Explore the cell grid, given a number of points to be sampled in each cell; the weights of the cell w...
Definition:
SimpleCellGrid.h:233
ExSample::SimpleCellGrid::adjustReferenceWeight
void adjustReferenceWeight(double w)
Adjust the reference weight.
Definition:
SimpleCellGrid.h:172
ExSample::SimpleCellGrid::weightInformation
std::vector< std::pair< Counter, Counter > > & weightInformation()
Access weight information for adaption steps.
Definition:
SimpleCellGrid.h:161
ExSample::SimpleCellGrid::fromXML
virtual void fromXML(const XML::Element &)
Fill CellGrid data from an XML element.
ExSample::SimpleCellGrid::makeInstance
virtual CellGrid * makeInstance() const
Produce a new instance of a cell grid.
ExSample::SimpleCellGrid::sample
double sample(RndGenerator &rnd, Function &f, std::vector< double > &p, bool unweight, bool adjustReference)
Sample a point and return its weight.
Definition:
SimpleCellGrid.h:291
ExSample::SimpleCellGrid::generate
std::pair< double, double > generate(RndGenerator &rnd, Function &f, std::vector< double > &p, const std::vector< bool > ¶meterFlags)
Sample a point and return its weight.
Definition:
SimpleCellGrid.h:336
ExSample::SimpleCellGrid::theReferenceWeight
double theReferenceWeight
The reference weight to be used for unweighting.
Definition:
SimpleCellGrid.h:369
ExSample::SimpleCellGrid::adapt
virtual void adapt(double gain, double epsilon, std::set< SimpleCellGrid * > &newCells)
Perform a default adaption step, splitting along the dimension which shows up the largest difference ...
ExSample::SimpleCellGrid::setWeights
virtual void setWeights()
Update the weights of the cells from information accumulated so far.
ExSample::SimpleCellGrid::SimpleCellGrid
SimpleCellGrid(const std::vector< double > &newLowerLeft, const std::vector< double > &newUpperRight, bool keepWeightInformation=true, double newWeight=0.0)
Construct given boundaries and a weight.
ExSample::SimpleCellGrid::toXML
virtual XML::Element toXML() const
Return an XML element for the data of this CellGrid.
ExSample::SimpleCellGrid::split
virtual void split(std::size_t newSplitDimension, double newSplitCoordinate)
Split this cell grid in the given dimension and coordinate, if it is a leaf.
ExSample::SimpleCellGrid::updateWeightInformation
virtual void updateWeightInformation(const std::vector< double > &p, double w)
Update the weight information for the given point.
ExSample::SimpleCellGrid::secondChild
const SimpleCellGrid & secondChild() const
Return the second child.
Definition:
SimpleCellGrid.h:81
ExSample::SimpleCellGrid::theWeightInformation
std::vector< std::pair< Counter, Counter > > theWeightInformation
Weight information for adaption steps.
Definition:
SimpleCellGrid.h:364
ExSample::SimpleCellGrid::selectCell
SimpleCellGrid * selectCell(RndGenerator &rnd)
Select a cell.
Definition:
SimpleCellGrid.h:269
XML::Element
Element represents a (tree of) XML elements.
Definition:
Element.h:56
sqr
constexpr auto sqr(const T &x) -> decltype(x *x)
ExSample::SimpleCellGrid::Counter
A simple counter to store information used for adaption.
Definition:
SimpleCellGrid.h:97
ExSample::SimpleCellGrid::Counter::averageWeight
double averageWeight() const
Return the average weight.
Definition:
SimpleCellGrid.h:140
ExSample::SimpleCellGrid::Counter::book
void book(double weight)
Book a point.
Definition:
SimpleCellGrid.h:130
ExSample::SimpleCellGrid::Counter::Counter
Counter()
Default constructor.
Definition:
SimpleCellGrid.h:102
ExSample::SimpleCellGrid::Counter::maxWeight
double maxWeight
The maximum weight.
Definition:
SimpleCellGrid.h:125
ExSample::SimpleCellGrid::Counter::sumOfWeights
double sumOfWeights
The sum of weights.
Definition:
SimpleCellGrid.h:115
ExSample::SimpleCellGrid::Counter::nPoints
double nPoints
The number of points.
Definition:
SimpleCellGrid.h:110
ExSample::SimpleCellGrid::Counter::sumOfSquaredWeights
double sumOfSquaredWeights
The sum of squared weights.
Definition:
SimpleCellGrid.h:120
ExSample::SimpleCellGrid::Counter::varianceOfAverage
double varianceOfAverage() const
Return the variance of the weights.
Definition:
SimpleCellGrid.h:145
Generated on Thu Jun 20 2024 17:50:53 for Herwig by
1.9.6