Sim.cpp

Go to the documentation of this file.

/*
 * Copyright 2011-2016 Universiteit Antwerpen
 *
 * Licensed under the EUPL, Version 1.1 or  as soon they will be approved by
 * the European Commission - subsequent versions of the EUPL (the "Licence");
 * You may not use this work except in compliance with the Licence.
 * You may obtain a copy of the Licence at: http://ec.europa.eu/idabc/eupl5
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the Licence is distributed on an "AS IS" basis,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the Licence for the specific language governing
 * permissions and limitations under the Licence.
 */
#include "Sim.h"

#include "CoreData.h"
#include "SimState.h"
#include "TimeSlicer.h"

#include "bio/Cell.h"
#include "bio/CBMBuilder.h"
#include "bio/MBMBuilder.h"
#include "bio/Mesh.h"
#include "bio/MeshCheck.h"
#include "bio/PBMBuilder.h"
#include "bio/Wall.h"
#include "math/RandomEngine.h"
#include "math/RandomEngineType.h"
#include "model/ComponentFactoryProxy.h"
#include "util/clock_man/TimeStamp.h"
#include "util/clock_man/Utils.h"
#include "util/misc/Exception.h"
#include "util/misc/log_debug.h"

#include <boost/optional.hpp>
#include <boost/property_tree/exceptions.hpp>
#include <boost/property_tree/xml_parser.hpp>

#include <algorithm>
#include <functional>
#include <iomanip>
#include <limits>
#include <sstream>
#include <string>
#include <tuple>

using namespace std;
using namespace std::chrono;
using namespace std::placeholders;
using namespace boost::property_tree;
using namespace boost::property_tree::xml_parser;
using boost::optional;
using namespace SimPT_Sim::ClockMan;
using namespace SimPT_Sim::Util;
using namespace SimPT_Sim;

