RelapSolver.cpp

Go to the documentation of this file.

 
 
#include "RelapSolver.h"
#include "../Global/fmodule.h"
#include "../Global/CheckedIO.h"
#include "../Global/SimulationContext.h"
#include "../Global/RestartTrim.h"
#include <cstring>
#include <cstdlib>
#include <iostream>
using namespace std;
 
// ---- Fortran entry points ----
extern "C" { void relap5_(); }
extern "C" { void ads_set_cache_mode_(int* mode); }
extern "C" { void ads_set_membuf_mode_(int mode); }
extern "C" { void ads_membuf_get_trimmed_(char** ptr, int* size); }
extern "C" { void ads_membuf_stats_(int* sz, int* dumps,
                                    long long* b_rstrec, long long* b_rnewp,
                                    long long* b_wrplid); }
extern "C" { void ads_membuf_get_(char** ptr, int* sz); }
extern "C" { void ads_membuf_set_(char* ptr, int sz); }
 
// ---- Static callback routing ----
// Set before relap5_() is called, cleared after it returns.
// Safe: RELAP5 is single-threaded; in MPI, each rank is a separate process.
static SolverCallback* s_active_callback = nullptr;
static ThermalHydraulicSolver* s_active_solver = nullptr;
 
SolverCallback* getActiveSolverCallback() { return s_active_callback; }
ThermalHydraulicSolver* getActiveSolver() { return s_active_solver; }
 
// ---- Sequence execution ----
 
void RelapSolver::runSequence(SolverCallback& cb) {
    s_active_callback = &cb;
    s_active_solver = this;
    relap5_();
    s_active_callback = nullptr;
    s_active_solver = nullptr;
}
 
// ---- State accessors (thin wrappers around r5pars_) ----
 
double RelapSolver::currentTime() const { return r5pars_.timesno; }
int RelapSolver::stopFlag() const { return r5pars_.csnostop; }
bool RelapSolver::isTerminateAll() const { return r5pars_.terminateall; }
void RelapSolver::setTerminateAll(bool v) { r5pars_.terminateall = v; }
 
// ---- Branch control ----
 
void RelapSolver::setBranchCount(int n) { r5pars_.cnbranch = n; }
void RelapSolver::setGenNode(bool v) { r5pars_.cgennode = v; }
 
// ---- Restart lifecycle ----
 
void RelapSolver::setNewRestart(bool v) { r5pars_.cnewrst = v; }
 
void RelapSolver::setInputDeck(const std::string& path) {
    strncpy(r5pars_.carg2, path.c_str(), sizeof(r5pars_.carg2) - 1);
    r5pars_.carg2[sizeof(r5pars_.carg2) - 1] = '\0';
    r5pars_.carg2len = static_cast<int>(path.size());
}
 
// ---- Checkpoint/restore ----
 
std::vector<char> RelapSolver::checkpoint() {
    char* buf_ptr = nullptr;
    int buf_sz = 0;
    // Environment variable ADS_MEMBUF_TRIM (read once per process): "1" makes
    // the capture return [file header + plot-id block][last restart dump]
    // instead of the whole accumulated buffer. A child inherits its parent's
    // entire buffer and appends to it, so an untrimmed capture carries every
    // dump in the lineage. The trimmed form is valid only if the buffer is
    // restored at the point it was captured. Experimental, off by default.
    static const bool trim = [] {
        const char* t = std::getenv("ADS_MEMBUF_TRIM");
        return t && *t == '1';
    }();
    if (trim) ads_membuf_get_trimmed_(&buf_ptr, &buf_sz);
    else      ads_membuf_get_(&buf_ptr, &buf_sz);
    // ADS_MEMBUF_TRACE=1: report what a transfer is actually made of. Reported
    // from C++ because RELAP rebinds Fortran unit 6, so writes from inside the
    // membuf module go somewhere unobservable.
    if (const char* t = std::getenv("ADS_MEMBUF_TRACE"); t && *t == '1') {
        int sz = 0, dumps = 0; long long r = 0, n = 0, w = 0;
        ads_membuf_stats_(&sz, &dumps, &r, &n, &w);
        std::cout << "[MEMBUF] capture " << buf_sz / 1048576.0 << " MB"
                  << "  dumps=" << dumps
                  << "  rstrec=" << r / 1048576.0 << " MB"
                  << "  rnewp=" << n / 1048576.0 << " MB"
                  << "  wrplid=" << w / 1048576.0 << " MB"
                  << "  (per dump " << (dumps ? r / 1048576.0 / dumps : 0.0) << " MB)"
                  << std::endl;
    }
    if (buf_sz > 0 && buf_ptr)
        return std::vector<char>(buf_ptr, buf_ptr + buf_sz);
    return {};
}
 
void RelapSolver::restore(const std::vector<char>& data) {
    if (data.empty()) return;
 
    // Splice the run-wide restart base back on if this buffer omits it.
    //
    // With ADS_RESTART_TRIM=1, MPI transfers ship [header][target] only: the
    // [header][base] prefix is byte-identical across lineages, so it is cached
    // per process instead of being sent with every branch. RELAP aborts in
    // fabend if the first restart record it sees is not a COMPLETE one
    // (rrestf.f), so the prefix has to come back before the restore.
    //
    // This lives HERE as well as in Client, because serial restores go
    // straight to solver.restore() from ADSIDAC.cpp, and an MPI checkpoint
    // resumed in serial mode carries base-less buffers.
    if (!adsidac::has_complete_restart(data)) {
        if (base_cache_.empty()) {
            std::vector<char> f;
            const std::string p = g_sim.SScheduler.getProjectDirectory() + "/restart_base.bin";
            if (read_file_verified(p, f) && !f.empty()) base_cache_ = std::move(f);
        }
        if (base_cache_.empty()) {
            std::cerr << "FATAL: restart buffer has no complete restart record and no "
                      << "cached base is available — restoring it would abend in RELAP."
                      << std::endl;
            std::abort();
        }
        auto spliced = adsidac::splice_base(base_cache_, data);
        if (spliced.empty()) {
            std::cerr << "FATAL: could not splice the restart base onto a buffer."
                      << std::endl;
            std::abort();
        }
        ads_membuf_set_(spliced.data(), static_cast<int>(spliced.size()));
        return;
    }
 
    // Self-contained: remember the base so later base-less buffers can be fixed
    // without touching the disk.
    if (base_cache_.empty()) base_cache_ = adsidac::extract_header_and_base(data);
    ads_membuf_set_(const_cast<char*>(data.data()), static_cast<int>(data.size()));
}
 
// ---- Configuration ----
 
void RelapSolver::setThermoCacheMode(int mode) { ads_set_cache_mode_(&mode); }
void RelapSolver::setMembufMode(int mode) { ads_set_membuf_mode_(mode); }