Metadata-Version: 2.4
Name: FLToptim
Version: 0.7.0
Summary: Fast Legendre-transform convex optimisation for energy-storage operation and spatial economic dispatch
Author-email: Robin Girard <robin.girard@minesparis.psl.eu>
License-Expression: MIT
Project-URL: Homepage, https://git.persee.minesparis.psl.eu/energy-alternatives/FLToptim
Project-URL: Source, https://git.persee.minesparis.psl.eu/energy-alternatives/FLToptim
Keywords: energy,storage,hydro,battery,dynamic programming,Legendre transform,convex optimisation,economic dispatch,parametric simplex,linear programming
Classifier: Development Status :: 4 - Beta
Classifier: Intended Audience :: Science/Research
Classifier: Operating System :: POSIX :: Linux
Classifier: Operating System :: MacOS
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.9
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Programming Language :: Python :: 3.13
Classifier: Programming Language :: C++
Classifier: Topic :: Scientific/Engineering
Requires-Python: >=3.9
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy
Requires-Dist: highspy
Provides-Extra: test
Requires-Dist: pytest; extra == "test"
Requires-Dist: scipy; extra == "test"
Dynamic: license-file

# FLToptim — Fast Legendre Transform Optim

Fast convex optimisation for energy-storage operation and spatial economic dispatch, built around the
**Legendre-transform dynamic program** for convex piecewise-linear (and quadratic) value functions.

## What's inside

- **`cplfunction` / `cpqfunction`** — the convex piecewise-linear/quadratic value-function DP
  (`OptimMargInt`). Each function is stored as a `map<position, slope-increment>`, so the
  inf-convolutions of the Legendre / dynamic-programming recursion stay near-linear. This is the
  original *dynprogstorage* core (Girard, Barbesant, Foucault, Kariniotakis, 2013).
- **`param_simplex`** — a bespoke parametric right-hand-side **dual simplex** that traces the exact
  injection-cost curve φ_n(y) of a small dense network LP in a single pass (one dual pivot per
  breakpoint), carrying its basis across hours.
- **`elec` / `mr_decompose`** — the **spatial-LP ↔ storage-DP decomposition** for perfect-foresight
  annual dispatch: alternate a per-hour spatial network LP with the per-node storage DP, passing each
  node's convex injection-cost *curve* (not a scalar price) so the scheme is curvature-damped and
  converges to the monolithic optimum. Handles reservoirs, batteries (efficiency kink) and coupled
  multi-energy (electricity + hydrogen) networks, with optional parallelism (hour-chunk curve build;
  Jacobi across storages).

## Install / build

```sh
pip install fltoptim          # from PyPI (builds from source on platforms without a wheel)
# or, from a checkout:
pip install -e .              # builds the Cython DP extension + the parametric-simplex shared library
# or
python setup.py build_ext --inplace
```

Building from source needs a C++17 compiler (Linux/macOS); `numpy` and `highspy` are the runtime
dependencies. The reference frontal solves in the decomposition use Gurobi when available (optional).

### DP value-function storage

The linear storage DP uses the sorted-array backend by default (`flat=True`, and
`dp_flat=True` in `elec.decompose` / `mr_decompose`). It runs the same algorithm and
produces schedules bit-for-bit identical to the historical `std::map` backend, while
measured annual workloads are about 4.5--6 times faster in the DP itself. Pass
`flat=False` or `dp_flat=False` to retain the historical backend for comparison or
diagnostics.

The array backend is the recommended choice while intermediate value functions contain
the usual handful to few dozen breakpoints. The map can theoretically become preferable
when functions retain hundreds or thousands of breakpoints and the recursion performs
many insertions in their middle, because array insertion is linear in the number of
breakpoints. No such crossover has been observed in the supported electricity or
multi-resource annual workloads; benchmark both backends before opting out rather than
using a fixed breakpoint threshold.

For the spatial decomposition, `stop_rel_cost=tol` optionally stops the outer
sweeps when two consecutive system costs differ by at most `tol` relatively,
after `min_sweeps` passes. The default remains a fixed number of sweeps because
larger networks can oscillate rather than converge monotonically. A tolerance of
`1e-6` reduced the measured 3-node annual case from eight to four passes, but did
not trigger within 30 passes on the 10-node case.

`reuse_basis=True` is retained as an experimental diagnostic: it carries the
parametric-simplex basis between outer sweeps and removes repeated HiGHS seeds.
Those seeds are already negligible, however, and degenerate optima can send the
outer iteration along a different trajectory. It is therefore disabled by
default and is not currently recommended as a performance option.

For multi-node cases with one storage direction per node, `shared_hourly_bases=True`
builds one optimal basis per hour at the common current dispatch and uses this
atlas to seed every nodal/storage trace. Complete injection-cost curves are still
retraced, so this is a seed/factorisation experiment rather than curve caching.
Annual measurements reduced trace time by roughly 29% at 3 nodes, 26% at 8 nodes
and 38% at 10 nodes. Degenerate hourly optima can nevertheless change the outer
Gauss--Seidel trajectory. The option is therefore off by default.

The implementation compares two atlas candidates at each hour: the basis propagated
from hour `t-1` and the basis of hour `t` from the preceding outer sweep.
`hourly_basis_selection` exposes five policies: `longitudinal`, `certain`
(switch only when the previous-hour candidate is infeasible and the previous-
sweep candidate is feasible), `predict` (fewest violated basic variables, then
least normalised violation), `canonical` (the predictor plus a lexicographic,
candidate-order-independent signature when both bases are feasible optima), and
`min_pivots` (run both reoptimisations). On the
first annual transition, primal admissibility gives a certain zero-pivot choice
for 55% / 26% / 20% of hours at 3 / 8 / 10 nodes when the longitudinal candidate
is already feasible; the reverse certain choice adds 15% / 14% / 16%. Among the
remaining neither-feasible cases, the violation predictor identifies the actual
minimum-pivot candidate 72% / 77% / 80% of the time. `canonical` is the default
selection inside the atlas: paired annual runs gave 0.386 s at
+0.00258% (3 nodes), 7.84 s at +0.06004% (8 nodes), and 11.17 s at +0.33686%
(10 nodes), versus 18.18 s at +0.32719% for the standard 10-node engine. Fewer
spatial pivots still do not mathematically guarantee a better outer Gauss--Seidel
trajectory; the canonical rule makes the observed degeneracy reproducible rather
than perturbing the economic objective.

Keep the atlas disabled for one-node problems, a single curve per hour,
ramp-aware or multi-resource decompositions, or whenever constructing it costs
more than the repeated reoptimisations it replaces. It is most useful once
several nodal storage curves are rebuilt at every outer sweep (the annual 8- and
10-node benchmarks are representative); measure both modes on smaller systems.

The test suite ships with the package:

```sh
pip install "fltoptim[test]"
python -m pytest --pyargs FLToptim.tests -m "not slow" -q
```

## Reference

R. Girard, V. Barbesant, F. Foucault, G. Kariniotakis, *Fast dynamic programming with application to
storage planning*, 2013. Please cite it if you use this software (see `CITATION.cff`).