namespace SimPT_Sim {

Sim::Sim()
        : m_is_stationary(false)
{
        m_cd.m_coupled_sim_transfer_data = make_shared<CoreData::TransferMapType>();
        m_cd.m_parameters    = make_shared<ptree>();
        m_cd.m_random_engine = make_shared<RandomEngine>();
        m_cd.m_time_data     = make_shared<TimeData>();
}

Sim& Sim::operator=(const Sim& other)
{
        if (this != &other) {
                //The mutex doesn't have to be copied.

                m_is_stationary   = other.m_is_stationary;
                m_project_name    = other.m_project_name;
                m_run_date        = other.m_run_date;
                m_timings         = other.m_timings;

                m_cd.m_coupled_sim_transfer_data   = make_shared<CoreData::TransferMapType>();
                m_cd.m_mesh                        = CBMBuilder().Build(*other.m_cd.m_mesh);
                m_cd.m_parameters                  = make_shared<ptree>(*other.m_cd.m_parameters);
                m_cd.m_time_data                   = other.m_cd.m_time_data;
                m_cd.m_random_engine               = make_shared<RandomEngine>(*other.m_cd.m_random_engine);
        }
        return *this;
}

Sim::~Sim()
{
}

const ptree& Sim::GetParameters() const
{
        return *m_cd.m_parameters;
}

string Sim::GetProjectName() const
{
        return m_project_name;
}

std::string Sim::GetRunDate() const
{
        return m_run_date;
}

int Sim::GetSimStep() const
{
        return m_cd.m_time_data->m_sim_step;
}

double Sim::GetSimTime() const
{
        return m_cd.m_time_data->m_sim_time;
}

SimState Sim::GetState() const
{
        SimState sim_state;
        sim_state.SetProjectName(m_project_name);                // Project name
        sim_state.SetTimeStep(GetSimStep());                     // Sim steps
        sim_state.SetTime(GetSimTime());                         // Sim time
        sim_state.SetParameters(*m_cd.m_parameters);             // Parameters
        sim_state.SetMeshState(m_cd.m_mesh->GetState());         // Mesh

        ptree re_state;
        RandomEngineType::Info re_info = m_cd.m_random_engine->GetInfo();
        re_state.put("type", re_info.type);                      // Random Engine
        re_state.put("seed", re_info.seed);                      // Random Engine
        re_state.put("state", re_info.state);                    // Random Engine
        sim_state.SetRandomEngineState(re_state);                // Random Engine

        return sim_state;
}

string Sim::GetStatusMessage() const
{
            stringstream ss;
            unsigned int  i = static_cast<unsigned int>(GetSimTime());
            ss << "Step: " << GetSimStep()
                    << " ; Time: " << i << " -- " << Utils::ToColonString(seconds(i))
                    << " ; Cells: " << m_cd.m_mesh->GetCells().size();
            return ss.str();
}

Sim::Timings Sim::GetTimings() const
{
        return m_timings;
}

void Sim::Initialize(const ptree& root_pt)
{
        ptree const& sim_pt = root_pt.get_child("vleaf2");

        // Initialize global info.
        m_project_name     = sim_pt.get<string>("project");
        m_run_date         = TimeStamp().ToString();
        m_cd.m_time_data   = std::make_shared<TimeData>();

        // Initilize CoreData: time info
        m_cd.m_time_data->m_sim_step = sim_pt.get<unsigned int>("sim_step", 0U);
        m_cd.m_time_data->m_sim_time = sim_pt.get<double>("sim_time");

        // Initilize CoreData: parameters
        Reinitialize(sim_pt.get_child("parameters"));

        // Initilize CoreData: build mesh.
        PBMBuilder build_director;
        m_cd.m_mesh = build_director.Build(sim_pt.get_child("mesh"));

        // Initilize CoreData: random engine.
        RandomEngineInitialize(sim_pt);
}

void Sim::Initialize(const SimState& sim_state)
{
        // Initialize global info.
        m_project_name  = sim_state.GetProjectName();
        m_run_date      = TimeStamp().ToString();

        // Initilize CoreData: time info
        m_cd.m_time_data                = make_shared<TimeData>();
        m_cd.m_time_data->m_sim_step    = sim_state.GetTimeStep();
        m_cd.m_time_data->m_sim_time    = sim_state.GetTime();

        // Initilize CoreData: parameters
        Reinitialize(sim_state.GetParameters());

        // Initilize CoreData: build mesh.
        MBMBuilder build_director;
        m_cd.m_mesh = build_director.Build(sim_state.GetMeshState());

        // Initilize CoreData: random engine.
        RandomEngineInitialize(sim_state);
}

bool Sim::IsAtTermination() const
{
        const double          sim_time      = GetSimTime();
        const unsigned int    sim_step      = GetSimStep();
        const unsigned int    cell_count    = m_cd.m_mesh->GetCells().size();
        const bool            is_stationary = IsStationary();

        const double max_sim_time
                = m_cd.m_parameters->get<double>("termination.max_sim_time", numeric_limits<double>::max());
        const unsigned int max_sim_step
                = m_cd.m_parameters->get<unsigned int>("termination.max_sim_steps", numeric_limits<unsigned int>::max());
        const unsigned int max_cell_count
                = m_cd.m_parameters->get<unsigned int>("termination.max_cell_count", numeric_limits<unsigned int>::max());
        const bool stationarity_check
                = m_cd.m_parameters->get<bool>("termination.stationarity_check");

        return (sim_time >= max_sim_time) || (sim_step >= max_sim_step)
                || (cell_count >= max_cell_count) || (stationarity_check && is_stationary);
}

bool Sim::IsStationary() const
{
        return m_is_stationary;
}

bool Sim::RandomEngineInitialize(const ptree& sim_pt)
{
        // A random engine as defined in RandomEngineType.h has 3 properties:
        // 'type', 'seed' and 'state'. The 'type' and 'seed' define the starting
        // condition of the engine. The 'state' defines the state of the engine at
        // some point in time.
        //
        // 'type' and 'seed' are saved in the parameters ptree.
        // 'type', 'seed' and 'state' are saved in the random_engine ptree
        //
        // Initialization logic:
        //   If re_state is empty, initialize a new engine from param's type / seed;
        //   Else
        //     If re_state and parameters.random_engine match on type / seed
        //       initialize from state info in re_state;
        //     Else
        //       initialize a new engine from param's type / seed (i.e: ignore re_state);

        bool status = false;      // Gets to be true is reinitialization succeeds.

        // Random engine info from parameters.random_engine (type & seed):
        const ptree& parameters = sim_pt.get_child("parameters");
        const ptree& re_params_pt = parameters.get_child("random_engine");

        // Optional random engine state info from random_engine_state:
        optional<const ptree&> re_state = sim_pt.get_child_optional("random_engine_state");

        if (!re_state) {
                status = m_cd.m_random_engine->Reinitialize(re_params_pt);
        } else {
                const auto params_type              = re_params_pt.get<string>("type");
                optional<unsigned int> params_seed  = re_params_pt.get_optional<unsigned int>("seed");
                const ptree& re_state_pt            = re_state.get();
                const auto state_type               = re_state_pt.get<string>("type");
                const auto state_seed               = re_state_pt.get<unsigned int>("seed");

                const bool types_match = (params_type == state_type);
                const bool seeds_match = params_seed ? (params_seed.get() == state_seed) : true;

                if (types_match && seeds_match) {
                        status = m_cd.m_random_engine->Reinitialize(re_state_pt);
                } else {
                        status = m_cd.m_random_engine->Reinitialize(re_params_pt);
                }
        }

        return status;
}

bool Sim::RandomEngineInitialize(const SimState& sim_state)
{
        // Logic for this method: see RandomEngineInitialize(const ptree&)

        bool status = false;      // Gets to be true is reinitialization succeeds.

        // Random engine info from parameters.random_engine (type & seed):
        const ptree& parameters    = sim_state.GetParameters();
        const ptree& re_params_pt  = parameters.get_child("random_engine");

        // Optional random engine state info from random_engine_state:
        const ptree& re_state_pt              = sim_state.GetRandomEngineState();
        optional<string> re_state_type        = re_state_pt.get_optional<string>("type");
        optional<unsigned int> re_state_seed  = re_state_pt.get_optional<unsigned int>("seed");
        optional<string> re_state_state       = re_state_pt.get_optional<string>("state");

        if (!re_state_type || !re_state_seed || !re_state_state) {
                status = m_cd.m_random_engine->Reinitialize(re_params_pt);
        } else {
                const auto params_type              = re_params_pt.get<string>("type");
                optional<unsigned int> params_seed  = re_params_pt.get_optional<unsigned int>("seed");
                const bool types_match = (params_type == re_state_type.get());
                const bool seeds_match = params_seed ? (params_seed.get() == re_state_seed.get()) : true;

                if (types_match && seeds_match) {
                        m_cd.m_random_engine->Reinitialize(re_state_pt);
                } else {
                        status = m_cd.m_random_engine->Reinitialize(re_params_pt);
                }
        }

        return status;
}

void Sim::Reinitialize(const ptree& p)
{
        lock_guard<mutex> timestep_guard(m_timestep_mutex);
        *m_cd.m_parameters = p;
}

void Sim::ResetTimings()
{
        m_timings.Clear();
}

void Sim::TimeSliceGo(const std::lock_guard<std::mutex>&, double time_slice, SimPhase phase)
{
        // Pre-condition.
        assert((m_time_evolver != nullptr) && ("TimeEvolver factory produces nullptr."));
        assert((m_cd.Check()) && "CoreData do no check out!");
        assert(MeshCheck(*m_cd.m_mesh).CheckAll() && "Mesh inconsistent!");

        // Effect time step, record timings, advance step counter ...
        const auto tup = m_time_evolver(time_slice, phase);
        m_timings.Merge(get<0>(tup));
        m_is_stationary = get<1>(tup);

        // Post-condition.
        assert(MeshCheck(*m_cd.m_mesh).CheckAll() && "Mesh inconsistent!");
}

void Sim::TimeSliceSetup(const std::lock_guard<std::mutex>&)
{
        // Pre-condition.
        assert((m_cd.Check()) && "CoreData do no check out!");

        // Create and initialize time evolver with current parameters.
        const auto f   = ComponentFactoryProxy::Create(m_cd.m_parameters->get<string>("model.group", ""));
        m_time_evolver = f->CreateTimeEvolver(m_cd);

        // Post condition.
        assert((m_time_evolver != nullptr) && ("Time evolver factory produces nullptr."));
}

void Sim::TimeSliceWrapup(const std::lock_guard<std::mutex>&)
{
        m_time_evolver = TimeEvolverComponent();
}

unique_ptr<TimeSlicer> Sim::TimeSlicing()
{
        return unique_ptr<TimeSlicer>(new TimeSlicer(m_timestep_mutex, this));
}

void Sim::TimeStep()
{
        auto slicer = TimeSlicing();
        const double time_step = m_cd.m_parameters->get<double>("model.time_step");
        slicer->Go(time_step);
        ++m_cd.m_time_data->m_sim_step;
}

ptree Sim::ToPtree() const
{
        ptree sim_pt;
        const string time_stamp = TimeStamp().ToString();
        const string mark = " Generated by SimPT processing at " + time_stamp;
        sim_pt.add("<xmlcomment>", mark);

        sim_pt.put("vleaf2.project",  GetProjectName());
        sim_pt.put("vleaf2.run_date", GetRunDate());
        sim_pt.put("vleaf2.sim_step", GetSimStep());
        sim_pt.put("vleaf2.sim_time", GetSimTime());
        sim_pt.put_child("vleaf2.parameters", GetParameters());
        sim_pt.put_child("vleaf2.mesh", m_cd.m_mesh->ToPtree());

        RandomEngineType::Info info = m_cd.m_random_engine->GetInfo();
        sim_pt.put("vleaf2.random_engine_state.type", info.type);
        sim_pt.put("vleaf2.random_engine_state.seed", info.seed);
        sim_pt.put("vleaf2.random_engine_state.state", info.state);

        return sim_pt;
}

} // namespace