Mcmc config

bindings/python/src/stochastic/CMakeLists.txt

@@ -21,7 +21,7 @@
 add_geode_python_binding(
     NAME "py_stochastic"
     SOURCES
-        "sampling/mcmc/helpers/fracture_simulation_runner.hpp"
+       # "sampling/mcmc/helpers/fracture_simulation_runner.hpp"
         "sampling/mcmc/helpers/simulation_monitor.hpp"
         "sampling/mcmc/helpers/simulation_printer.hpp"
         "sampling/mcmc/simulation_runner.hpp"

bindings/python/src/stochastic/sampling/mcmc/helpers/fracture_simulation_runner.hpp

@@ -27,7 +27,7 @@
 {
     void define_fracture_simulation( pybind11::module& module )
     {
         using namespace geode;
 
         pybind11::class_< FractureSetDescription >(
             module, "FractureSetDescription" )
@@ -88,7 +88,7 @@
             //                pybind11::arg( "engine" ), pybind11::arg( "steps"
             //                ), "Run simulation for a fixed number of steps." )
             .def( "run",
-                static_cast< StatisticsMonitor ( FractureSimulationRunner::* )(
+                static_cast< StatisticsTracker ( FractureSimulationRunner::* )(
                     RandomEngine&, const SimulationConfigurator& ) >(
                     &FractureSimulationRunner::run ),
                 pybind11::arg( "engine" ), pybind11::arg( "config" ),

bindings/python/src/stochastic/sampling/mcmc/helpers/simulation_monitor.hpp

@@ -21,30 +21,30 @@
  *
  */
 
-#include <geode/stochastic/sampling/mcmc/helpers/simulation_monitor.hpp>
+#include <geode/stochastic/inference/statistics_tracker.hpp>
 
 namespace geode
 {
     void define_simulation_monitor( pybind11::module &module )
     {
-        pybind11::class_< geode::StatisticsMonitor >(
-            module, "StatisticsMonitor" )
+        pybind11::class_< geode::StatisticsTracker >(
+            module, "StatisticsTracker" )
             .def( pybind11::init< geode::index_t >(),
                 pybind11::arg( "nb_energy_terms" ),
-                "Create a StatisticsMonitor for a given number of energy "
+                "Create a StatisticsTracker for a given number of energy "
                 "terms" )
-            .def( "add_realization", &geode::StatisticsMonitor::add_realization,
+            .def( "add_realization", &geode::StatisticsTracker::add_realization,
                 pybind11::arg( "values" ),
                 "Add a realization (vector of doubles) to update statistics" )
-            .def( "statiscal_count", &geode::StatisticsMonitor::statiscal_count,
+            .def( "statiscal_count", &geode::StatisticsTracker::statiscal_count,
                 "Return the number of realizations added" )
-            .def_property_readonly( "means", &geode::StatisticsMonitor::means,
+            .def_property_readonly( "means", &geode::StatisticsTracker::means,
                 "Return the computed mean values for each energy term" )
             .def_property_readonly( "variances",
-                &geode::StatisticsMonitor::variances,
+                &geode::StatisticsTracker::variances,
                 "Return the computed variances for each energy term" )
-            .def( "__repr__", []( const geode::StatisticsMonitor &self ) {
-                return "<StatisticsMonitor count="
+            .def( "__repr__", []( const geode::StatisticsTracker &self ) {
+                return "<StatisticsTracker count="
                        + std::to_string( self.statiscal_count() ) + ">";
             } );
     }

bindings/python/src/stochastic/sampling/mcmc/helpers/simulation_printer.hpp

@@ -64,7 +64,7 @@ namespace geode
         //                &SimulationPrinter::print_statistics_summary,
         //                pybind11::arg( "monitor" ),
         //                pybind11::arg( "energy_term_names" ) = "",
-        //                "Print statistics summary from a StatisticsMonitor."
+        //                "Print statistics summary from a StatisticsTracker."
         //                );
     }
 } // namespace geode

bindings/python/src/stochastic/stochastic.cpp

@@ -26,8 +26,8 @@
 
 #include "sampling/direct/double_sampler.hpp"
 
-#include "sampling/mcmc/helpers/fracture_simulation_runner.hpp"
-#include "sampling/mcmc/helpers/simulation_monitor.hpp"
+// #include "sampling/mcmc/helpers/fracture_simulation_runner.hpp"
+// #include "sampling/mcmc/helpers/simulation_monitor.hpp"
 #include "sampling/mcmc/helpers/simulation_printer.hpp"
 #include "sampling/mcmc/simulation_runner.hpp"
 
@@ -49,8 +49,8 @@ PYBIND11_MODULE( opengeode_stochastic_py_stochastic, module )
     geode::define_random_engine( module );
     geode::define_double_sampler( module );
 
-    geode::define_simulation_monitor( module );
+    // geode::define_simulation_monitor( module );
     geode::define_simulation_printer( module );
     geode::define_simulation_runner( module );
-    geode::define_fracture_simulation( module );
+    // geode::define_fracture_simulation( module );
 }

bindings/python/tests/stochastic/CMakeLists.txt

@@ -18,8 +18,8 @@
 # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
 # SOFTWARE.
 
-add_geode_python_test(
-    SOURCE "test-py-mh-fractures.py"
-    DEPENDENCIES
-        ${PROJECT_NAME}::py_stochastic
-)
+#add_geode_python_test(
+#    SOURCE "test-py-mh-fractures.py"
+#    DEPENDENCIES
+#        ${PROJECT_NAME}::py_stochastic
+#)

include/geode/stochastic/inference/statistics_tools.hpp

@@ -0,0 +1,76 @@
+/*
+ * Copyright (c) 2019 - 2026 Geode-solutions
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+ * SOFTWARE.
+ *
+ */
+#pragma once
+
+#include <absl/container/flat_hash_map.h>
+
+#include <geode/stochastic/inference/statistics_tracker.hpp>
+#include <geode/stochastic/inference/target_statistics.hpp>
+
+namespace geode
+{
+    namespace statistics
+    {
+        template < typename ObjectType >
+        void validate( const StatisticsTracker< ObjectType >& tracker,
+            const TargetStatistics< ObjectType >& targets )
+        {
+            const auto& model = targets.model();
+
+            for( const auto& term_uuid : targets.active_terms() )
+            {
+                const auto mean = tracker.mean( term_uuid );
+                const auto target = targets.target( term_uuid );
+
+                const auto rel_error =
+                    std::fabs( mean - target )
+                    / ( std::fabs( target ) + geode::GLOBAL_EPSILON );
+
+                OpenGeodeStochasticStochasticException::check_exception(
+                    rel_error < targets.tolerance( term_uuid ), nullptr,
+                    OpenGeodeException::TYPE::result,
+                    "[StatisticsValidator] Failure for term ",
+                    model.term_name( term_uuid ), "\n  mean    = ", mean,
+                    "\n  target  = ", target, "\n  error   = ", rel_error,
+                    "\n  tol     = ", targets.tolerance( term_uuid ) );
+            }
+        }
+
+        template < typename ObjectType >
+        double quadratic_loss( const StatisticsTracker< ObjectType >& tracker,
+            const TargetStatistics< ObjectType >& targets )
+        {
+            double loss = 0.0;
+
+            for( const auto& term_uuid : targets.active_terms() )
+            {
+                const auto diff =
+                    tracker.mean( term_uuid ) - targets.value( term_uuid );
+
+                loss += diff * diff;
+            }
+
+            return loss;
+        }
+    } // namespace statistics
+} // namespace geode

include/geode/stochastic/inference/statistics_tracker.hpp

@@ -0,0 +1,86 @@
+#pragma once
+
+// #include <absl/container/flat_hash_map.h>
+#include <geode/basic/common.hpp>
+#include <geode/basic/uuid.hpp>
+
+#include <geode/stochastic/common.hpp>
+#include <geode/stochastic/models/model.hpp>
+
+namespace geode
+{
+    template < typename ObjectType >
+    class StatisticsTracker
+    {
+    public:
+        StatisticsTracker( const Model< ObjectType >& model ) : model_{ model }
+        {
+            means_.resize( model.nb_terms(), 0.0 );
+            m2_.resize( model.nb_terms(), 0.0 );
+        }
+
+        [[nodiscard]] index_t statiscal_count() const
+        {
+            return count_;
+        }
+
+        void add_realization( const std::vector< double >& values )
+        {
+            ++count_;
+            for( const auto value_id : geode::Range{ values.size() } )
+            {
+                auto& value = values[value_id];
+                auto& mean = means_[value_id];
+                auto& sum_of_squares = m2_[value_id];
+
+                const auto delta = value - mean;
+                mean += delta / count_;
+                const auto delta2 = value - mean;
+                sum_of_squares += delta * delta2;
+            }
+        }
+
+        [[nodiscard]] double mean( const uuid& term_uuid ) const
+        {
+            return means_[model_.term_index( term_uuid )];
+        }
+
+        [[nodiscard]] const std::vector< double >& means() const
+        {
+            return means_;
+        }
+
+        [[nodiscard]] double variance( const uuid& term_uuid ) const
+        {
+            return variance( model_.term_index( term_uuid ) );
+        }
+
+        [[nodiscard]] std::vector< double > variances() const
+        {
+            std::vector< double > variances;
+            variances.reserve( model_.nb_terms() );
+            for( const auto variance_id : geode::Range{ model_.nb_terms() } )
+            {
+                variances.emplace_back( this->variance( variance_id ) );
+            }
+            return variances;
+        }
+
+    private:
+        [[nodiscard]] double variance( index_t term_index ) const
+        {
+            if( count_ < 2 )
+            {
+                return 0.0;
+            }
+            return m2_[term_index] / ( count_ - 1 );
+        }
+
+    private:
+        const Model< ObjectType >& model_;
+
+        std::vector< double > means_;
+        std::vector< double > m2_;
+        index_t count_{ 0 };
+    };
+} // namespace geode

include/geode/stochastic/inference/target_statistics.hpp

@@ -0,0 +1,100 @@
+/*
+ * Copyright (c) 2019 - 2026 Geode-solutions
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+ * SOFTWARE.
+ *
+ */
+#pragma once
+#include <geode/basic/uuid.hpp>
+
+namespace geode
+{
+    struct TargetStatisticConfig
+    {
+        std::string term_name;
+        double value;
+        double tolerance;
+    };
+
+    template < typename ObjectType >
+    class TargetStatistics
+    {
+    public:
+        explicit TargetStatistics( const Model< ObjectType >& model )
+            : model_( model )
+        {
+            values_.resize( model.nb_terms(), 0.0 );
+            tolerances_.resize( model.nb_terms(), 0.0 );
+            active_.resize( model.nb_terms(), false );
+        }
+
+        const Model< ObjectType >& model() const
+        {
+            return model_;
+        }
+
+        void set_target( const TargetStatisticConfig& statistic )
+        {
+            const auto term_uuid =
+                model_.terms().get_term_uuid( statistic.term_name );
+            const auto idx = model_.term_index( term_uuid );
+
+            values_[idx] = statistic.value;
+            tolerances_[idx] = statistic.tolerance;
+            active_[idx] = true;
+        }
+
+        bool has_target( const uuid& term_uuid ) const
+        {
+            return active_[model_.term_index( term_uuid )];
+        }
+
+        double target( const uuid& term_uuid ) const
+        {
+            return values_[model_.term_index( term_uuid )];
+        }
+
+        double tolerance( const uuid& term_uuid ) const
+        {
+            return tolerances_[model_.term_index( term_uuid )];
+        }
+
+        std::vector< uuid > active_terms() const
+        {
+            std::vector< uuid > active_terms_uuid;
+
+            for( const auto& term : model_.terms().energy_terms() )
+            {
+                const auto& id = term->id();
+                if( active_[model_.term_index( id )] )
+                {
+                    active_terms_uuid.push_back( id );
+                }
+            }
+            return active_terms_uuid;
+        }
+
+    private:
+        const Model< ObjectType >& model_;
+
+        std::vector< double > values_;
+        std::vector< double > tolerances_;
+        std::vector< bool > active_;
+    };
+} // namespace geode