Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Node difference for abscissa in Primal/Dual Integrals #268

Open
wants to merge 2 commits into
base: latest
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion libecole/include/ecole/reward/bound-integral.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@ template <Bound bound> class ECOLE_EXPORT BoundIntegral : public RewardFunction
public:
using BoundFunction = std::function<std::tuple<Reward, Reward>(scip::Model& model)>;

ECOLE_EXPORT BoundIntegral(bool wall_ = false, const BoundFunction& bound_function_ = {});
ECOLE_EXPORT BoundIntegral(bool wall_ = false, bool use_nnodes_ = false, const BoundFunction& bound_function_ = {});
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let's change use_nnodes from a boolean to an enumeration like abscissa:

enum struct Abscissa { wall_time, n_nodes };

And

BoundIntegral(bool wall_ = false, Abscissa abscissa = Abscissa::wall_time, const BoundFunction& bound_function_ = {});


ECOLE_EXPORT void before_reset(scip::Model& model) override;
ECOLE_EXPORT Reward extract(scip::Model& model, bool done = false) override;
Expand All @@ -26,6 +26,7 @@ template <Bound bound> class ECOLE_EXPORT BoundIntegral : public RewardFunction
Reward initial_dual_bound = 0.0;
Reward offset = 0.0;
bool wall = false;
bool use_nnodes;
};

using PrimalIntegral = BoundIntegral<Bound::primal>;
Expand Down
115 changes: 92 additions & 23 deletions libecole/src/reward/bound-integral.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -25,15 +25,23 @@ class IntegralEventHandler : public ::scip::ObjEventhdlr {
inline static auto constexpr base_name = "ecole::reward::IntegralEventHandler";
inline static auto integral_reward_function_counter = 0;

IntegralEventHandler(SCIP* scip, bool wall_, bool extract_primal_, bool extract_dual_, const char* name_) :
IntegralEventHandler(
SCIP* scip,
bool wall_,
bool use_nnodes_,
bool extract_primal_,
bool extract_dual_,
const char* name_) :
ObjEventhdlr(scip, name_, "Event handler for primal and dual integrals"),
wall{wall_},
use_nnodes{use_nnodes_},
extract_primal{extract_primal_},
extract_dual{extract_dual_} {}

~IntegralEventHandler() override = default;

[[nodiscard]] std::vector<std::chrono::nanoseconds> const& get_times() const noexcept { return times; }
[[nodiscard]] std::vector<SCIP_Longint> const& get_nodes() const noexcept { return nodes; }
[[nodiscard]] std::vector<SCIP_Real> const& get_primal_bounds() const noexcept { return primal_bounds; }
[[nodiscard]] std::vector<SCIP_Real> const& get_dual_bounds() const noexcept { return dual_bounds; }

Expand All @@ -50,9 +58,11 @@ class IntegralEventHandler : public ::scip::ObjEventhdlr {

private:
bool wall;
bool use_nnodes;
bool extract_primal;
bool extract_dual;
std::vector<std::chrono::nanoseconds> times;
std::vector<SCIP_Longint> nodes;
Comment on lines 64 to +65
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Then what about having a single vector here like std::vector<double> abscissa instead?
Also related to #198 ?

std::vector<SCIP_Real> primal_bounds;
std::vector<SCIP_Real> dual_bounds;
};
Expand Down Expand Up @@ -127,6 +137,23 @@ auto get_dual_bound(SCIP* scip) {
}
}

/* Get the number of nodes in the branch and bound tree */
auto get_nnodes(SCIP* scip) {
switch (SCIPgetStage(scip)) {
case SCIP_STAGE_TRANSFORMED:
case SCIP_STAGE_INITPRESOLVE:
case SCIP_STAGE_PRESOLVING:
case SCIP_STAGE_EXITPRESOLVE:
case SCIP_STAGE_PRESOLVED:
case SCIP_STAGE_INITSOLVE:
case SCIP_STAGE_SOLVING:
case SCIP_STAGE_SOLVED:
return SCIPgetNNodes(scip);
default:
return static_cast<SCIP_Longint>(0);
}
}

auto time_now(bool wall) -> std::chrono::nanoseconds {
if (wall) {
return std::chrono::steady_clock::now().time_since_epoch();
Expand Down Expand Up @@ -157,7 +184,11 @@ void IntegralEventHandler::extract_metrics(SCIP* scip, SCIP_EVENTTYPE event_type
dual_bounds.push_back(dual_bounds.back());
}
}
times.push_back(time_now(wall));
if (use_nnodes) {
nodes.push_back(get_nnodes(scip));
} else {
times.push_back(time_now(wall));
}
Comment on lines -160 to +191
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let's abstract that in a function/method, something like get_abscissa.
And elsewhere...

}

