diff --git a/benchmarks/Memory/main.cpp b/benchmarks/Memory/main.cpp new file mode 100644 index 00000000..4f468acd --- /dev/null +++ b/benchmarks/Memory/main.cpp @@ -0,0 +1,151 @@ +// ReSharper disable CppUnusedIncludeDirective +#include +#include +#include +#include +#include + +#include // Required for parallel_setup + +#include "fourdst/logging/logging.h" +#include "fourdst/atomic/species.h" +#include "fourdst/composition/utils.h" + +#include "quill/Logger.h" +#include "quill/Backend.h" +#include "CLI/CLI.hpp" + +#include + +#include "gridfire/gridfire.h" +#include "fourdst/composition/composition.h" +#include "gridfire/utils/gf_omp.h" + +#include +#include +#include + +static std::atomic g_allocated_bytes{0}; + +void* operator new(std::size_t size) { + g_allocated_bytes += size; + if (void* ptr = std::malloc(size)) { + return ptr; + } + throw std::bad_alloc(); +} + +void operator delete(void* ptr, std::size_t size) { + g_allocated_bytes -= size; + std::free(ptr); +} + +void operator delete(void* ptr) { + std::free(ptr); +} + +struct MemoryScopeTracker { + size_t start_bytes; + MemoryScopeTracker() : start_bytes(g_allocated_bytes.load()) {} + size_t bytes_allocated() const { + return g_allocated_bytes.load() - start_bytes; + } + + void reset_tracking() { + start_bytes = 0; + g_allocated_bytes = 0; + } +}; + + +static std::terminate_handler g_previousHandler = nullptr; +void quill_terminate_handler(); + +gridfire::NetIn init(const double temp, const double rho, const double tMax) { + std::setlocale(LC_ALL, ""); + g_previousHandler = std::set_terminate(quill_terminate_handler); + quill::Logger* logger = fourdst::logging::LogManager::getInstance().getLogger("log"); + logger->set_log_level(quill::LogLevel::Info); + + using namespace gridfire; + const std::vector X = {0.7081145999999999, 2.94e-5, 0.276, 0.003, 0.0011, 9.62e-3, 1.62e-3, 5.16e-4}; + const std::vector symbols = {"H-1", "He-3", "He-4", "C-12", "N-14", "O-16", "Ne-20", "Mg-24"}; + + + const fourdst::composition::Composition composition = fourdst::composition::buildCompositionFromMassFractions(symbols, X); + + NetIn netIn; + netIn.composition = composition; + netIn.temperature = temp; + netIn.density = rho; + netIn.energy = 0; + + netIn.tMax = tMax; + netIn.dt0 = 1e-12; + + return netIn; +} + +void quill_terminate_handler() +{ + quill::Backend::stop(); + if (g_previousHandler) + g_previousHandler(); + else + std::abort(); +} + + +int main(int argc, char* argv[]) { + using namespace gridfire; + + double temp = 1.5e7; + double rho = 1.6e2; + double tMax = 3e17; + std::string output_filename = "gf_mem.csv"; + + CLI::App app("GridFire Memory Benchmarks"); + app.add_option("--temperature", temp, "Temperature in degrees")->default_val(std::format("{:5.2E}", temp)); + app.add_option("--density", rho, "Density in Kg")->default_val(std::format("{:5.2E}", rho)); + app.add_option("--tmax", tMax, "Maximum time in seconds")->default_val(std::format("{:5.2E}", tMax)); + app.add_option("--output", output_filename, "Output filename for intermediate results")->default_val("gf_mem.csv"); + + CLI11_PARSE(app, argc, argv); + + const NetIn netIn = init(temp, rho, tMax); + + std::unique_ptr engine; + + int prev_reactions = 0; + int prev_species = 0; + engine = std::make_unique(netIn.composition, 1); + + MemoryScopeTracker tracker; + std::ofstream mem_file(output_filename, std::ios::out); + mem_file << "depth,species,reactions,engine_memory_bytes,solver_memory_bytes\n"; + for (int depth = 1; depth <= 100; depth++) { + tracker.reset_tracking(); + engine = std::make_unique(netIn.composition, depth); + auto blob = engine->constructStateBlob(); + size_t engine_usage = tracker.bytes_allocated(); + + size_t current_num_species = engine->getNetworkSpecies(*blob).size(); + size_t current_num_reactions = engine->getNetworkReactions(*blob).size(); + + if (prev_reactions == current_num_reactions && prev_species == current_num_species) { + std::println("Found end of useful graph traversal at a depth of {}", depth); + break; + } + tracker.reset_tracking(); + const solver::PointSolver localSolver(*engine); + solver::PointSolverContext solverCtx(*blob); + size_t solver_usage = tracker.bytes_allocated(); + + mem_file << std::format("{},{},{},{},{}\n", depth, current_num_species, current_num_reactions, engine_usage, solver_usage); + + prev_reactions = current_num_reactions; + prev_species = current_num_species; + } + mem_file.close(); + std::println("Memory benchmarks results written to {}", output_filename); +} diff --git a/benchmarks/Memory/meson.build b/benchmarks/Memory/meson.build new file mode 100644 index 00000000..2550fad6 --- /dev/null +++ b/benchmarks/Memory/meson.build @@ -0,0 +1,5 @@ +executable( + 'gf_bench_memory', + 'main.cpp', + dependencies: [gridfire_dep, cli11_dep], +) diff --git a/benchmarks/SingleZoneSolver/main.cpp b/benchmarks/SingleZoneSolver/main.cpp index d450a66b..5d10f6df 100644 --- a/benchmarks/SingleZoneSolver/main.cpp +++ b/benchmarks/SingleZoneSolver/main.cpp @@ -98,7 +98,7 @@ int main() { metadata["Compiler"] = "Unknown"; #endif - metadata["Threads"] = omp_get_max_threads(); + // metadata["Threads"] = omp_get_max_threads(); metadata["Runs"] = runs; metadata["Temperature"] = temp; metadata["Density"] = rho; diff --git a/benchmarks/Timing/main.cpp b/benchmarks/Timing/main.cpp new file mode 100644 index 00000000..4584d7ad --- /dev/null +++ b/benchmarks/Timing/main.cpp @@ -0,0 +1,147 @@ +// ReSharper disable CppUnusedIncludeDirective +#include +#include +#include +#include +#include + +#include "fourdst/logging/logging.h" +#include "fourdst/atomic/species.h" +#include "fourdst/composition/utils.h" + +#include "quill/Logger.h" +#include "quill/Backend.h" +#include "CLI/CLI.hpp" + +#include + +#include "gridfire/gridfire.h" +#include "fourdst/composition/composition.h" +#include "gridfire/utils/gf_omp.h" + +#include +#include +#include + +static std::terminate_handler g_previousHandler = nullptr; +void quill_terminate_handler(); + +gridfire::NetIn init(const double temp, const double rho, const double tMax) { + std::setlocale(LC_ALL, ""); + g_previousHandler = std::set_terminate(quill_terminate_handler); + quill::Logger* logger = fourdst::logging::LogManager::getInstance().getLogger("log"); + logger->set_log_level(quill::LogLevel::Info); + + using namespace gridfire; + const std::vector X = {0.7081145999999999, 2.94e-5, 0.276, 0.003, 0.0011, 9.62e-3, 1.62e-3, 5.16e-4}; + const std::vector symbols = {"H-1", "He-3", "He-4", "C-12", "N-14", "O-16", "Ne-20", "Mg-24"}; + + + const fourdst::composition::Composition composition = fourdst::composition::buildCompositionFromMassFractions(symbols, X); + + NetIn netIn; + netIn.composition = composition; + netIn.temperature = temp; + netIn.density = rho; + netIn.energy = 0; + + netIn.tMax = tMax; + netIn.dt0 = 1e-12; + + return netIn; +} + + +void quill_terminate_handler() +{ + quill::Backend::stop(); + if (g_previousHandler) + g_previousHandler(); + else + std::abort(); +} + +int main(int argc, char* argv[]) { + using namespace gridfire; + + double temp = 1.5e7; + double rho = 1.6e2; + double tMax = 3e17; + + std::string output_filename = "gridfire_timings.csv"; + + CLI::App app("GridFire Timeing Benchmarks"); + app.add_option("--temperature", temp, "Temperature in degrees")->default_val(std::format("{:5.2E}", temp)); + app.add_option("--density", rho, "Density in Kg")->default_val(std::format("{:5.2E}", rho)); + app.add_option("--tmax", tMax, "Maximum time in seconds")->default_val(std::format("{:5.2E}", tMax)); + app.add_option("--output", output_filename, "Output filename for intermediate results")->default_val("gridfire_timings.csv"); + + CLI11_PARSE(app, argc, argv); + + const NetIn netIn = init(temp, rho, tMax); + + std::unique_ptr engine; + + + struct TimingInfo { + double depth; + int num_reactions; + int num_species; + double timing; + }; + + std::vector timings; + + int prev_reactions = 0; + int prev_species = 0; + engine = std::make_unique(netIn.composition, 1); + for (int depth = 1; depth <= 100; depth++) { + engine = std::make_unique(netIn.composition, depth); + auto blob = engine->constructStateBlob(); + + TimingInfo info; + info.depth = depth; + info.num_species = engine->getNetworkSpecies(*blob).size(); + info.num_reactions = engine->getNetworkReactions(*blob).size(); + + if (prev_reactions == info.num_reactions && prev_species == info.num_species) { + std::println("Found end of useful graph traversal at a depth of {}", depth); + break; + } + const solver::PointSolver localSolver(*engine); + solver::PointSolverContext solverCtx(*blob); + + + solverCtx.stdout_logging = true; + + try { + auto start = std::chrono::high_resolution_clock::now(); + auto result = localSolver.evaluate(solverCtx, netIn, false, false); + auto end = std::chrono::high_resolution_clock::now(); + + double ns = std::chrono::duration(end - start).count(); + + info.timing = ns; + + prev_reactions = info.num_reactions; + prev_species = info.num_species; + + timings.push_back(info); + } catch (gridfire::exceptions::CVODESolverFailureError& e) { + continue; + } + } + + std::ofstream csvFile(output_filename, std::ios::out); + csvFile << std::format("# Temperature (K): {}", temp); + csvFile << std::format("# Density: {}", rho); + csvFile << std::format("# TMax: {}", tMax); + csvFile << "depth,reactions,species,time\n"; + for (const auto& [depth, numReactions, numSpecies, ns]: timings) { + std::string line = std::format("{},{},{},{}\n", depth, numReactions, numSpecies, ns); + csvFile << line; + } + csvFile.close(); + + std::println("Timeing Benchmarks results written to {}", output_filename); +} diff --git a/benchmarks/Timing/meson.build b/benchmarks/Timing/meson.build new file mode 100644 index 00000000..fd26557b --- /dev/null +++ b/benchmarks/Timing/meson.build @@ -0,0 +1,5 @@ +executable( + 'gf_bench_timeing', + 'main.cpp', + dependencies: [gridfire_dep, cli11_dep], +) diff --git a/benchmarks/meson.build b/benchmarks/meson.build index 52a7e63b..236c79dc 100644 --- a/benchmarks/meson.build +++ b/benchmarks/meson.build @@ -1,3 +1,5 @@ if get_option('build_benchmarks') subdir('SingleZoneSolver') + subdir('Memory') + subdir('Timing') endif \ No newline at end of file diff --git a/src/include/gridfire/io/generative/mesa.h b/src/include/gridfire/io/generative/mesa.h new file mode 100644 index 00000000..e69de29b diff --git a/src/lib/io/generative/mesa.cpp b/src/lib/io/generative/mesa.cpp new file mode 100644 index 00000000..e69de29b diff --git a/tests/graphnet_sandbox/main.cpp b/tests/graphnet_sandbox/main.cpp deleted file mode 100644 index 3de8db86..00000000 --- a/tests/graphnet_sandbox/main.cpp +++ /dev/null @@ -1,249 +0,0 @@ -// ReSharper disable CppUnusedIncludeDirective -#include -#include -#include -#include -#include - -#include // Required for parallel_setup - -#include "fourdst/logging/logging.h" -#include "fourdst/atomic/species.h" -#include "fourdst/composition/utils.h" - -#include "quill/Logger.h" -#include "quill/Backend.h" -#include "CLI/CLI.hpp" - -#include - -#include "gridfire/gridfire.h" -#include "fourdst/composition/composition.h" -#include "gridfire/utils/gf_omp.h" - - -static std::terminate_handler g_previousHandler = nullptr; -static std::vector>>> g_callbackHistory; -static bool s_wrote_abundance_history = false; -void quill_terminate_handler(); - -gridfire::NetIn init(const double temp, const double rho, const double tMax) { - std::setlocale(LC_ALL, ""); - g_previousHandler = std::set_terminate(quill_terminate_handler); - quill::Logger* logger = fourdst::logging::LogManager::getInstance().getLogger("log"); - logger->set_log_level(quill::LogLevel::Info); - - using namespace gridfire; - const std::vector X = {0.7081145999999999, 2.94e-5, 0.276, 0.003, 0.0011, 9.62e-3, 1.62e-3, 5.16e-4}; - const std::vector symbols = {"H-1", "He-3", "He-4", "C-12", "N-14", "O-16", "Ne-20", "Mg-24"}; - - - const fourdst::composition::Composition composition = fourdst::composition::buildCompositionFromMassFractions(symbols, X); - - NetIn netIn; - netIn.composition = composition; - netIn.temperature = temp; - netIn.density = rho; - netIn.energy = 0; - - netIn.tMax = tMax; - netIn.dt0 = 1e-12; - - return netIn; -} - -void log_results(const gridfire::NetOut& netOut, const gridfire::NetIn& netIn) { - std::vector logSpecies = { - fourdst::atomic::H_1, - fourdst::atomic::He_3, - fourdst::atomic::He_4, - fourdst::atomic::C_12, - fourdst::atomic::N_14, - fourdst::atomic::O_16, - fourdst::atomic::Ne_20, - fourdst::atomic::Mg_24 - }; - - std::vector initial; - std::vector final; - std::vector delta; - std::vector fractional; - for (const auto& species : logSpecies) { - double initial_X = netIn.composition.getMassFraction(species); - double final_X = netOut.composition.getMassFraction(species); - double delta_X = final_X - initial_X; - double fractionalChange = (delta_X) / initial_X * 100.0; - - initial.push_back(initial_X); - final.push_back(final_X); - delta.push_back(delta_X); - fractional.push_back(fractionalChange); - } - - initial.push_back(0.0); // Placeholder for energy - final.push_back(netOut.energy); - delta.push_back(netOut.energy); - fractional.push_back(0.0); // Placeholder for energy - - initial.push_back(0.0); - final.push_back(netOut.dEps_dT); - delta.push_back(netOut.dEps_dT); - fractional.push_back(0.0); - - initial.push_back(0.0); - final.push_back(netOut.dEps_dRho); - delta.push_back(netOut.dEps_dRho); - fractional.push_back(0.0); - - initial.push_back(0.0); - final.push_back(netOut.specific_neutrino_energy_loss); - delta.push_back(netOut.specific_neutrino_energy_loss); - fractional.push_back(0.0); - - initial.push_back(0.0); - final.push_back(netOut.specific_neutrino_flux); - delta.push_back(netOut.specific_neutrino_flux); - fractional.push_back(0.0); - - initial.push_back(netIn.composition.getMeanParticleMass()); - final.push_back(netOut.composition.getMeanParticleMass()); - delta.push_back(final.back() - initial.back()); - fractional.push_back((final.back() - initial.back()) / initial.back() * 100.0); - - std::vector rowLabels = [&]() -> std::vector { - std::vector labels; - for (const auto& species : logSpecies) { - labels.emplace_back(species.name()); - } - labels.emplace_back("ε"); - labels.emplace_back("dε/dT"); - labels.emplace_back("dε/dρ"); - labels.emplace_back("Eν"); - labels.emplace_back("Fν"); - labels.emplace_back("<μ>"); - return labels; - }(); - - - gridfire::utils::Column paramCol("Parameter", rowLabels); - gridfire::utils::Column initialCol("Initial", initial); - gridfire::utils::Column finalCol ("Final", final); - gridfire::utils::Column deltaCol ("δ", delta); - gridfire::utils::Column percentCol("% Change", fractional); - - std::vector> columns; - columns.push_back(std::make_unique>(paramCol)); - columns.push_back(std::make_unique>(initialCol)); - columns.push_back(std::make_unique>(finalCol)); - columns.push_back(std::make_unique>(deltaCol)); - columns.push_back(std::make_unique>(percentCol)); - - - gridfire::utils::print_table("Simulation Results", columns); -} - - -void record_abundance_history_callback(const gridfire::solver::PointSolverTimestepContext& ctx) { - s_wrote_abundance_history = true; - const auto& engine = ctx.engine; - // std::unordered_map> abundances; - std::vector Y; - for (const auto& species : engine.getNetworkSpecies(ctx.state_ctx)) { - const size_t sid = engine.getSpeciesIndex(ctx.state_ctx, species); - double y = N_VGetArrayPointer(ctx.state)[sid]; - Y.push_back(y > 0.0 ? y : 0.0); // Regularize tiny negative abundances to zero - } - - fourdst::composition::Composition comp(engine.getNetworkSpecies(ctx.state_ctx), Y); - - - std::unordered_map> abundances; - for (const auto& sp : comp | std::views::keys) { - abundances.emplace(std::string(sp.name()), std::make_pair(sp.mass(), comp.getMolarAbundance(sp))); - } - g_callbackHistory.emplace_back(ctx.t, abundances); -} - - -void save_callback_data(const std::string_view filename) { - std::set unique_species; - for (const auto &abundances: g_callbackHistory | std::views::values) { - for (const auto &species_name: abundances | std::views::keys) { - unique_species.insert(species_name); - } - } - std::ofstream csvFile(filename.data(), std::ios::out); - csvFile << "t,"; - - size_t i = 0; - for (const auto& species_name : unique_species) { - csvFile << species_name; - if (i < unique_species.size() - 1) { - csvFile << ","; - } - i++; - } - - csvFile << "\n"; - - for (const auto& [time, data] : g_callbackHistory) { - csvFile << time << ","; - size_t j = 0; - for (const auto& species_name : unique_species) { - if (!data.contains(species_name)) { - csvFile << "0.0"; - } else { - csvFile << data.at(species_name).second; - } - if (j < unique_species.size() - 1) { - csvFile << ","; - } - ++j; - } - csvFile << "\n"; - } - - csvFile.close(); -} - -void log_callback_data(const double temp) { - if (s_wrote_abundance_history) { - std::cout << "Saving abundance history to abundance_history.csv" << std::endl; - save_callback_data("abundance_history_" + std::to_string(temp) + ".csv"); - } - -} - -void quill_terminate_handler() -{ - log_callback_data(1.5e7); - quill::Backend::stop(); - if (g_previousHandler) - g_previousHandler(); - else - std::abort(); -} - -void callback_main(const gridfire::solver::PointSolverTimestepContext& ctx) { - record_abundance_history_callback(ctx); -} - -int main() { - using namespace gridfire; - - constexpr double temp = 1.5e7; - constexpr double rho = 1.6e2; - constexpr double tMax = 3e17; - - const NetIn netIn = init(temp, rho, tMax); - - - auto engine = engine::GraphEngine(netIn.composition, engine::NetworkBuildDepth::Full); - auto blob = engine.constructStateBlob(); - - const solver::PointSolver localSolver(engine); - solver::PointSolverContext solverCtx(*blob); - - auto result = localSolver.evaluate(solverCtx, netIn, false, false); - std::cout << result << std::endl; -} diff --git a/validation/vv/figures.ipynb b/validation/vv/figures.ipynb new file mode 100644 index 00000000..54f657bb --- /dev/null +++ b/validation/vv/figures.ipynb @@ -0,0 +1,37 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "id": "initial_id", + "metadata": { + "collapsed": true + }, + "outputs": [], + "source": [ + "" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 2 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython2", + "version": "2.7.6" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}