Skip to content

Commit 47cd7b9

Browse files
feat(EnergyTerm): add uuids for tracability
1 parent e582470 commit 47cd7b9

15 files changed

Lines changed: 430 additions & 210 deletions

include/geode/stochastic/sampling/direct/object_set_sampler/point_set_sampler.hpp

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -32,11 +32,10 @@
3232
namespace geode
3333
{
3434
template < index_t dimension >
35-
class UniformPointObjectSetSampler
36-
: public ObjectSetSampler< Point< dimension > >
35+
class UniformPointSetSampler : public ObjectSetSampler< Point< dimension > >
3736
{
3837
public:
39-
UniformPointObjectSetSampler(
38+
UniformPointSetSampler(
4039
const BoundingBox< dimension >& box, uuid subset_id )
4140
: ObjectSetSampler< Point< dimension > >{ subset_id }, box_( box )
4241
{

include/geode/stochastic/sampling/mcmc/metropolis_hasting_sampler.hpp

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -265,7 +265,8 @@ namespace geode
265265
"[MH] Change proposal has no index" );
266266
const auto delta_log_energy = energy_.delta_log_energy_change(
267267
state, proposal.old_object_id.value(),
268-
proposal.new_object.value().first );
268+
proposal.new_object.value().first,
269+
proposal.new_object.value().second );
269270
// should we test that objects are in the same group?
270271
// should be ensured by the dynamic
271272
return accept_or_reject( proposal, state, engine, delta_log_energy,

include/geode/stochastic/sampling/mcmc/models/components/energy_term.hpp

Lines changed: 93 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -26,92 +26,148 @@
2626
#include <geode/stochastic/spatial/object_set.hpp>
2727
#include <optional>
2828

29-
namespace
29+
namespace geode
3030
{
31-
struct NegLogParam
31+
namespace detail
3232
{
33-
explicit NegLogParam( double param )
33+
struct EnergyScale
3434
{
35-
OPENGEODE_EXCEPTION( param >= 0.,
36-
"[Gibbs energy term] - The model parameter "
37-
"cannot be smaller than 0." );
38-
if( param >= geode::GLOBAL_EPSILON )
35+
explicit EnergyScale( double param )
3936
{
40-
value = -std::log( param );
37+
OPENGEODE_EXCEPTION( param >= 0.,
38+
"[Gibbs energy term] - The model parameter "
39+
"cannot be negative." );
40+
41+
if( param >= geode::GLOBAL_EPSILON )
42+
{
43+
value = -std::log( param ); // store log-space parameter
44+
}
45+
// else value = std::nullopt → special case: param == 0
4146
}
42-
}
4347

44-
double scale( double multiplier ) const
45-
{
46-
if( value )
48+
/// Compute energy contribution for a given statistic multiplier
49+
double contribution( double multiplier ) const
4750
{
48-
return value.value() * multiplier;
51+
if( value )
52+
{
53+
return value.value() * multiplier;
54+
}
55+
// Hard constraint: param == 0
56+
return ( multiplier > 0 )
57+
? std::numeric_limits< double >::infinity()
58+
: 0.0;
4959
}
50-
return ( multiplier > 0 )
51-
? std::numeric_limits< double >::infinity()
52-
: 0.0;
53-
}
5460

55-
double param() const
56-
{
57-
if( value )
61+
/// Return original parameter (gamma)
62+
double parameter() const
5863
{
59-
return std::exp( -value.value() );
64+
if( value )
65+
{
66+
return std::exp( -value.value() );
67+
}
68+
return 0.;
6069
}
61-
return 0.;
62-
}
63-
std::optional< double > value;
64-
};
65-
} // namespace
70+
71+
private:
72+
std::optional< double > value; // empty if param == 0 (hardcore)
73+
};
74+
} // namespace detail
75+
} // namespace geode
6676

6777
namespace geode
6878
{
6979
template < typename Type >
7080
class EnergyTerm
7181
{
7282
public:
73-
explicit EnergyTerm( double parameter )
74-
: neg_log_parameter_( parameter )
83+
explicit EnergyTerm( std::string_view name, double param )
84+
: name_{ name }, energy_scale_{ param }
85+
{
86+
}
87+
88+
explicit EnergyTerm( std::string_view name,
89+
double param,
90+
const uuid& targeted_subset_id )
91+
: name_{ name },
92+
energy_scale_{ param },
93+
targeted_subset_id_{ targeted_subset_id }
7594
{
7695
}
7796

7897
virtual ~EnergyTerm() = default;
7998

80-
geode::uuid id() const
99+
const uuid& id() const
81100
{
82-
return id_;
101+
return energy_term_id_;
83102
}
84103

85-
const std::string& name() const
104+
std::string_view name() const
86105
{
87106
return name_;
88107
}
89108

90109
double parameter() const
91110
{
92-
return neg_log_parameter_.param();
111+
return energy_scale_.parameter();
112+
}
113+
114+
std::optional< uuid > targeted_subset_id() const
115+
{
116+
return targeted_subset_id_;
117+
}
118+
119+
/// Energy contribution for a given statistic multiplier
120+
double contribution( double multiplier ) const
121+
{
122+
return energy_scale_.contribution( multiplier );
93123
}
94124

95125
virtual double total_log( const ObjectSet< Type >& state ) const = 0;
96126

97127
virtual double delta_log_add( const ObjectSet< Type >& state,
98-
const Type& object,
99-
const uuid subset_id ) const = 0;
128+
const Type& new_object,
129+
const uuid& new_object_subset_id ) const = 0;
100130

101131
virtual double delta_log_remove(
102132
const ObjectSet< Type >& state, ObjectId object_id ) const = 0;
103133

104134
virtual double delta_log_change( const ObjectSet< Type >& state,
105135
ObjectId old_object_id,
106-
const Type& new_object ) const = 0;
136+
const Type& new_object,
137+
const uuid& new_object_subset_id ) const = 0;
107138