void IntegralEventHandler::clear_bounds() {
Expand All @@ -171,9 +202,16 @@ void IntegralEventHandler::clear_bounds() {
primal_bounds.clear();
primal_bounds.push_back(last_primal);
}
auto last_time = times.back();
times.clear();
times.push_back(last_time);

if (use_nnodes) {
auto last_node = nodes.back();
nodes.clear();
nodes.push_back(last_node);
} else {
auto last_time = times.back();
times.clear();
times.push_back(last_time);
}
}

/*************************************
Expand All @@ -183,17 +221,24 @@ void IntegralEventHandler::clear_bounds() {
auto compute_dual_integral(
std::vector<SCIP_Real> const& dual_bounds,
std::vector<std::chrono::nanoseconds> const& times,
std::vector<SCIP_Longint> const& nodes,
bool use_nnodes,
SCIP_Real const offset,
SCIP_Real const initial_dual_bound,
SCIP_Objsense obj_sense) {
SCIP_Real dual_integral = 0.0;
for (std::size_t i = 0; i < dual_bounds.size() - 1; ++i) {
auto const time_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
double metric_diff;
if (use_nnodes) {
metric_diff = static_cast<double>(nodes[i + 1] - nodes[i]);
} else {
metric_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
}
auto const dual_bound = dual_bounds[i];
if (obj_sense == SCIP_OBJSENSE_MINIMIZE) {
dual_integral += (offset - std::max(dual_bound, initial_dual_bound)) * time_diff;
dual_integral += (offset - std::max(dual_bound, initial_dual_bound)) * metric_diff;
} else {
dual_integral += -(offset - std::min(dual_bound, initial_dual_bound)) * time_diff;
dual_integral += -(offset - std::min(dual_bound, initial_dual_bound)) * metric_diff;
}
}
return dual_integral;
Expand All @@ -202,17 +247,24 @@ auto compute_dual_integral(
auto compute_primal_integral(
std::vector<SCIP_Real> const& primal_bounds,
std::vector<std::chrono::nanoseconds> const& times,
std::vector<SCIP_Longint> const& nodes,
bool use_nnodes,
SCIP_Real const offset,
SCIP_Real const initial_primal_bound,
SCIP_Objsense obj_sense) {
SCIP_Real primal_integral = 0.0;
for (std::size_t i = 0; i < primal_bounds.size() - 1; ++i) {
auto const time_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
double metric_diff;
if (use_nnodes) {
metric_diff = static_cast<double>(nodes[i + 1] - nodes[i]);
} else {
metric_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
}
auto const primal_bound = primal_bounds[i];
if (obj_sense == SCIP_OBJSENSE_MINIMIZE) {
primal_integral += -(offset - std::min(primal_bound, initial_primal_bound)) * time_diff;
primal_integral += -(offset - std::min(primal_bound, initial_primal_bound)) * metric_diff;
} else {
primal_integral += (offset - std::max(primal_bound, initial_primal_bound)) * time_diff;
primal_integral += (offset - std::max(primal_bound, initial_primal_bound)) * metric_diff;
}
}
return primal_integral;
Expand All @@ -222,21 +274,28 @@ auto compute_primal_dual_integral(
std::vector<SCIP_Real> const& primal_bounds,
std::vector<SCIP_Real> const& dual_bounds,
std::vector<std::chrono::nanoseconds> const& times,
std::vector<SCIP_Longint> const& nodes,
bool use_nnodes,
SCIP_Real const initial_primal_bound,
SCIP_Real const initial_dual_bound,
SCIP_Objsense obj_sense) {
SCIP_Real primal_dual_integral = 0.0;

for (std::size_t i = 0; i < primal_bounds.size() - 1; ++i) {
auto const time_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
double metric_diff;
if (use_nnodes) {
metric_diff = static_cast<double>(nodes[i + 1] - nodes[i]);
} else {
metric_diff = std::chrono::duration<double>(times[i + 1] - times[i]).count();
}
auto const dual_bound = dual_bounds[i];
auto const primal_bound = primal_bounds[i];
if (obj_sense == SCIP_OBJSENSE_MINIMIZE) {
primal_dual_integral +=
-(std::max(dual_bound, initial_dual_bound) - std::min(primal_bound, initial_primal_bound)) * time_diff;
-(std::max(dual_bound, initial_dual_bound) - std::min(primal_bound, initial_primal_bound)) * metric_diff;
} else {
primal_dual_integral +=
(std::min(dual_bound, initial_dual_bound) - std::max(primal_bound, initial_primal_bound)) * time_diff;
(std::min(dual_bound, initial_dual_bound) - std::max(primal_bound, initial_primal_bound)) * metric_diff;
}
}
return primal_dual_integral;
Expand All @@ -252,8 +311,15 @@ auto get_eventhdlr(scip::Model& model, const char* name) -> auto& {
}

/** Add the integral event handler to the model. */
void add_eventhdlr(scip::Model& model, bool wall, bool extract_primal, bool extract_dual, const char* name) {
auto handler = std::make_unique<IntegralEventHandler>(model.get_scip_ptr(), wall, extract_primal, extract_dual, name);
void add_eventhdlr(
scip::Model& model,
bool wall,
bool use_nnodes,
bool extract_primal,
bool extract_dual,
const char* name) {
auto handler =
std::make_unique<IntegralEventHandler>(model.get_scip_ptr(), wall, use_nnodes, extract_primal, extract_dual, name);
scip::call(SCIPincludeObjEventhdlr, model.get_scip_ptr(), handler.get(), true);
// NOLINTNEXTLINE memory ownership is passed to SCIP
handler.release();
Expand Down Expand Up @@ -287,7 +353,8 @@ auto default_primal_dual_bound_function(scip::Model& model) -> std::tuple<Reward
} // namespace

template <Bound bound>
ecole::reward::BoundIntegral<bound>::BoundIntegral(bool wall_, const BoundFunction& bound_function_) : wall{wall_} {
ecole::reward::BoundIntegral<bound>::BoundIntegral(bool wall_, bool use_nnodes_, const BoundFunction& bound_function_) :
wall{wall_}, use_nnodes{use_nnodes_} {
if constexpr (bound == Bound::dual) {
bound_function = bound_function_ ? bound_function_ : default_dual_bound_function;
} else if constexpr (bound == Bound::primal) {
Expand All @@ -306,13 +373,13 @@ template <Bound bound> void BoundIntegral<bound>::before_reset(scip::Model& mode
// Initalize bounds and event handler
if constexpr (bound == Bound::dual) {
std::tie(offset, initial_dual_bound) = bound_function(model);
add_eventhdlr(model, wall, false, true, name.c_str());
add_eventhdlr(model, wall, use_nnodes, false, true, name.c_str());
} else if constexpr (bound == Bound::primal) {
std::tie(offset, initial_primal_bound) = bound_function(model);
add_eventhdlr(model, wall, true, false, name.c_str());
add_eventhdlr(model, wall, use_nnodes, true, false, name.c_str());
} else if constexpr (bound == Bound::primal_dual) {
std::tie(initial_primal_bound, initial_dual_bound) = bound_function(model);
add_eventhdlr(model, wall, true, true, name.c_str());
add_eventhdlr(model, wall, use_nnodes, true, true, name.c_str());
}

// Extract metrics before resetting to get initial reference point
Expand All @@ -327,17 +394,19 @@ template <Bound bound> Reward BoundIntegral<bound>::extract(scip::Model& model,
auto const& dual_bounds = handler.get_dual_bounds();
auto const& primal_bounds = handler.get_primal_bounds();
auto const& times = handler.get_times();
auto const& nodes = handler.get_nodes();
auto const obj_sense = SCIPgetObjsense(model.get_scip_ptr());

// Compute primal integral and difference
SCIP_Real integral = 0.;
if constexpr (bound == Bound::dual) {
integral = compute_dual_integral(dual_bounds, times, offset, initial_dual_bound, obj_sense);
integral = compute_dual_integral(dual_bounds, times, nodes, use_nnodes, offset, initial_dual_bound, obj_sense);
} else if constexpr (bound == Bound::primal) {
integral = compute_primal_integral(primal_bounds, times, offset, initial_primal_bound, obj_sense);
integral =
compute_primal_integral(primal_bounds, times, nodes, use_nnodes, offset, initial_primal_bound, obj_sense);
} else if constexpr (bound == Bound::primal_dual) {
integral = compute_primal_dual_integral(
primal_bounds, dual_bounds, times, initial_primal_bound, initial_dual_bound, obj_sense);
primal_bounds, dual_bounds, times, nodes, use_nnodes, initial_primal_bound, initial_dual_bound, obj_sense);
}

// Reset arrays for storing bounds
Expand Down
18 changes: 15 additions & 3 deletions python/src/ecole/core/reward.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -179,8 +179,9 @@ void bind_submodule(py::module_ const& m) {
it includes time spent in :py:meth:`~ecole.environment.Environment.reset` and time spent waiting on the agent.
)");
dualintegral.def(
py::init<bool, DualIntegral::BoundFunction>(),
py::init<bool, bool, DualIntegral::BoundFunction>(),
py::arg("wall") = false,
py::arg("use_nnodes") = false,
py::arg("bound_function") = DualIntegral::BoundFunction{},

R"(
Expand All @@ -190,6 +191,9 @@ void bind_submodule(py::module_ const& m) {
----------
wall :
If true, the wall time will be used. If False (default), the process time will be used.
use_nnodes :
If true, the integral will be computed with respect to the number of nodes. Otherwise
the integral is computed with respect to time.
bound_function :
A function which takes an ecole model and returns a tuple of an initial dual bound and the offset
to compute the dual bound with respect to. Values should be ordered as (offset, initial_dual_bound).
Expand All @@ -212,8 +216,9 @@ void bind_submodule(py::module_ const& m) {
it includes time spent in :py:meth:`~ecole.environment.Environment.reset` and time spent waiting on the agent.
)");
primalintegral.def(
py::init<bool, PrimalIntegral::BoundFunction>(),
py::init<bool, bool, PrimalIntegral::BoundFunction>(),
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For the enum, we have a binding function that will make it implicitly convertible from strings.

It will need to be moved inside the extension-helper library (and you might need to rebase first), and be renamed in something more sensible. Let me know if you need help with it.

py::arg("wall") = false,
py::arg("use_nnodes") = false,
py::arg("bound_function") = PrimalIntegral::BoundFunction{},
R"(
Create a PrimalIntegral reward function.
Expand All @@ -222,6 +227,9 @@ void bind_submodule(py::module_ const& m) {
----------
wall :
If true, the wall time will be used. If False (default), the process time will be used.
use_nnodes :
If true, the integral will be computed with respect to the number of nodes. Otherwise
the integral is computed with respect to time.
bound_function :
A function which takes an ecole model and returns a tuple of an initial primal bound and the offset
to compute the primal bound with respect to. Values should be ordered as (offset, initial_primal_bound).
Expand All @@ -243,8 +251,9 @@ void bind_submodule(py::module_ const& m) {
it includes time spent in :py:meth:`~ecole.environment.Environment.reset` and time spent waiting on the agent.
)");
primaldualintegral.def(
py::init<bool, PrimalDualIntegral::BoundFunction>(),
py::init<bool, bool, PrimalDualIntegral::BoundFunction>(),
py::arg("wall") = false,
py::arg("use_nnodes") = false,
py::arg("bound_function") = PrimalDualIntegral::BoundFunction{},
R"(
Create a PrimalDualIntegral reward function.
Expand All @@ -253,6 +262,9 @@ void bind_submodule(py::module_ const& m) {
----------
wall :
If true, the wall time will be used. If False (default), the process time will be used.
use_nnodes :
If true, the integral will be computed with respect to the number of nodes. Otherwise
the integral is computed with respect to time.
bound_function :
A function which takes an ecole model and returns a tuple of an initial primal bound and dual bound.
Values should be ordered as (initial_primal_bound, initial_dual_bound). The default function returns
Expand Down
3 changes: 3 additions & 0 deletions python/tests/test_reward.py
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,9 @@ def pytest_generate_tests(metafunc):
ecole.reward.PrimalIntegral(bound_function=lambda x: (0.0, 0.0)),
ecole.reward.DualIntegral(bound_function=lambda x: (0.0, 0.0)),
ecole.reward.PrimalDualIntegral(bound_function=lambda x: (0.0, 0.0)),
ecole.reward.PrimalIntegral(use_nnodes=True, bound_function=lambda x: (0.0, 0.0)),
ecole.reward.DualIntegral(use_nnodes=True, bound_function=lambda x: (0.0, 0.0)),
ecole.reward.PrimalDualIntegral(use_nnodes=True, bound_function=lambda x: (0.0, 0.0)),
)
metafunc.parametrize("reward_function", all_reward_functions)

Expand Down