108139
virtual double statistic( const ObjectSet< Type >& state ) const = 0;
109140

110141
protected:
111-
NegLogParam neg_log_parameter_;
142+
bool is_targeted_subset( const uuid& subset_id ) const
143+
{
144+
return !targeted_subset_id_ || subset_id == *targeted_subset_id_;
145+
}
146+
147+
template < typename Func >
148+
void for_each_targeted_object(
149+
const ObjectSet< Type >& state, Func&& do_apply ) const
150+
{
151+
if( targeted_subset_id_ )
152+
{
153+
for( const auto id : geode::Range{
154+
state.nb_objects_in_subset( *targeted_subset_id_ ) } )
155+
{
156+
do_apply( { id, *targeted_subset_id_ } );
157+
}
158+
return;
159+
}
160+
for( const auto& id : state.get_all_object() )
161+
{
162+
do_apply( id );
163+
}
164+
}
112165

113166
private:
114-
geode::uuid id_;
115167
std::string name_;
168+
detail::EnergyScale energy_scale_;
169+
170+
std::optional< uuid > targeted_subset_id_{};
171+
uuid energy_term_id_{};
116172
};
117173
} // namespace geode

include/geode/stochastic/sampling/mcmc/models/components/intensity_term.hpp

Lines changed: 36 additions & 28 deletions
Original file line numberDiff line numberDiff line change
@@ -28,57 +28,65 @@
2828

2929
namespace geode
3030
{
31+
/// ObjectCountTerm
3132
template < typename Type >
3233
class IntensityTerm : public EnergyTerm< Type >
3334
{
3435
public:
35-
explicit IntensityTerm( double lambda, uuid subset_id )
36-
: EnergyTerm< Type >( lambda ), subset_id_{ subset_id }
36+
explicit IntensityTerm( std::string_view name, double lambda )
37+
: EnergyTerm< Type >( name, lambda )
3738
{
3839
}
3940

40-
double total_log( const ObjectSet< Type >& state ) const final
41+
explicit IntensityTerm(
42+
std::string_view name, double lambda, const uuid& subset_id )
43+
: EnergyTerm< Type >( name, lambda, subset_id )
4144
{
42-
const auto n = static_cast< double >( number_of_objects( state ) );
43-
return this->neg_log_parameter_.scale( n );
4445
}
4546

46-
double delta_log_add( const ObjectSet< Type >& state,
47-
const Type& sample,
48-
const uuid target_subset_id ) const final
47+
double total_log( const ObjectSet< Type >& state ) const override
4948
{
50-
return target_subset_id == subset_id_
51-
? this->neg_log_parameter_.scale( 1. )
52-
: 0.;
49+
const auto n = this->statistic( state );
50+
return this->contribution( n );
5351
}
5452

55-
double delta_log_remove(
56-
const ObjectSet< Type >& state, ObjectId object_id ) const final
53+
double delta_log_add( const ObjectSet< Type >& /*state*/,
54+
const Type& /*new_object*/,
55+
const uuid& new_object_subset_id ) const override
5756
{
58-
return object_id.subset == subset_id_
59-
? this->neg_log_parameter_.scale( -1. )
60-
: 0.;
57+
if( !this->is_targeted_subset( new_object_subset_id ) )
58+
{
59+
return 0.0;
60+
}
61+
return this->contribution( 1.0 );
6162
}
6263

63-
double delta_log_change( const ObjectSet< Type >& state,
64-
ObjectId old_object_id,
65-
const Type& new_sample ) const final
64+
double delta_log_remove( const ObjectSet< Type >& /*state*/,
65+
ObjectId object_id ) const override
6666
{
67-
return 0.0;
67+
if( !this->is_targeted_subset( object_id.subset ) )
68+
{
69+
return 0.0;
70+
}
71+
return this->contribution( -1.0 );
6872
}
6973

70-
double statistic( const ObjectSet< Type >& state ) const final
74+
double delta_log_change( const ObjectSet< Type >& /*state*/,
75+
ObjectId /*old_object_id*/,
76+
const Type& /*new_object*/,
77+
const uuid& /*new_object_subset_id*/ ) const override
7178
{
72-
return static_cast< double >( number_of_objects( state ) );
79+
return 0.0;
7380
}
7481

75-
private:
76-
index_t number_of_objects( const ObjectSet< Type >& state ) const
82+
double statistic( const ObjectSet< Type >& state ) const override
7783
{
78-
return state.nb_objects_in_subset( subset_id_ );
84+
if( this->targeted_subset_id() )
85+
{
86+
return static_cast< double >( state.nb_objects_in_subset(
87+
this->targeted_subset_id().value() ) );
88+
}
89+
return static_cast< double >( state.nb_objects() );
7990
}
80-
81-
private:
82-
uuid subset_id_;
8391
};
8492
} // namespace geode

0 commit comments

Comments
 (0)