# /for-agents # For agents The whole manual is at [`/llms-full.txt`](/llms-full.txt); the index is at [`/llms.txt`](/llms.txt). ## The mental model **A variable is a dimension.** `m.var("x", (P, W))` occupies a block of the model's single column space. A member's column is computed from its multi-index and is not stored. A variable over millions of columns costs its members, not its columns. **An expression is symbolic.** `cost[P, W] * x[P, W]` contains references, not arrays. Writing it allocates nothing. It becomes matrix entries when a constraint is materialised. **A constraint is an array.** It is a `nimblend` array over its free sets crossed with the column space. There is no assembly step: the array is the matrix. **A definition is a model without its data.** `Definition` provides the vocabulary a model is written in, `set`, `param`, `var`, `constraint` and `set_objective`, over symbols declared with no members and no values. `explain()` reports what it declares. `build(data)` binds a copy and returns a `Model`. One definition builds a model for each dataset it is given. **A built model is inspected through three methods.** `explain()` reports what it built. `row(name, **coords)` computes one row: the entries `assemble` writes for it, and no other row. `absent(name)` reports which coordinates were dropped from a constraint and by which rule. **A session keeps the solver open.** `model.session()` assembles once and keeps the solver's model. `diagnose()` queries the solved instance for the conflicting rows of an infeasible model, or for the ray of an unbounded one. `available()` lists the installed adapters. `capabilities(name)` reports what one adapter supports and which capabilities it rejects together. A model with integer columns has no duals: a mixed-integer model's duals are not its relaxation's. **`nimblend` is the layer below.** Its vocabulary is dimensions, labels, entries and alignment, and it contains no optimization term. Import from `nimblend` itself, never from `nimblend.sparse` or another submodule, and never read an array's `.index` or `.data` or a domain's `.codes`. Each of the three has a reader above it: `coordinates()`, `values()`, `positions_of_coordinates()` and `as_coord()`. Do not build an index matrix either: a domain returns the array over its own members through `array(values)` and `identity(into, coord, start)`. ## The public surface | From | Names | | --- | --- | | `nimopt` | `COLUMN`, `ROW`, `Absence`, `Alias`, `Assembled`, `Coefficient`, `Constraint`, `Definition`, `Diagnosis`, `Explanation`, `Expression`, `Model`, `Option`, `Param`, `Piecewise`, `Relation`, `Row`, `Session`, `Set`, `Solution`, `Sum`, `Term`, `Variable`, `available`, `capabilities`, `dumps`, `load`, `loads`, `options`, `product`, `save`, `subset`, `subset_of` | | `nimblend` | `Array`, `DenseArray`, `Domain`, `EntryBuffer`, `SparseArray`, `combined_dims`, `from_long`, `from_dense`, `is_canonical`, `sum_arrays`, `StoredCoord`, `ProductCoord`, `SubsetCoord` | **A coefficient composes.** A coefficient is a parameter read at its sets or an arithmetic combination of such readings: `price[G, T] / eta[G, T]` is a coefficient written before any data exists, read at its sets like a parameter, and evaluated once when the matrix is built. `+`, `-`, `*`, `/` and a power by a number combine coefficients. An expression also contains a constant: `x + 1 <= 5` produces the row `x <= 4`. ## What goes wrong **A chained comparison.** `0 <= expr <= 10` raises `TypeError`. Python evaluates it as two comparisons joined by `and` and keeps only the second. A relation has no truth value, and the chained form raises instead of dropping the first bound. Write each bound as its own constraint. **A sum over a lag.** `Sum(T - 1, ...)` raises: a sum runs over a set's members. Put the lag on the variable reference, `x[T - 1]`. **The built-in `sum` over a set's members.** `sum(x[S, t] for t in members)` returns the correct expression at a cost: it produces one term per member, where `Sum(T, x[S, T])` produces one term and reduces a dimension. The terms concatenate pairwise and each materialises its own block. Building a model that way runs 24 times slower at 25 members and 275 times slower at 400, and the factor grows with the member count. Use the built-in `sum` for a short list of distinct expressions and `Sum` for a set's members. **A right-hand side over the wrong dimensions.** A constraint's right-hand side is a parameter over exactly its free dimensions. The error message gives both. **Reading values from a model that did not solve.** `objective` and `primal` raise where `feasible` is False. They raise at status `unbounded` and `unbounded_or_infeasible` whatever `feasible` reports. `bound` and `gap` are `None` at those two statuses. `dual` raises where `status` is not `"optimal"`. A solve stopped at a limit reports `feasible` True where the solver found a point, with `bound` and `gap` beside it. Read `status` first. **A domain over a definition's sets.** `product((B, T))` needs each set's coordinate, and a declared set has none. In a definition, give `where=`, `over=` and `subset=` as a tuple of its sets or as one of its parameters, whose coefficients are the coordinates. **A row that is not there.** `row()` raises for a coordinate at which the constraint has no row. `absent()` reports which rule dropped it: a coefficient absent inside a sum removes a **term** and keeps the row; a term absent along a **free** dimension removes the **row**. **Reading a MILP's duals.** A model with integer columns has no duals. `dual()` raises; it does not return the relaxation's duals. Read `primal`. **A conflict HiGHS cannot prove.** HiGHS computes its conflict over the linear relaxation. A model infeasible only through its integrality produces no conflict, and `diagnose()` raises. The Gurobi conflict covers the integrality. **Two operands that share no dimension.** Every binary operator combines two dimensioned operands only where they share a dimension. The rule applies to a coefficient multiplied by a variable and to two coefficients alike. Frames sharing no dimension raise `ValueError`; their combination would be an outer product. A number has no dimension and scales every entry. **A division by zero.** A divisor that is zero raises `ZeroDivisionError` with the coordinate, for a Python number, a NumPy scalar and a coefficient with a zero at one coordinate alike. Handle the divisor before it is passed to an expression. **A derived coefficient read at the wrong sets.** A combination is read at its sets as a parameter is, and the reading is checked against the dimensions it has. `unit_cost[T, G]` raises where it is written and reports `('G', 'T')`. **Bypassing the `nimblend` interface.** A test fails on an import from a `nimblend` submodule. It also fails on a read of an array's `.index` or `.data` or a domain's `.codes`, and on a module of the package that assembles an index matrix of its own. ## Every error, and where it is shown The prose above covers the common mistakes. The table lists every error the documentation demonstrates. Each row is executed to produce the message beside it. | Raises | Message | Shown at | | --- | --- | --- | | `ValueError` | the upper bound 'cap' has no value at member ('b',) of variable 'x'; give the bound a value at every member of the variable | [/guides/bounds-from-parameters](/guides/bounds-from-parameters) | | `ValueError` | variable 'x' is declared over ('G',) and is not over ['W']; its upper bound 'cap' is declared over ('W',) | [/guides/bounds-from-parameters](/guides/bounds-from-parameters) | | `TypeError` | a coefficient is a parameter; build one with `Param.from_dense` or `Param.from_long` and read it at its sets. A product of two expressions is not linear. | [/guides/coefficient-arithmetic](/guides/coefficient-arithmetic) | | `ValueError` | coefficient (fuel_price / efficiency) is declared over ('G', 'T'); got ('T', 'G') | [/guides/coefficient-arithmetic](/guides/coefficient-arithmetic) | | `ZeroDivisionError` | divisor holed is zero at 1 coordinate(s), first at {'G': 'base', 'T': 1}; remove the zeros or divide by another parameter | [/guides/coefficient-arithmetic](/guides/coefficient-arithmetic) | | `ValueError` | frames ('G',) and ('T',) share no dimension; pass operands that share a dimension | [/guides/coefficient-arithmetic](/guides/coefficient-arithmetic) | | `ValueError` | constraint 'capacity' has free dimensions ('P',); its condition is over ('W',) | [/guides/conditions](/guides/conditions) | | `ValueError` | constraint 'capacity' is given over= and where= together; pass one of them | [/guides/conditions](/guides/conditions) | | `ValueError` | variable 'x' is read at member 't9' of dimension 'T'; read it at a member that set contains | [/guides/fixed-members](/guides/fixed-members) | | `ValueError` | a lag is a whole number of members; got 1.7 | [/guides/lags](/guides/lags) | | `ValueError` | a sum is over the members of ['T'] and takes the set, not a lag of it; write the lag at the variable's reference | [/guides/lags](/guides/lags) [/reference/expression](/reference/expression) | | `ValueError` | parameter 'rate' is read at a lag ['T']; write the lag at the variable's reference | [/guides/lags](/guides/lags) | | `ValueError` | piecewise 'fuel' has points that are not convex, required by sign '>=' at {'G': 'a'}; use method 'incremental' | [/guides/piecewise](/guides/piecewise) | | `ValueError` | 'max(gen[G, T]) <= 10': the syntax supports one call; write Sum | [/guides/saving-and-loading](/guides/saving-and-loading) | | `ValueError` | member '2030-01-01T00:30' does not convert exactly to datetime64[h] at dimension 'T' of variable 'gen'; write a member in the unit of that dimension | [/guides/saving-and-loading](/guides/saving-and-loading) | | `ValueError` | capital does not fall from base to what follows it; pass a capital cost that falls across the merit order | [/models/expansion](/models/expansion) | | `ValueError` | absence is 'unknown' and the array has no value at 3 of 4 coordinates; pass fill= to to_dense() | [/nimblend/arrays](/nimblend/arrays) | | `ValueError` | frames ('P',) and ('Q',) share no dimension; pass operands that share a dimension | [/nimblend/arrays](/nimblend/arrays) | | `ValueError` | label columns have lengths {'t': 2} and the value column has length 1; pass columns of equal length | [/nimblend/arrays](/nimblend/arrays) | | `ValueError` | 3 member(s) numbered from 4 end at position 6, and dimension 'k' has extent 6; pass a smaller start or a larger coord | [/nimblend/domains](/nimblend/domains) | | `ValueError` | a domain of 3 member(s) requires values of shape (3,); got shape (2,) | [/nimblend/domains](/nimblend/domains) | | `ValueError` | constraint 'supply' has free dimensions ('P',); its right-hand side 'demand' is over ('W',) | [/reference/constraint](/reference/constraint) [/tutorial/constraints](/tutorial/constraints) | | `ValueError` | data does not cover ['S']; add an entry for each | [/reference/definition](/reference/definition) | | `ValueError` | parameter 'S' is already declared as a set; declare another name | [/reference/definition](/reference/definition) | | `TypeError` | a relation already has one bound; compare the expression again in its own constraint | [/reference/expression](/reference/expression) | | `TypeError` | a relation has no truth value; write each bound in its own constraint | [/reference/expression](/reference/expression) [/tutorial/constraints](/tutorial/constraints) | | `TypeError` | an LP has no row for a strict inequality; write `<=` or `>=`, and reduce with `Sum` in place of `min` or `max` | [/reference/expression](/reference/expression) | | `TypeError` | an expression has no absolute value: expressions are linear; bound the expression with two rows, or reduce it with `Sum` over its sets | [/reference/expression](/reference/expression) | | `TypeError` | an expression is reduced over the sets it is summed across; specify them with `Sum(I, J, expression)` | [/reference/expression](/reference/expression) | | `TypeError` | cannot divide by an expression: expressions are linear; declare the reciprocal as a coefficient the variable multiplies | [/reference/expression](/reference/expression) | | `TypeError` | cannot raise an expression to a power: expressions are linear; raise a coefficient to the power and multiply it by a variable | [/reference/expression](/reference/expression) | | `ValueError` | term 'x' already sums over ['T']; sum over each dimension once | [/reference/expression](/reference/expression) | | `ValueError` | constraint 'cap' gives where= a domain with no name; declare its members as a parameter and refer to that parameter | [/reference/files](/reference/files) | | `ValueError` | parameter 'c' is given columns ['value', 'S']; a table lists the dimensions then value: ['S', 'value'] | [/reference/files](/reference/files) | | `ValueError` | variable 'x' contains the unknown key 'bound'; write only 'sets', 'subset', 'lower', 'upper', 'integer' | [/reference/files](/reference/files) | | `ValueError` | constraint 'cap' has no row at {'P': 'p3'}; read `absent('cap')` for the rule that dropped it | [/reference/inspection](/reference/inspection) | | `ValueError` | parameter 'cost': label columns have lengths {'P': 1, 'W': 2} and the value column has length 2; pass columns of equal length | [/reference/param](/reference/param) | | `TypeError` | parameter 'price' is over ('G',) and expresses no coefficient until it is read; read it at its sets as price[G] | [/reference/param](/reference/param) | | `ValueError` | status is 'infeasible' and the solver reports no feasible point; read `status` before reading values | [/reference/solution](/reference/solution) [/tutorial/solving](/tutorial/solving) | | `ValueError` | model 'm' has integer columns and 'highs' reports no duals for it; read primal values only | [/reference/solvers](/reference/solvers) | | `TypeError` | parameter 'supply' is over ('P',) and expresses no coefficient until it is read; read it at its sets as supply[P] | [/tutorial/constraints](/tutorial/constraints) | | `ValueError` | absence is 'unknown' and the array has no value at 1 of 6 coordinates; pass fill= to to_dense() | [/tutorial/reading-the-answer](/tutorial/reading-the-answer) | | `ValueError` | parameter 'cost' is over sets of shape (2, 3); got values of shape (2, 2) | [/tutorial/sets-and-parameters](/tutorial/sets-and-parameters) | --- # / # nimopt `nimopt` is a Python library for building linear and mixed-integer programs. A model is declared symbolically over named index sets, in the form of parameters, variables and constraints. The declaration is expanded into a coefficient matrix at assembly or at solve. Solutions are returned as arrays over the same index sets. A primal value is read by label, not by column position. `nimopt` is built on `nimblend`, a labeled sparse N-dimensional array library. `nimblend` contains no optimization vocabulary and does not import `nimopt`. It is documented in [its own section](/nimblend). ## Design **A variable is a dimension.** A constraint is an array indexed over its free sets crossed with the model's column space, with the coefficients as values. There is no assembly step converting the model into a matrix: the array is the matrix. **Absence is distinct from zero.** An entry is either stored or absent, and every array declares what absence means: `"empty"` for a coordinate that contributes nothing, `"unknown"` for one that was never modeled. A missing result is never counted as zero. Division by an absent value raises an error; it does not return an infinity. **A subset determines the columns.** A variable declared over a subset of a set product has one column per member of the subset and none for the rest. The full product is never materialised, at declaration or after it. **Expressions are symbolic.** An expression contains references to variables and parameters, not their values. `cost[P, W] * x[P, W]` is the same expression over a million routes and over six. The values are read when the matrix is built. **Dropped rows are reported.** A row whose terms have no value at some coordinate is dropped; it is not written incompletely. `absent()` lists every dropped row with the rule that dropped it. `row()` returns one row of the assembled matrix in the form passed to the solver. ## Install ```bash pip install "nimopt[highs]" ``` The extra installs `nimopt` and its dependency `nimblend` from PyPI, with HiGHS as the solver backend. HiGHS is the default solver, and `[highs]` installs it. `[gurobi]` and `[mosek]` add those adapters instead, `[bench]` adds the comparison suite and `[dev]` the test and lint tooling. `available()` reports the solvers whose backend can be imported in the current environment. `capabilities(name)` reports what one adapter supports, whether or not its backend is installed. Development installs come from a checkout. `nimblend` is a dependency and installs first, from its clone. `nimopt` then installs from its own root: ```bash pip install /path/to/nimblend pip install ".[highs]" ``` ## A first model A transport problem: two plants with limited supply ship to three warehouses with fixed demand, and the objective is total shipping cost. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(solution.status, solution.objective) print(solution.primal("x").to_dense()) ``` Output: ```text optimal 135.0 [[20. 0. 10.] [ 0. 15. 5.]] ``` The primal values are returned as a 2 by 3 array over plants and warehouses, in the order the sets declare their members. ## Declaring before the data exists A `Definition` declares the same model without binding data. Its sets and parameters are declared by name and its constraints use the same expression syntax. `explain()` reports the whole declaration before any value is read. ```python import numpy as np from nimopt import Definition, Sum d = Definition("transport") P, W = d.set("P"), d.set("W") cost = d.param("cost", (P, W)) supply = d.param("supply", (P,)) demand = d.param("demand", (W,)) x = d.var("x", (P, W)) d.constraint("supply", Sum(W, x[P, W]) <= supply[P]) d.constraint("demand", Sum(P, x[P, W]) >= demand[W]) d.set_objective(Sum(P, W, cost[P, W] * x[P, W])) print(d.explain()) ``` Output: ```text transport min not built sets P · W parameters cost (P,W) · supply (P) · demand (W) variables x (P×W) [0.0, inf] constraint supply (P) Sum(W, x[P, W]) <= supply[P] constraint demand (W) Sum(P, x[P, W]) >= demand[W] objective min Sum(P, W, cost[P, W] * x[P, W]) ``` A definition is copied before it is bound. One definition builds a model for each dataset it is given, and no build changes the definition. The built model is inspected the same way. `row()` computes one row of the assembled matrix in the form passed to the solver. ```python import numpy as np from nimopt import Definition, Sum d = Definition("transport") P, W = d.set("P"), d.set("W") cost = d.param("cost", (P, W)) supply = d.param("supply", (P,)) demand = d.param("demand", (W,)) x = d.var("x", (P, W)) d.constraint("supply", Sum(W, x[P, W]) <= supply[P]) d.constraint("demand", Sum(P, x[P, W]) >= demand[W]) d.set_objective(Sum(P, W, cost[P, W] * x[P, W])) data = { "P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris", "rome"]), "cost": np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]]), "supply": np.array([30.0, 25.0]), "demand": np.array([20.0, 15.0, 15.0]), } m = d.build(data) print(m) print(m.row("demand", W="paris")) print(m.absent("demand")) ``` Output: ```text Model('transport', 1 variables, 6 columns, 5 rows) demand[W='paris'] row 3 1·x[lisbon,paris] + 1·x[porto,paris] >= 15 demand 3 of 3 rows stated by terms ``` ## A variable over a subset Where a variable spans an arc list instead of a full product, it has one column per arc and the product is never built. A thousand plants each serving three warehouses is three thousand columns, not a million. ```python import numpy as np from nimopt import Model, Set, subset P = Set("P", np.array([f"p{i}" for i in range(1000)])) W = Set("W", np.array([f"w{i}" for i in range(1000)])) served = np.array([f"w{(i * 7 + k) % 1000}" for i in range(1000) for k in range(3)]) arcs = subset((P, W), {"P": np.repeat(P.labels, 3), "W": served}) m = Model("transport") x = m.var("x", (P, W), subset=arcs) print(f"{m.n_columns} columns over a product of {len(P) * len(W)}") ``` Output: ```text 3000 columns over a product of 1000000 ``` ## Features - Sets, aliases, subsets and set products as the index structure of every declaration - Parameters from dense arrays or long-form columns, broadcast where a parameter is narrower than the variable it multiplies - Composable coefficients: a parameter read at its sets, or an arithmetic of parameters written before data exists - Conditions on a sum and on a constraint, lags that drop or wrap at the ends of a set, and members fixed at a label - Per-column bounds from a parameter, and variables declared over a subset of a set product - A `Definition` written before data exists and built against any number of datasets - `explain()` on a definition or a built model, `row()` into the assembled matrix, and `absent()` reporting dropped rows and the rule that dropped each - A `Session` that keeps the solved instance open, and `diagnose()` reporting the conflicting rows of an infeasible model or the ray of an unbounded one - Primals and duals returned over their index sets, with absence distinct from zero - Continuous and integer columns, solved through HiGHS, Gurobi or Mosek behind one adapter contract, with `capabilities()` reporting what each adapter supports and which capabilities it rejects together - One option vocabulary translated into each solver's own option names, with a time limit written the same way for every solver - A model written to and read back from YAML, with its data inline or in a sidecar - A corpus of worked models under `nimopt.models`, each with its formulation, its inputs at any size, and an objective computed by arithmetic instead of by a solver ## Performance Where a variable's columns are a subset of a set product, not materialising the product saves memory and time. On a transport model of 400 000 arcs over a 20 000 000-cell product, `nimopt` builds the matrix in 72.9 MB of resident memory against linopy's 1 682.9 MB. The build takes 292.9 ms against 844.5 ms. Where nothing is sparse, the alignment work costs time and saves nothing. On a fully dense temporally coupled model at 2 111 080 rows, the comparison reverses: linopy builds the matrix three times faster, for seven percent more resident memory. The benchmark suite measures both models. The [benchmark page](/explanation/what-the-numbers-measure) gives each figure, the baseline it is measured against, and what it does not claim. ## Documentation - [Get started](/get-started): installation, and the transport model above solved and read back. - [Vocabulary](/vocabulary): the terms used throughout the documentation. - [Tutorial](/tutorial/sets-and-parameters): the transport model built in six steps, one concept per page. - [Playground](/playground): every example runs in the browser and can be edited. - [For agents](/for-agents): the mental model, the public surface and the failure modes on one page. Every Python example on this site is executed by the test suite and shows the output it produced. ## License MIT. See `LICENSE`. ## Citing The package includes a `CITATION.cff`. Cite it by author, name and version: > Gaete-Morales, Carlos. *nimopt* (version 0.20260923.0). MIT. ## Contributing Issues and patches are welcome once the repositories are published. Until then, the most useful contribution is a model that does not fit: the formulations that are awkward to write determine the next features. --- # /get-started # Get started ## Install ```bash pip install nimopt ``` `nimblend` is installed as a dependency. HiGHS is the default solver. ## A transport model Two plants, Lisbon and Porto, ship to three warehouses, Berlin, Paris and Rome. Plant `p` has supply `s[p]`, warehouse `w` has demand `d[w]`, and one unit shipped on route `(p, w)` costs `c[p, w]`. The decision variable `x[p, w]` is the quantity shipped on each route. ```text minimize Σ_{p,w} c[p,w] · x[p,w] subject to Σ_w x[p,w] ≤ s[p] for each plant p Σ_p x[p,w] ≥ d[w] for each warehouse w x[p,w] ≥ 0 ``` In `nimopt` the sets index every declaration and the parameters contain the data. `m.var` declares the decision variable, `m.constraint` adds each constraint family under a name, and `set_objective` sets the objective function. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(solution.status, solution.objective) ``` Output: ```text optimal 135.0 ``` The solver reports an optimal solution with objective 135. ## Reading the solution `primal("x")` returns the shipments as an array indexed over the sets `x` was declared on. `dual("demand")` returns the dual value of each demand row: the change in the objective per unit increase in that warehouse's demand. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() shipped = solution.primal("x") print(shipped.dims) print(shipped.to_dense()) print(solution.dual("demand").to_dense()) ``` Output: ```text ('P', 'W') [[20. 0. 10.] [ 0. 15. 5.]] [3. 1. 6.] ``` Rows are plants and columns are warehouses. Lisbon ships 20 to Berlin and 10 to Rome; Porto ships 15 to Paris and 5 to Rome. The duals of the demand rows are 3, 1 and 6: the marginal cost of one additional unit at each warehouse. Every code block in this documentation is self-contained. The second block therefore repeats the model. Each block runs in a Python session as it is, or opens in the playground with "Run this example". ## Next - [Vocabulary](/vocabulary) defines the terms used throughout: set, member, frame, row, absence, and others. - The [tutorial](/tutorial/sets-and-parameters) builds this model one concept per page. - The [guides](/guides/subsets) cover sparse networks, time lags, conditions on rows, and bounds from data. - [Explanation](/explanation/a-variable-is-a-dimension) covers the design and its costs. --- # /vocabulary # Vocabulary Terms used throughout the documentation, each defined once. Examples refer to the transport model: plants `P = {lisbon, porto}` ship to warehouses `W = {berlin, paris, rome}`. ## Index sets **Set.** A named index dimension with labels. `Set("P", np.array(["lisbon", "porto"]))` is the set of plants. Parameters, variables and constraints are indexed over sets, and solution values are returned over the same sets. **Member.** One element of a set. `"lisbon"` is a member of `P`. **Label.** The name of a member, a string or a number. Labels are the caller-facing identifiers; integer positions are used internally. **Set product.** The Cartesian product of several sets. `P × W` has six members, `("lisbon", "berlin")`, `("lisbon", "paris")` and so on. Variables and parameters are indexed over set products. **Subset.** An explicit list of members of a set product. `subset((P, W), {"P": ..., "W": ...})` lists the routes that exist. A variable over a subset has a column per listed member and none for the rest. **Alias.** A second name for a set, sharing its labels. It allows a parameter or a constraint to relate a set to itself, such as a flow between two nodes of one node set. **Domain.** A set of coordinates over some dimensions. `product((P, W))` and `subset(...)` return one. `subset=`, `where=` and `over=` take a domain. ## Data and decisions **Parameter.** Data indexed over a set product: one value per member. `cost` is indexed over `(P, W)`; `supply` over `P`. A parameter has no column in the matrix. **Coefficient.** The multiplier of a variable in a row. A parameter indexed at its sets, `cost[P, W]`, is a coefficient, and so is an arithmetic combination of such readings, `price[G, T] / eta[G, T]`. **Variable.** A decision variable. `m.var("x", (P, W))` declares one decision per route. Values are read after a solve with `primal("x")`. **Column.** One decision in the coefficient matrix. Each member of a variable is one column. Column indices are computed from member positions and are never assigned by the caller. **Bound.** The interval a column may take values in. The default lower bound is 0 and the default upper bound is infinity. ## Expressions **Expression.** A linear combination of variables, such as `Sum(W, x[P, W])`. Writing an expression records its structure and computes nothing. Values are read when the matrix is built. **Term.** One component of an expression: one variable, an optional coefficient, and the sets summed over. `cost[P, W] * x[P, W]` is one term. **Frame.** The dimensions an expression is still indexed over, also called its free dimensions. `x[P, W]` has frame `(P, W)`; `Sum(W, x[P, W])` has frame `(P,)`. An empty frame is a scalar. **Sum.** Summation over the members of the named sets. The summed sets are removed from the frame. **Lag.** A reference to the previous or next member of a set. `x[T - 1]` references the previous period. A lag either drops the row with no predecessor or, with `T.cyclic`, wraps to the last member. **Fixed member.** A label in place of a set in a reference, `x[G, "t0"]`. It selects that member and removes the set from the frame. ## Constraints and the matrix **Relation.** An expression compared with `<=`, `>=` or `==` to a right-hand side. `Sum(W, x[P, W]) <= supply[P]` is a relation. It becomes of the model when passed to `m.constraint`. **Constraint.** A relation added to the model under a name. It produces one row per member of its expression's frame. **Row.** One inequality or equality of the coefficient matrix. The supply constraint over two plants produces two rows. **Right-hand side.** The scalar or parameter on the other side of the relation. A scalar applies to every row. A parameter is indexed over exactly the frame of the constraint. Any other frame raises `ValueError`. Each row then has its own value. **Objective.** A scalar expression, one with an empty frame, that the solver minimizes or maximizes. `Sum(P, W, cost[P, W] * x[P, W])` is the total shipping cost. **Sense.** The optimization direction, `"min"` or `"max"`, set once on the `Model`. **Materialise.** Evaluate a parameter or an expression into an array of values. Materialisation runs when the matrix is built, not when the expression is written. **Assemble.** Build the coefficient matrix from every constraint. `solve()` assembles before calling the solver. `assemble()` returns the matrix without solving. **Nonzero.** One stored coefficient of the matrix. `nnz` is the count. ## Solutions **Solution.** The return value of `solve()`: a status, an objective value, and primal and dual values. **Status.** The outcome the solver reported: `optimal`, `infeasible`, `unbounded`, or a limit reached. Values are defined only for `optimal`. **Primal.** The value of a variable in the solution, returned over the sets it was declared on. **Dual.** The dual value of a constraint, also called the shadow price: the change in the objective per unit change in that row's right-hand side. Returned over the constraint's frame. **Absence.** A coordinate at which an array has no value, as distinct from a stored zero. Every array declares the meaning of absence: `"empty"` for a coordinate that contributes nothing, used by parameters, or `"unknown"` for one that was never modeled, used by solutions. **Session.** A solver instance kept open on one assembled model. A caller queries it after the solve, for the conflicting rows of an infeasible model. **Definition.** A model written before its data exists, in the same vocabulary. `build(data)` produces a `Model` for one dataset. --- # /reference/constraint # Constraint ## `Constraint` Returned by `Model.constraint`. Rows over an expression's frame, bounded by a right-hand side. ``` Model.constraint(name, relation, where=None, over=None) ``` | Argument | Meaning | | --- | --- | | `name` | the name `Solution.dual` reads it back by | | `relation` | an expression, a sense and a right-hand side | | `where` | a domain intersecting the rows | | `over` | the rows, given explicitly | | Member | Returns | | --- | --- | | `n_rows` | the number of rows it produces | | `nnz` | the number of coefficients they contain | | `labels_at(positions)` | the labels of the rows at those positions, per free dimension | | `position_of(coords)` | the position of the row at a coordinate, one label per free dimension | | `row_at(position)` | the columns, the coefficients and the lower and upper bound of the row at that position | | `row_of(name)` on the `Assembled` | the position of those rows in the matrix | A row derived from the terms exists where every term has a value and the right-hand side has a value. A coefficient absent inside a sum removes a term and keeps the row. A term absent along a free dimension removes the row: a row missing one of its terms would express a constraint that was not written. `over=` gives the rows explicitly instead, and a term covering some of them contributes where it has values. A condition given with `where=` intersects the row domain, and a row outside the condition is not produced. The expression is symbolic, and the constraint stores the term list and no block. The expression is materialised once to compute its shape and once to write it, and it stores nothing between the two. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) m = Model("transport") x = m.var("x", (P, W)) rows = m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) print(rows.n_rows, rows.nnz) print(m.assemble().row_of("supply")) ``` Output: ```text 2 6 slice(0, 2, None) ``` The right-hand side is a number, applied to every row, or a parameter over exactly the free dimensions of the constraint. A parameter gives each row its own value. A parameter over other dimensions raises `ValueError`, and the message gives both sets of dimensions. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= demand[W]) ``` Raises ValueError: ```text ValueError: constraint 'supply' has free dimensions ('P',); its right-hand side 'demand' is over ('W',) ``` --- # /reference/definition # Definition ## `Definition` ``` Definition(name="definition", sense="min") ``` A definition declares the sets, parameters and variables a model is written from, and its constraints, in the expression syntax a model uses. It contains no data. A set declared here identifies a dimension and has no members, and a parameter identifies a shape and has no values. An expression contains references, not arrays. The free dimensions of an equation and its sense are read from the relation, and neither is declared beside it. `sense` is `"min"` or `"max"`, set once here. | Member | Returns | | --- | --- | | `set(name)` | a declared `Set`, whose members arrive with the data | | `alias(name, base)` | a declared `Alias` over one of this definition's sets | | `param(name, sets)` | a declared `Param`, whose values arrive with the data | | `var(name, sets, subset=None, lower=0.0, upper=inf, integer=False)` | a declared `Variable` | | `constraint(name, relation, where=None, over=None)` | nothing; registers the constraint | | `piecewise(name, x, x_points, y, y_points, sign, method, active=None, relaxed=False, where=None)` | a `Piecewise`; `build` generates its declarations and checks its breakpoints | | `build(data)` | a `Model` over the declarations, bound to `data` | | `explain()` | an `Explanation` of what is declared | | `set_objective(expression)` | nothing; sets the objective | | `sense` | `"min"` or `"max"`, as declared | | `sets`, `aliases`, `parameters`, `variables`, `constraints`, `piecewise_declarations` | the registries, keyed by name | ```python from nimopt import Definition, Sum d = Definition("dispatch", sense="min") snapshot = d.set("snapshot") generator = d.set("generator") p_max = d.param("p_max", (generator,)) load = d.param("load", (snapshot,)) cost = d.param("cost", (generator,)) p = d.var("p", (snapshot, generator), lower=0.0, upper=p_max) d.constraint("balance", Sum(generator, p[snapshot, generator]) == load[snapshot]) d.set_objective(Sum(snapshot, generator, cost[generator] * p[snapshot, generator])) print(list(d.sets), list(d.parameters)) print(d.constraints["balance"][0].expression.frame) print(list(d.variables), d) ``` Output: ```text ['snapshot', 'generator'] ['p_max', 'load', 'cost'] ('snapshot',) ['p'] Definition('dispatch', 1 variables, 1 constraints) ``` ## One namespace for sets and parameters Sets and parameters share one key space. The data a definition is built from is keyed by declared name, and one key identifies one symbol. Declaring a parameter under the name of a set raises `ValueError`. ```python raises=ValueError from nimopt import Definition d = Definition("d") S = d.set("S") d.param("S", (S,)) ``` Raises ValueError: ```text ValueError: parameter 'S' is already declared as a set; declare another name ``` Equations are in no data mapping. A constraint may therefore take the name of the parameter that bounds it. ```python from nimopt import Definition, Sum d = Definition("d") S = d.set("S") supply = d.param("supply", (S,)) one = d.param("one", (S,)) x = d.var("x", (S,)) d.constraint("supply", Sum(S, one[S] * x[S]) <= supply[S]) print(list(d.parameters), list(d.constraints)) ``` Output: ```text ['supply', 'one'] ['supply'] ``` ## An alias in a definition `alias(name, base)` declares a second name for one of the sets of the definition. A model relates a set to itself through an alias. The alias has no data of its own and reads the labels bound to its base set. `build` takes members for the set and none for the alias, and an alias in `data` raises `ValueError`. ```python import numpy as np from nimopt import Definition, Sum d = Definition("network", sense="min") N = d.set("N") NP = d.alias("NP", N) limit = d.param("limit", (N, NP)) flow = d.var("flow", (N, NP), lower=0.0) d.constraint("cap", flow[N, NP] <= limit[N, NP]) d.set_objective(Sum(N, NP, limit[N, NP] * flow[N, NP])) m = d.build({"N": np.array(["a", "b"]), "limit": np.ones((2, 2))}) print(m.n_columns, m.n_rows) ``` Output: ```text 4 4 ``` ## Domains in a definition A `Domain` resolves labels through the coordinate of each set, and a declared set has none. `where=` and `over=` on `constraint`, and `subset=` on `var`, therefore take a tuple of the sets of the definition, meaning their full product. They also take one of its parameters, whose coefficients are the coordinates. Both forms resolve to the same domain. A model and a definition declare a sparse variable or an explicit row domain the same way. ## Building `build(data)` copies the declaration graph, binds the copy, numbers the columns and returns a `Model`. `data` maps the name of a declared set to its members and the name of a declared parameter to its values. The definition is unchanged, and it builds one model per dataset it is given. The values of a parameter are given dense over its product, as an array of one value per cell. They are also given long over its entries, as one mapping of label columns and one value column. The long form gives a parameter with coefficients at some coordinates and none at the rest. A variable declared with `subset=` that parameter takes its members from those coordinates. ```python import numpy as np from nimopt import Definition, Sum d = Definition("transport", sense="min") P, W = d.set("P"), d.set("W") cost = d.param("cost", (P, W)) supply = d.param("supply", (P,)) demand = d.param("demand", (W,)) flow = d.var("flow", (P, W), subset=cost, lower=0.0) d.constraint("supply", Sum(W, cost[P, W] * flow[P, W]) <= supply[P]) d.constraint("demand", Sum(P, cost[P, W] * flow[P, W]) >= demand[W]) d.set_objective(Sum(P, W, cost[P, W] * flow[P, W])) m = d.build( { "P": np.array(["p1", "p2"]), "W": np.array(["w1", "w2"]), "cost": ( {"P": np.array(["p1", "p1", "p2"]), "W": np.array(["w1", "w2", "w1"])}, np.array([1.0, 2.0, 3.0]), ), "supply": np.array([3.0, 3.0]), "demand": np.array([1.0, 1.0]), } ) # three arcs, so three columns rather than the four the product would span print(m.n_columns, m.n_rows) print(m.solve().status) ``` Output: ```text 3 4 optimal ``` Data that omits a declaration, or contains a key the definition never declared, raises `ValueError` before anything is bound. ```python raises=ValueError from nimopt import Definition d = Definition("d") d.set("S") d.build({}) ``` Raises ValueError: ```text ValueError: data does not cover ['S']; add an entry for each ``` --- # /reference/explanation # Explanation ## `Explanation` Returned by `Definition.explain` and `Model.explain`. A frozen record of every declaration and what it built. The shapes it is made of are frozen too. A reader takes a field and parses no rendered text. | Field | Contains | | --- | --- | | `name`, `sense` | the name of the declaration and the direction it optimizes | | `built` | whether counts are facts about data or absent | | `sets` | one `SetShape` per dimension | | `parameters` | one `ParamShape` per parameter | | `variables` | one `VariableShape` per variable | | `constraints` | one `ConstraintShape` per equation | | `piecewise` | one `PiecewiseShape` per piecewise declaration | | `objective` | the objective expression, or `None` | | `columns`, `rows`, `nonzeros` | the model's shape, or `None` | A count is `None` where nothing is bound. It is never zero. A count of zero is a value a caller acts on, and a declaration with no data reports no count. | Shape | Fields | | --- | --- | | `SetShape` | `name`, `size` | | `ParamShape` | `name`, `dims`, `entries` | | `VariableShape` | `name`, `dims`, `members`, `columns`, `lower`, `upper`, `integer` | | `ConstraintShape` | `name`, `free`, `sense`, `rows`, `nonzeros`, `relation` | | `PiecewiseShape` | `name`, `free`, `method`, `sign`, `breakpoints`, `generated` | `VariableShape.members` identifies the parameter a sparse variable took its members from, and is `None` for one over the full product. Columns are absent until data binds. Without that field a sparse declaration and a dense one read identically. `ConstraintShape.free` and `.sense` are read off the relation. Neither is declared beside it. An expression contains references and reports both. ```python from nimopt import Definition, Sum d = Definition("transport", sense="min") P, W = d.set("P"), d.set("W") cost = d.param("cost", (P, W)) supply = d.param("supply", (P,)) flow = d.var("flow", (P, W), subset=cost, lower=0.0) d.constraint("supply", Sum(W, cost[P, W] * flow[P, W]) <= supply[P]) d.set_objective(Sum(P, W, cost[P, W] * flow[P, W])) e = d.explain() print(e.built, e.columns, e.variables[0].members) print(e.constraints[0].free, e.constraints[0].sense) print(e) ``` Output: ```text False None cost ('P',) <= transport min not built sets P · W parameters cost (P,W) · supply (P) variables flow (P×W) over cost [0.0, inf] constraint supply (P) Sum(W, cost[P, W] * flow[P, W]) <= supply[P] objective min Sum(P, W, cost[P, W] * flow[P, W]) ``` `PiecewiseShape.generated` lists the variables and constraints a model generated for the declaration. A definition generates none, and the tuple is empty. --- # /reference/expression # Expressions ## `Term` One variable, an optional coefficient, the dimensions summed over, and a scale factor. A term describes a block of coefficients and contains references, not arrays. Writing it allocates nothing: an expression over a million columns costs the same as one over ten. | Member | Returns | | --- | --- | | `free_dims` | the dimensions it is still indexed over | | `carried_dims` | every dimension it has | | `with_coefficient(coefficient)` | the term, scaled by a parameter | | `summing(dims)` | the term, reduced over those dimensions | | `scaled(by)` | the term, multiplied by a number | | `restricted_to(domain)` | the term, over those members only | A caller builds terms through the operators rather than these members: `cost[P, W] * x[P, W]` gives a coefficient, `Sum` gives the reduction, and `-` gives the scale. ## `Expression` A list of terms and the frame they share. The frame is the union of the terms' free dimensions, ordered by the term that introduces each. A term narrower than the frame is broadcast over it when the expression is materialised. | Member | Returns | | --- | --- | | `terms` | the terms it contains | | `frame` | the dimensions it is indexed over | | `coords` | the coordinates of that frame | | `materialise()` | its coefficients as a `nimblend` array | | `materialise_at(at)` | its coefficients at one coordinate of the frame, as a `nimblend` array over the column dimension | ```python import numpy as np from nimopt import Model, Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) m = Model("transport") x = m.var("x", (P, W)) y = m.var("y", (P, W)) combined = cost[P, W] * x[P, W] - y[P, W] print(combined.frame) print(len(combined.terms)) print(combined.terms[0].free_dims) ``` Output: ```text ('P', 'W') 2 ('P', 'W') ``` ## `Sum` ``` Sum(I, J, ..., expression, where=None) ``` The expression reduced over the named sets. Each set named leaves the frame. `where=` takes a domain and restricts the entries of each term before the reduction. The sum then runs over the coordinates given, not over every coordinate of the product. ```python import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) print(x[P, W].frame) print(Sum(W, x[P, W]).frame) print(Sum(P, W, x[P, W]).frame) ``` Output: ```text ('P', 'W') ('P',) () ``` A sum is over the members of a set and takes the set, not a lag of it. Write the lag at the variable reference. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (T,)) Sum(T - 1, x[T]) ``` Raises ValueError: ```text ValueError: a sum is over the members of ['T'] and takes the set, not a lag of it; write the lag at the variable's reference ``` ## `Relation` An expression, a sense and a right-hand side, produced by comparing an expression with `<=`, `>=` or `==`. `Model.constraint` turns one into a constraint. ```python import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) bounded = Sum(W, x[P, W]) <= 30.0 print(type(bounded).__name__, bounded.sense) ``` Output: ```text Relation <= ``` A relation has no truth value. Python evaluates `0 <= expr <= 10` as two comparisons joined by `and` and keeps only the second. The chained form raises, and the first bound is not dropped. ```python raises=TypeError import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) 0.0 <= Sum(W, x[P, W]) <= 10.0 ``` Raises TypeError: ```text TypeError: a relation has no truth value; write each bound in its own constraint ``` ## Forms that are not linear `nimopt` expresses linear terms. Each form below raises where it is written, and the message gives the form to write instead. A variable raised to a power is not linear; a coefficient takes the power and a variable multiplies it. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) x[T] ** 2 ``` Raises TypeError: ```text TypeError: cannot raise an expression to a power: expressions are linear; raise a coefficient to the power and multiply it by a variable ``` A variable in a denominator is not linear either; the reciprocal is written as a coefficient the variable multiplies. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) 1.0 / x[T] ``` Raises TypeError: ```text TypeError: cannot divide by an expression: expressions are linear; declare the reciprocal as a coefficient the variable multiplies ``` The absolute value of an expression is not linear. A magnitude is written with two rows bounding the expression, and a reduction with `Sum` over its sets. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) abs(x[T]) ``` Raises TypeError: ```text TypeError: an expression has no absolute value: expressions are linear; bound the expression with two rows, or reduce it with `Sum` over its sets ``` An LP has no row for a strict inequality. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) x[T] < 5.0 ``` Raises TypeError: ```text TypeError: an LP has no row for a strict inequality; write `<=` or `>=`, and reduce with `Sum` in place of `min` or `max` ``` The built-in `sum` of expressions with no set to reduce over calls `Expression.sum`. That call raises, and it returns no unreduced expression. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) x[T].sum() ``` Raises TypeError: ```text TypeError: an expression is reduced over the sets it is summed across; specify them with `Sum(I, J, expression)` ``` A relation expresses one bound. Comparing it a second time raises, and the first bound is not dropped. ```python raises=TypeError import numpy as np from nimopt import Model, Set T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) (x[T] <= 5.0) >= 1.0 ``` Raises TypeError: ```text TypeError: a relation already has one bound; compare the expression again in its own constraint ``` A dimension is reduced once. A second reduction over the same dimension raises. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum T = Set("T", np.arange(3)) m = Model("m") x = m.var("x", (T,)) Sum(T, Sum(T, x[T])) ``` Raises ValueError: ```text ValueError: term 'x' already sums over ['T']; sum over each dimension once ``` --- # /reference/files # Files ## `dumps`, `load`, `loads`, `save` | Function | Does | | --- | --- | | `dumps(what, inline=False, instructions=False, version=4)` | returns the text of a definition's file, or of a model's file; a model's data is included with `inline=True` alone; a definition with `inline=True` raises `ValueError` | | `load(path, data=None)` | reads a file; returns a `Definition`, or a `Model` where the file contains data or `data=` gives it | | `loads(text, data=None)` | the same over text; a sidecar name in text raises `ValueError`. Text has no directory | | `save(what, path, inline=False, instructions=False, version=4)` | writes a definition's file, or a model's with an `.npz` beside it, or one file with an inline block when `inline=True`; `instructions=True` writes the comment block that describes the format at the top of the file; `version` is `4` or `3` | `data=` is the mapping `build` takes or the path of an `.npz`. A file that contains data and a `data=` together raises `ValueError`. `dumps` returns the text `save` writes, without a sidecar line. Only `save` writes a sidecar and the line that identifies it. `version=4` writes each piecewise declaration under `piecewise`, and omits the sets, parameters, variables and constraints it generated. `version=3` writes a model's generated declarations as ordinary declarations, and omits the piecewise declarations and any parameter only they read. A model loaded from that file contains the same rows and no piecewise declaration. A definition with a piecewise declaration raises `ValueError` at version 3. Any other version raises `ValueError`. With `instructions=True`, every writer puts a fixed comment block at the top of the text. The block describes the format: the keys and their order, the defaults, the rules that determine which rows a constraint has, and the expression syntax. The text is the same in every file and describes the format, not the model. A reader with one file interprets it without this package. The block is a YAML comment: a file with it and a file without it load to the same model. ## The file | Key | Contains | | --- | --- | | `version` | `4`, or `3` where the caller asks for it; files of version `2` and `3` also load; any other value raises `ValueError` and reports the versions this reader accepts | | `name`, `sense` | the model's | | `sets` | a list of names | | `aliases` | each alias to its base set; absent where the model declares none | | `parameters` | each name to its dimensions | | `variables` | each name to `sets`, and to `subset`, `lower`, `upper`, `integer` where they differ from no subset, `0`, infinity and `false` | | `constraints` | each name to `relation`, and to `where` or `over` where given | | `piecewise` | each name to `x`, `x_points`, `y`, `y_points`, `sign`, `method`, and to `active`, `relaxed` and `where` where given; version 4 only; absent where the model declares none | | `objective` | the objective expression; absent where the model declares none | | `data` | an inline mapping, or the name of an `.npz` beside the file | `subset`, `where` and `over` take a parameter's name, meaning the coordinates it contains, or a list of set names, meaning their full product. A symbol's name is a Python identifier other than `Sum`. A dimension listed in `parameters`, in a variable's `sets`, or in a list of set names is a set or an alias. An alias is declared after its base set. The expressions are written as they are typed in Python and are read back through the same operators. The file is a fixed point: reading it and writing it again gives the same text. ```python from nimopt import Definition, Sum, dumps, loads d = Definition("d") S = d.set("S") c = d.param("c", (S,)) x = d.var("x", (S,), integer=True) d.constraint("cap", 2 * c[S] * x[S] - 1 <= 5) text = dumps(d) print(text) print(dumps(loads(text)) == text) ``` Output: ```text version: 4 name: d sense: min sets: [S] parameters: c: [S] variables: x: sets: [S] integer: true constraints: cap: relation: (c[S] * 2) * x[S] - 1 <= 5 True ``` A key the format does not define raises `ValueError`, at the top level and inside an entry. ```python raises=ValueError from nimopt import loads loads( "version: 3\nname: d\nsense: min\nsets: [S]\n" "variables:\n x: {sets: [S], bound: 1}\n" ) ``` Raises ValueError: ```text ValueError: variable 'x' contains the unknown key 'bound'; write only 'sets', 'subset', 'lower', 'upper', 'integer' ``` ## Data The data in a file is the mapping `build` takes, in three shapes. | Shape | Inline | In the `.npz` | | --- | --- | --- | | a set's members | a list | a one-dimensional label array | | a set of `datetime64` or `timedelta64` members | `dtype` and `members` | a one-dimensional label array | | a dense parameter | nested lists in row-major order | its grid | | a long parameter | `columns`, the dimensions then `value`, and `rows` | a structured array with one field per dimension and `value` | A parameter is written dense where its array covers its full product and long otherwise. The `.npz` is read with `allow_pickle=False`. An array of object dtype raises `ValueError` at save and reports the symbol: the container would pickle it, and the reader rejects a pickled array. ```python raises=ValueError from nimopt import loads loads( "version: 3\nname: d\nsense: min\nsets: [S]\nparameters:\n c: [S]\n" "data:\n S: [a]\n c:\n columns: [value, S]\n rows:\n - [1.0, a]\n" ) ``` Raises ValueError: ```text ValueError: parameter 'c' is given columns ['value', 'S']; a table lists the dimensions then value: ['S', 'value'] ``` ## Datetime members A set whose members are `datetime64` or `timedelta64` is written inline as a mapping of `dtype` and `members`, instead of a list. ```yaml data: T: dtype: datetime64[s] members: ['2030-01-01T00:00:00', '2030-01-01T01:00:00'] ``` A `datetime64` member is written as its ISO 8601 string. A `timedelta64` member is written as its integer count of the unit in the `dtype`. A label column of a long table is written in the same text and takes no marker: the column belongs to a set, and the reader converts it to that set's dtype. The `.npz` stores the dtype of every array and uses no separate form. A member fixed in a relation is written as text. A `datetime64` member is its quoted ISO 8601 string, such as `x['2030-01-01T00:00:00']`. A `timedelta64` member is a quoted count and numpy unit code, such as `x['3 h']`. Every member is converted to the dtype of its set. A string is parsed as ISO 8601. A `datetime.datetime`, a `datetime.date` and a `datetime64` of another unit are converted. An integer against a `timedelta64` set is a count of that set's own unit. A conversion that is not exact raises `ValueError`: `'2030-01-01T00:30'` against a set in hours raises instead of truncating to the hour. A member specifies no time zone. `datetime64` represents no offset, and a conversion to UTC would move the member. A string with an offset or a trailing `Z` raises `ValueError`, and a `datetime.datetime` with a `tzinfo` raises the same error. Write the naive form, `'2030-01-01T00:00:00'`. A member outside the range its dtype represents raises `ValueError` and reports the first and the last member of that range. A `NaT` member raises `ValueError`. ## What raises before anything is written | Written | Reason | | --- | --- | | a model whose `subset`, `where` or `over` is a domain with no name | declare the members as a parameter | | two parameter objects or two set objects sharing a name in one model | the file keys a symbol by name | | a symbol whose name is not an identifier, or is `Sum` | the expression syntax cannot address it | | an array of object dtype | the container would pickle it | ```python raises=ValueError import numpy as np from nimopt import Model, Set, dumps, subset P = Set("P", np.array(["a", "b"])) m = Model("m") x = m.var("x", (P,)) m.constraint("cap", x[P] <= 1.0, where=subset((P,), {"P": np.array(["a"])})) dumps(m) ``` Raises ValueError: ```text ValueError: constraint 'cap' gives where= a domain with no name; declare its members as a parameter and refer to that parameter ``` --- # /reference/inspection # Inspecting a built model ## `Row` Returned by `Model.row(name, **coords)`. One row as the assembled matrix stores it: the columns, the coefficients and the bounds that are passed to the solver. `Model.row` computes the entries of that row alone. It builds no matrix and no other row. `coords` gives one label per free dimension. A constraint with no free dimension has one row, and `Model.row(name)` returns it. | Field | Contains | | --- | --- | | `constraint` | the equation this row belongs to | | `coordinate` | the row's own coordinate, per free dimension | | `index` | the solver's own row number | | `terms` | one `RowTerm` per coefficient | | `sense`, `lower`, `upper` | read from the row's bounds | | `RowTerm` field | Contains | | --- | --- | | `column` | the solver's own column number | | `variable` | the variable that column belongs to | | `coordinate` | that column's coordinate, per dimension | | `coefficient` | the value in the matrix | A variable occupies a contiguous range of the column space from its `start`. A column resolves to its variable through that range, and to a coordinate through the numbering rule of that variable. `sense` is read from the bounds: equal bounds are `==`, an infinite lower bound is `<=`, an infinite upper bound is `>=`. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, cost[P, W] * x[P, W]) <= supply[P]) print(m.row("supply", P="porto")) ``` Output: ```text supply[P='porto'] row 1 3·x[porto,berlin] + 1·x[porto,paris] + 6·x[porto,rome] <= 25 ``` A coordinate at which the constraint has no row raises `ValueError`; the message points to the function that reports why it is missing. `row` and `absent` raise `KeyError` for a name that is not a declared constraint, and the message lists the declared constraints. `row` raises `ValueError` for a label that is not a member of its dimension, and the message identifies the dimension. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set P = Set("P", np.array(["p1", "p2", "p3"])) m = Model("m") x = m.var("x", (P,), upper=5.0) one = Param.from_dense("one", (P,), np.ones(3)) rhs = Param.from_long("rhs", (P,), {"P": np.array(["p1", "p2"])}, np.ones(2)) m.constraint("cap", one[P] * x[P] <= rhs[P]) m.row("cap", P="p3") ``` Raises ValueError: ```text ValueError: constraint 'cap' has no row at {'P': 'p3'}; read `absent('cap')` for the rule that dropped it ``` ## `Absence` Returned by `Model.absent(name)`. What a constraint set out to produce, what it produced, and which coordinates were dropped. | Field | Contains | | --- | --- | | `constraint` | the equation this is about | | `stated_by` | `"terms"` where the rows are derived, `"over"` where given explicitly | | `expected`, `standing` | rows expected, rows kept | | `dropped_rows` | one `DroppedRow(coordinate, rule, detail)` per row lost | | `dropped_terms` | one `DroppedTerm(coordinate, variable, rule, detail)` per term lost | `expected - len(dropped_rows) == standing`. | `dropped_rows` rule | Meaning | | --- | --- | | `term-does-not-reach` | a term has no value at that coordinate; the row would express a constraint that was not written | | `where` | the condition excludes it | | `absent-rhs` | the right-hand side has no value there | | `dropped_terms` rule | Meaning | | --- | --- | | `absent-coefficient` | a coefficient absent inside a sum; the row is kept with one term fewer | The two rules differ in what they remove. A coefficient absent inside a sum removes a **term** and keeps the row. A term absent along a **free** dimension removes the **row**. Under `over=` the rows are given explicitly. Nothing is dropped, and a right-hand side that omits one raises instead. An empty `dropped_rows` beside `stated_by="over"` follows from that rule. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["p1", "p2"])) W = Set("W", np.array(["w1", "w2", "w3"])) m = Model("t") flow = m.var("flow", (P, W)) cost = Param.from_long( "cost", (P, W), {"P": np.array(["p1", "p1", "p2"]), "W": np.array(["w1", "w2", "w1"])}, np.array([1.0, 2.0, 3.0]), ) supply = Param.from_dense("supply", (P,), np.array([3.0, 3.0])) m.constraint("supply", Sum(W, cost[P, W] * flow[P, W]) <= supply[P]) print(m.absent("supply")) ``` Output: ```text supply 2 of 2 rows stated by terms term absent P='p1', W='w3' flow absent-coefficient (cost) term absent P='p2', W='w2' flow absent-coefficient (cost) term absent P='p2', W='w3' flow absent-coefficient (cost) ``` --- # /reference/model # Model ## `Model` ``` Model(name="model", sense="min") ``` A model contains one column space, the constraints declared against it, and an objective. `name` labels it and is otherwise unused. `sense` is `"min"` or `"max"`, set once here. Any other value raises `ValueError`. | Member | Returns | | --- | --- | | `var(name, sets, subset=None, lower=0.0, upper=inf, integer=False)` | a `Variable` occupying the next range of columns | | `constraint(name, relation, where=None, over=None)` | a `Constraint` occupying the next range of rows | | `piecewise(name, x, x_points, y, y_points, sign, method, active=None, relaxed=False, where=None)` | a `Piecewise`; declares the variables and constraints of its method | | `set_objective(expression)` | nothing; sets the objective | | `sense` | `"min"` or `"max"`, as declared | | `solve(solver="highs", options=None)` | a `Solution` | | `assemble()` | an `Assembled`: the matrix, with no solver involved | | `n_columns`, `n_rows`, `nnz` | the shape declared so far | | `column_bounds()` | the lower and upper bound vectors, in column order | | `integrality()` | one flag per column | | `objective_coefficients()` | one coefficient per column | | `explain()` | an `Explanation` of what the model built | | `piecewise_declarations` | the piecewise declarations, keyed by name | | `objective` | the objective expression, or `None` | Declaring costs shapes, not blocks: `n_rows` and `nnz` are known when a constraint is added, and no matrix exists until `assemble` or `solve`. ```python import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= 30.0) m.constraint("total", Sum(P, W, x[P, W]) <= 100.0) m.set_objective(Sum(P, W, x[P, W])) print(m.n_columns, m.n_rows, m.nnz) print(m.objective_coefficients()) ``` Output: ```text 6 3 12 [1. 1. 1. 1. 1. 1.] ``` ## `Piecewise` ``` Model.piecewise(name, x, x_points, y, y_points, sign, method, active=None, relaxed=False, where=None) Definition.piecewise(name, x, x_points, y, y_points, sign, method, active=None, relaxed=False, where=None) ``` A piecewise-linear relation of the expression `y` to the expression `x`. `x` is on the curve through `x_points` and `y_points`. `sign` compares `y` with the curve: `"=="`, `"<="` or `">="`. The two points are parameters read at their sets. Each is over some or all of the sets of `x` and over one breakpoint set, the one set `x` is not over. An entity lists its first breakpoints, and its last breakpoints may be absent. An entity with no breakpoint has no generated rows and no generated columns. `where` restricts the declaration to some entities: a parameter, a tuple of sets or a domain over the sets of `x_points` other than the breakpoint set. The breakpoint checks, the generated columns and the generated rows cover the entities at its coordinates. `x`, `y` and `active` are compared at those coordinates only. | `method` | Generates | Requires | | --- | --- | --- | | `"incremental"` | per segment, one continuous and one integer column and their rows | breakpoints strictly increasing or strictly decreasing | | `"tangent"` | one row per segment, and two rows that keep `x` between the first and the last breakpoint | points convex under `>=`, concave under `<=`; no `active`; no `==`; no constant in `x` | | `"auto"` | the declarations of `"tangent"` where its requirements hold at every entity, and of `"incremental"` otherwise | breakpoints strictly increasing or strictly decreasing | `formulation` reports the method a model generates. For `"auto"` it is `None` on a definition. A definition reserves the generated names of both methods for an `"auto"` declaration. `active` is a binary variable over the sets of `x`, or a sum of them. Where it is 0, `x` is 0 and `y` is compared with 0. A term that is scaled or bounded outside 0 and 1 raises ValueError, and so does a continuous term under the default. With `where`, the bounds are read at its coordinates. `relaxed=True` accepts a continuous `active` between 0 and 1 and scales the curve by its value, which is the linear relaxation of the switch. `relaxed=True` with no `active` raises ValueError. `Model.piecewise` generates the declarations at once. `Definition.piecewise` stores the declaration, and `build` generates them. A generated name is `name`, an underscore and a suffix: | `method` | Sets | Parameters | Variables | Constraints | | --- | --- | --- | --- | --- | | `"incremental"` | `segment` | `members`, `x_step`, `y_step`, `x_first`, `y_first` | `fill`, `order` | `x`, `y`, `order_bound`, `fill_order`, `order_link`, `active` | | `"tangent"` | `segment` | `slope`, `intercept`, `x_low`, `x_high` | none | `tangent`, `x_min`, `x_max` | `{name}_active` exists only where `active` is given. The members of `{name}_segment` are the breakpoint set's members without the first. A segment is identified by its end breakpoint. | Member | Contains | | --- | --- | | `name`, `x`, `x_points`, `y`, `y_points`, `sign`, `method`, `active`, `relaxed`, `where` | the arguments | | `breakpoints` | the name of the breakpoint set | | `formulation` | `"incremental"` or `"tangent"`: the method a model generates; `None` for `"auto"` before a model generates it | | `entity` | the dimensions of `x_points` other than the breakpoint set | | `where_domain()` | the domain of `where` over `entity`, or `None` | | `names()` | the names the declaration generates, keyed by `"sets"`, `"parameters"`, `"variables"` and `"constraints"`; both methods' names for an `"auto"` declaration before a model generates it | | `generated` | the names a model generated, keyed the same way; empty on a definition | | `generated_names()` | every generated name, as a frozenset | An argument error raises when the declaration is made. `TypeError` is raised for an `x`, `y` or `active` that is not an expression, and for points that are not a parameter read at its sets. `ValueError` is raised for a name that is not a Python identifier, an unknown `method` or `sign`, expressions over different sets, points without exactly one breakpoint set, points over different sets, `"tangent"` with `"=="`, with `active` or with a constant in `x`, an `active` with a constant, a `where` of another type or over other sets than the entity sets of `x_points`, and a generated name the model or definition declares. A breakpoint error raises `ValueError` when the data is bound, before any declaration, and identifies the first entity at fault: points with no breakpoint, points present at different breakpoints, an entity with one breakpoint, an absent breakpoint before a present one, a value that is not finite, breakpoints that are not strictly monotonic, and, for `"tangent"`, points whose curvature does not match `sign`. ## `Assembled` The model's matrix in CSR form, returned by `assemble`. `indices` and `values` are views of the one buffer the model allocated; only `indptr` is built. | Member | Returns | | --- | --- | | `indptr`, `indices`, `values` | the matrix in CSR form | | `n_rows`, `n_cols` | its shape | | `row_lower`, `row_upper` | one bound per row | | `col_lower`, `col_upper`, `col_cost`, `integrality` | one entry per column | | `row_of(name)` | a constraint's rows, as a slice | | `to_dense()` | the matrix as an ndarray | ```python import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= 30.0) assembled = m.assemble() print(assembled.n_rows, assembled.n_cols) print(assembled.indptr) print(assembled.row_of("supply")) print(assembled.to_dense()) ``` Output: ```text 2 6 [0 3 6] slice(0, 2, None) [[1. 1. 1. 0. 0. 0.] [0. 0. 0. 1. 1. 1.]] ``` `to_dense` is for a small model. A model of any size is read through `row_of` and the CSR arrays. ## What a model built `explain()` reports every declaration with the count it built, and has `built=True`. It returns the record type a `Definition` returns with every count absent, and one reader covers both. A model contains variables and constraints. Its sets and parameters are collected from them, in order of first appearance. A dimension introduced by a coefficient belongs to no variable and is found through the parameter that has it. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, cost[P, W] * x[P, W]) <= supply[P]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) print(m.explain()) ``` Output: ```text transport min 6 columns · 2 rows · 6 nonzeros sets P 2 · W 3 parameters cost (P,W) 6 · supply (P) 2 variables x (P×W) 6 cols [0.0, inf] constraint supply (P) Sum(W, cost[P, W] * x[P, W]) <= supply[P] 2 rows 6 nz objective min Sum(P, W, cost[P, W] * x[P, W]) ``` --- # /reference/param # Param ## `Param` Coefficients over a set product. A parameter is data, not a model object: it has no columns and produces no rows. It supplies a term's coefficient and a constraint's right-hand side. The array of a parameter declares `absence="empty"`. A coordinate it does not have contributes no coefficient. Absence is the additive identity a sum requires. ### `Param.from_dense(name, sets, values)` Every cell of `values` as a coefficient. `values.shape` must equal the sizes of `sets`, in order; a mismatch raises `ValueError` with both shapes. ```python import numpy as np from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) print(cost.dims, cost.nnz) ``` Output: ```text ('P', 'W') 6 ``` ### `Param.from_long(name, sets, columns, values)` Coefficients from one label column per set and one value column. `columns` is a mapping keyed by set name. Each column and `values` are read in parallel and have the same length. ```python import numpy as np from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_long( "cost", (P, W), {"P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris"])}, np.array([2.0, 1.0]), ) print(cost.nnz) ``` Output: ```text 2 ``` A label column of a different length raises `ValueError`; the message gives the parameter, the column, its length and the value column's. ```python raises=ValueError import numpy as np from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) Param.from_long( "cost", (P, W), {"P": np.array(["lisbon"]), "W": np.array(["berlin", "paris"])}, np.array([2.0, 1.0]), ) ``` Raises ValueError: ```text ValueError: parameter 'cost': label columns have lengths {'P': 1, 'W': 2} and the value column has length 2; pass columns of equal length ``` ### `Param.from_positions(name, sets, index, values)` Coefficients from an index matrix and one value column. `index` has one row per set, and each column identifies one member by its position in each set. A caller that already has positions resolves no labels. The columns are in any order, and each member appears once. ```python import numpy as np from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_positions("cost", (P, W), [[1, 0], [1, 0]], [1.0, 2.0]) print(cost.array.to_dense()) ``` Output: ```text [[2. 0. 0.] [0. 1. 0.]] ``` A position outside its set raises `ValueError`, as does a member given twice. ### `Param.from_array(name, sets, array)` Coefficients from a nimblend array over the sets' names, in any order. The labels are members of the sets, in any order and extent. The parameter has an entry where the array has one, under absence `'empty'` and `'unknown'` alike. A label outside its set, other dimension names and a set with no members raise `ValueError`; an object that is not a nimblend array raises `TypeError`. ```python import numpy as np import nimblend as nb from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_array( "cost", (P, W), nb.from_dense( np.array([[1.0]]), {"P": np.array(["porto"]), "W": np.array(["paris"])} ), ) print(cost.array.to_dense()) ``` Output: ```text [[0. 0. 0.] [0. 1. 0.]] ``` ### Members | Member | Returns | | --- | --- | | `dims` | the sets it is indexed over | | `nnz` | the number of coefficients | | `materialise()` | the coefficients as a `nimblend` array | | `param[sets]` | a reference, with the sets given checked against `dims` | A label in place of a set fixes that dimension at one member: the coefficients at that member are read and the dimension leaves the reference. ```python import numpy as np from nimopt import Param, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) print(cost[P, W].dims) print(cost[P, "berlin"].dims) ``` Output: ```text ('P', 'W') ('P',) ``` ## `Coefficient` What a term reads as its coefficient: a `name` to report, the `dims` it is indexed over, the array it `materialise()`s to, and a reading at its sets. A parameter read at its sets is a coefficient, and so is an arithmetic combination of coefficients. One interface therefore covers both. `+`, `-`, `*`, `/` and a power by a number combine coefficients. The combination is symbolic: it contains references, derives its dimensions from its operands, and is evaluated once, when the term it multiplies is materialised. It can therefore be written in a definition before any data exists. ```python import numpy as np from nimopt import Param, Set G = Set("G", np.array(["base", "peak"])) T = Set("T", np.arange(3)) price = Param.from_dense("fuel_price", (G, T), np.full((2, 3), 30.0)) eta = Param.from_dense("efficiency", (G, T), np.array([[0.5] * 3, [0.4] * 3])) unit_cost = price[G, T] / eta[G, T] print(unit_cost.name, unit_cost.dims) print(unit_cost[G, T].materialise().to_dense()[:, 0]) ``` Output: ```text (fuel_price / efficiency) ('G', 'T') [60. 75.] ``` A parameter has no arithmetic of its own. It is read at its sets, and the references combine. ```python raises=TypeError import numpy as np from nimopt import Param, Set G = Set("G", np.array(["a", "b"])) price = Param.from_dense("price", (G,), np.array([1.0, 2.0])) eta = Param.from_dense("eta", (G,), np.array([0.5, 0.4])) price / eta ``` Raises TypeError: ```text TypeError: parameter 'price' is over ('G',) and expresses no coefficient until it is read; read it at its sets as price[G] ``` --- # /reference/sets # Sets and domains ## `Set` ``` Set(name, labels) ``` A named dimension with labels. `labels` is an array; `name` is what every reference to the dimension uses. | Member | Returns | | --- | --- | | `name`, `labels` | the declared name and labels | | `len(set)` | the number of members | | `position_of(labels)` | the position of each label given | | `coord` | the coordinate the labels resolve through | | `cyclic` | the same set, with a lag that wraps instead of dropping | | `set - 1` | the set lagged, dropping the members outside the set | ```python import numpy as np from nimopt import Set T = Set("T", np.array(["t0", "t1", "t2"])) print(T.name, len(T)) print(T.labels) print(T.position_of(np.array(["t2", "t0"]))) ``` Output: ```text T 3 ['t0' 't1' 't2'] [2 0] ``` ## `Alias` ``` Alias(name, set) ``` A second name for a set, sharing its labels and its coordinate. A parameter over a set and its alias is an ordinary two-dimensional array. A model relates a set to itself without declaring a second set. No labels are copied: the alias uses the coordinate the set already built. ```python import numpy as np from nimopt import Alias, Param, Set N = Set("N", np.array(["a", "b"])) M = Alias("M", N) flow = Param.from_dense("flow", (N, M), np.array([[0.0, 1.0], [1.0, 0.0]])) print(flow.dims) print(flow.materialise().to_dense()) ``` Output: ```text ('N', 'M') [[0. 1.] [1. 0.]] ``` ## `product` ``` product(sets) ``` Every member of a set product, as a domain. Passed to `over=`, it declares the rows of a constraint explicitly, for a constraint whose terms each cover some of its rows. A factor is a set or a domain: a set contributes every member, and a domain its own members. `product((subset((G,), columns), T))` is the members of `G` in `columns` crossed with every member of `T`. Two factors over one dimension raise `ValueError`. ## `subset` ``` subset(sets, columns) ``` The members of a set product a model uses, given by label. `columns` contains one label column per set, keyed by the name of the set. The columns are read in parallel: the k-th entry of each column belongs to the same member. The result is a list of members, not a cross product. ## `subset_of` ``` subset_of(sets, index) ``` The same, given by position. Each column of `index` is one member. A caller with positions passes them directly and builds no labels to resolve back. ```python import numpy as np from nimopt import Set, product, subset, subset_of P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) print(product((P, W)).size) print( subset( (P, W), {"P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris"])} ).size ) print(subset_of((P, W), np.array([[0, 1], [0, 1]])).size) ``` Output: ```text 6 2 2 ``` The product has six members. Both subsets have two, `lisbon` with `berlin` and `porto` with `paris`. The columns are read in parallel. --- # /reference/solution # Solution ## `Solution` Returned by `Model.solve`. Primal and dual values, returned over the sets they were declared over. | Member | Returns | | --- | --- | | `status` | the outcome the solver reported | | `feasible` | whether the solver reports a primal-feasible point | | `objective` | the objective value of that point | | `bound` | the bound on the optimal objective the solver proved, or `None` | | `gap` | the relative distance from the objective to the bound, or `None` | | `primal(name)` | the named variable's values over its own sets | | `dual(name, kind=None)` | a constraint's duals over its free sets, or a variable's reduced costs over its own sets | | `has_duals` | whether `dual` returns values: status `optimal`, and a solver that reports duals | `status` and `feasible` are readable whatever the solver reported. `objective` and `primal` raise `ValueError` where `feasible` is False. They raise at status `unbounded` and `unbounded_or_infeasible` whatever `feasible` reports. An unbounded model has no optimal value, and `bound` and `gap` are `None` there. A solve stopped at a limit reports `feasible` True where the solver found a point, and those reads then return it. `dual` raises `ValueError` where `status` is not `optimal`, and for a model with integer columns. Read `status` first, or `has_duals` before `dual`. `primal` and `dual` raise `KeyError` for a name the model does not declare. `primal` takes a variable, and its message reports a constraint name as one `dual` reads. `dual` takes either, and its message lists the declared constraints and variables. Both raise `KeyError` for a name declared after the solve. A model declares its constraints and its variables in two registries, so one name identifies one of each. `dual` raises `ValueError` for such a name and reads it under `kind="constraint"` or `kind="variable"`. Any other `kind` raises `ValueError`. `dual` returns a reduced cost for a variable: its objective coefficient less the duals of the rows it appears in, weighted by its coefficients in them, in the model's own objective under either sense. nimopt derives the value from the row duals the solver reports, so the convention does not vary by solver. The values follow the dual solution the solver returns. A degenerate model has more than one such solution, and two solvers can report different reduced costs for it. A reduced cost follows the variable's members by the rule `primal` follows: a `DenseArray` over a full product, a `SparseArray` over a subset. `bound` is a lower bound on the optimal objective under sense `min` and an upper bound under sense `max`. It is `None` where the solver reports none. For a model without integer columns it is the objective at status `optimal` and `None` at any other status. `bound` is readable at every status, and the solvers report none at status `unbounded` and `unbounded_or_infeasible`. `gap` is `abs(objective - bound) / abs(objective)`. It is `None` where `feasible` is False, where `bound` is `None`, and at status `unbounded` and `unbounded_or_infeasible`. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(solution.status) print(solution.objective) print(solution.primal("x").to_dense()) print(solution.dual("demand").to_dense()) ``` Output: ```text optimal 135.0 [[20. 0. 10.] [ 0. 15. 5.]] [3. 1. 6.] ``` Reading a value where the solver reports no feasible point raises `ValueError`; the message gives the status. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("infeasible") x = m.var("x", (P, W)) m.constraint("floor", Sum(W, x[P, W]) >= 10.0) m.constraint("ceiling", Sum(W, x[P, W]) <= 1.0) m.set_objective(Sum(P, W, x[P, W])) m.solve().objective ``` Raises ValueError: ```text ValueError: status is 'infeasible' and the solver reports no feasible point; read `status` before reading values ``` ## The array type of a value A variable over a full product has a value at every cell of its frame. The solver returns those values in column order, and they reshape into a `DenseArray` with no index built. A variable over a subset has values at its members alone. A dense frame would be the grid the declaration avoids, and those values remain a `SparseArray`. A dual follows the rows of its constraint by the same rule. Every array declares `absence="unknown"`. A coordinate the model does not have has no value, and combining the results of two models adds no zero for it. --- # /reference/solvers # Solvers ## `available` and `capabilities` `available()` lists every adapter whose backend can be imported in the current environment, with the capabilities each declares. `capabilities(name)` reports for an adapter whether or not its backend is installed. A descriptor describes the adapter as shipped. A caller reads one to choose what to install. A descriptor describes the adapter, not the library behind it: a solver feature the adapter does not call is `absent`. What `available()` lists depends on the machine; what `capabilities(name)` reports does not. ```python from nimopt import available, capabilities print("highs" in available()) print(capabilities("highs")) print(capabilities("gurobi")) print(capabilities("mosek")) ``` Output: ```text True highs integrality native · duals native · conflict native · ray native rejects duals+integrality gurobi integrality native · duals native · conflict native · ray native rejects duals+integrality mosek integrality native · duals native · conflict absent · ray native rejects duals+integrality ``` ## `Capabilities` | Member | Returns | | --- | --- | | `solver` | the adapter's name | | `support` | one of `"native"` or `"absent"` per capability | | `rejected` | the pairs this adapter rejects together | | `supports(capability)` | whether the adapter handles it at all | | `rejects(one, other)` | whether it rejects the two together | The capabilities are `integrality`, `duals`, `conflict` and `ray`. A flat set is insufficient: a solver can support two and reject their combination. Every adapter rejects `integrality` with `duals`. A mixed-integer model's duals are not the relaxation's. A model with integer columns has no duals at all, and `Solution.dual` raises. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set, Sum T = Set("T", np.arange(2)) one = Param.from_dense("one", (T,), np.ones(2)) m = Model("m") x = m.var("x", (T,), upper=3.0, integer=True) m.constraint("cap", one[T] * x[T] <= 2.0) m.set_objective(Sum(T, one[T] * x[T])) m.solve().dual("cap") ``` Raises ValueError: ```text ValueError: model 'm' has integer columns and 'highs' reports no duals for it; read primal values only ``` ## `Session` Returned by `Model.session(solver="highs", options=None)`. The model of one solver, opened on one assembled model and kept open. A solve passes the matrix to the solver, and a later query reads the same solved instance. | Member | Returns | | --- | --- | | `assembled` | the matrix the session was opened on | | `solver`, `capabilities` | which adapter, and what it can do | | `status` | what the last solve reported, or `None` before one | | `solve()` | a `Solution` | | `close()` | releases the solver's model | `Model.solve()` opens a session, solves and closes it. A caller who needs only a solution needs no session. The session assembles the matrix when it opens. `solve()` and `diagnose()` raise `ValueError` where the model declares a variable or a constraint, or sets its objective, after that. Open a new session on the changed model. ```python import numpy as np from nimopt import Model, Param, Set, Sum SNAP = Set("snapshot", np.arange(3)) GEN = Set("generator", np.array(["wind", "gas"])) p_max = Param.from_dense("p_max", (GEN,), np.array([10.0, 20.0])) load = Param.from_dense("load", (SNAP,), np.array([25.0, 20.0, 5.0])) cost = Param.from_dense("cost", (GEN,), np.array([1.0, 5.0])) m = Model("dispatch", sense="min") p = m.var("p", (SNAP, GEN), lower=0.0, upper=p_max) m.constraint("balance", Sum(GEN, p[SNAP, GEN]) == load[SNAP]) m.set_objective(Sum(SNAP, GEN, cost[GEN] * p[SNAP, GEN])) with m.session() as session: solution = session.solve() print(solution.status, solution.objective) print(session.status) ``` Output: ```text optimal 150.0 optimal ``` ## `Diagnosis` Returned by `Session.diagnose()`. Why a model did not solve, as the rows and columns that explain it. | Member | Returns | | --- | --- | | `status`, `solver` | what the solve reported, and which adapter | | `method` | `"native"` where the solver computed the conflict | | `conflict` | one `Row` per conflicting row, or `None` where the model is not infeasible | | `columns` | one `ColumnBound` per column of the conflict, with the bounds the model declares | | `ray` | one `RayTerm` per column the ray moves, or `None` where the model is not unbounded | `conflict` contains the same `Row` that `Model.row` returns, and a conflicting row reads in one format. Two solvers may return different irreducible sets. Removing either set makes the model feasible. The two sets need not be equal. A conflict is computed on the solved backend, and the session keeps that backend available. ```python import numpy as np from nimopt import Model, Param, Set, Sum SNAP = Set("snapshot", np.arange(3)) GEN = Set("generator", np.array(["wind", "gas"])) p_max = Param.from_dense("p_max", (GEN,), np.array([10.0, 20.0])) load = Param.from_dense("load", (SNAP,), np.array([25.0, 100.0, 5.0])) cost = Param.from_dense("cost", (GEN,), np.array([1.0, 5.0])) m = Model("dispatch", sense="min") p = m.var("p", (SNAP, GEN), lower=0.0, upper=p_max) m.constraint("balance", Sum(GEN, p[SNAP, GEN]) == load[SNAP]) m.set_objective(Sum(SNAP, GEN, cost[GEN] * p[SNAP, GEN])) with m.session() as session: print(session.solve().status) print(session.diagnose()) ``` Output: ```text infeasible infeasible highs conflict native balance[snapshot=1] row 1 1·p[1,wind] + 1·p[1,gas] == 100 bound p[snapshot=1, generator='wind'] [0, 10] bound p[snapshot=1, generator='gas'] [0, 20] ``` HiGHS computes its conflict over the model's linear relaxation. A model that is feasible as an LP and infeasible only through its integrality therefore yields no conflict, and the adapter raises. It reports no row it did not prove. Gurobi's conflict covers the integrality. Mosek's adapter computes no conflict. A session on it raises, and the message refers to `capabilities("mosek")`. ## `options` and `Option` `options()` lists every option a caller can set, under the names of `nimopt`. An option outside the list raises `ValueError`. No option is passed to a solver that ignores it, and an unrecognized name stops the solve. `options(solver)` lists the same options with the name and the values of that solver. A caller follows those names into the documentation of the solver. ```python from nimopt import options for option in options("highs"): print(f"{option.name:<16} {option.native}") ``` Output: ```text time_limit time_limit iteration_limit simplex_iteration_limit node_limit mip_max_nodes mip_gap mip_rel_gap mip_abs_gap mip_abs_gap feasibility_tol primal_feasibility_tolerance optimality_tol dual_feasibility_tolerance threads threads seed random_seed log output_flag presolve presolve method solver newton_system hipo_system crossover run_crossover pdlp_tol pdlp_optimality_tolerance ``` An `Option` has a `name`, the `kind` it takes, what it `does`, and its `choices` where it takes one of a set. Read for a solver it also has `native` and `native_choices`. | Option | Takes | Does | `highs` | `gurobi` | `mosek` | | --- | --- | --- | --- | --- | --- | | `time_limit` | float | seconds the solver may run for | `time_limit` | `TimeLimit` | `optimizer_max_time` | | `iteration_limit` | int | simplex iterations the solver may take | `simplex_iteration_limit` | `IterationLimit` | `sim_max_iterations` | | `node_limit` | int | branch-and-bound nodes the solver may explore | `mip_max_nodes` | `NodeLimit` | `mio_max_num_branches` | | `mip_gap` | float | relative gap at which a mixed-integer solve stops | `mip_rel_gap` | `MIPGap` | `mio_tol_rel_gap` | | `mip_abs_gap` | float | absolute gap at which a mixed-integer solve stops | `mip_abs_gap` | `MIPGapAbs` | `mio_tol_abs_gap` | | `feasibility_tol` | float | how far a primal solution may miss a row | `primal_feasibility_tolerance` | `FeasibilityTol` | `basis_tol_x` | | `optimality_tol` | float | how far a dual solution may miss a bound | `dual_feasibility_tolerance` | `OptimalityTol` | `basis_tol_s` | | `threads` | int | threads the solver may use; 0 leaves it the choice | `threads` | `Threads` | `num_threads` | | `seed` | int | the seed the solver randomizes from | `random_seed` | `Seed` | `mio_seed` | | `log` | bool | whether the solver writes its own iteration log | `output_flag` | `OutputFlag` | `log` | | `presolve` | `off` / `choose` / `on` | how hard the solver presolves | `presolve` | `Presolve` | `presolve_use` | | `method` | `choose` / `simplex` / `barrier` / `hipo` / `pdlp` | the algorithm the solver runs | `solver` | `Method` | `optimizer` | | `newton_system` | `choose` / `augmented` / `normaleq` | the Newton system an interior point method factorizes | `hipo_system` | not supported | not supported | | `crossover` | `choose` / `off` / `on` | whether an interior point is moved to a vertex after the solve | `run_crossover` | `Crossover` | `intpnt_basis` | | `pdlp_tol` | float | relative tolerance at which the first-order method stops | `pdlp_optimality_tolerance` | not supported | not supported | A choice each solver writes differently is given once and translated. The value a caller writes has one meaning for every solver. No solver supports every option or every choice. `newton_system` and `pdlp_tol` are specific to HiGHS, and so are `hipo` and `pdlp` under `method`. Passing one of them to Gurobi or Mosek raises, and the message identifies it. No solver runs a different algorithm in its place. Mosek runs only its mixed-integer optimizer on a model with integer columns. `method` is `choose` there, and any other choice raises. The guide on [interior point and first-order methods](/guides/highs-methods) gives the memory each method requires. It also gives the installation of a HiGHS with HiPO and a GPU. ## Progress reporting `build()`, `assemble()`, `session()` and `solve()` take `progress=`. `progress=True` draws a report in a terminal and nothing where output is redirected. A reporter given by the caller is used as passed, and a notebook or another interface writes through it. A reporter implements `start(total, what)`, `step(done, what)` and `done()`. That is the whole contract. The report covers building. Its resolution is the structure of the model. The measuring pass counts constraints and the writing pass counts nonzeros. A model with one constraint of one term reports one step and no fraction. `log=True` requests the log of the solver. Building finishes before a solver starts, and the report and the log do not interleave. --- # /reference/variable # Variable ## `Variable` Returned by `Model.var`. A variable over a set product, or over a subset of one. ``` Model.var(name, sets, subset=None, lower=0.0, upper=inf, integer=False) ``` | Argument | Meaning | | --- | --- | | `name` | the name `Solution.primal` reads it back by | | `sets` | the dimensions it is declared over | | `subset` | the members it has; the full product when omitted | | `lower`, `upper` | the bound every one of its columns takes | | `integer` | whether its columns are integral | The variable's columns are a virtual coordinate: a member's column is computed from its multi-index by stride arithmetic for a full product, or is its rank among the codes of a subset. Nothing stores a column index. A variable over millions of columns therefore costs only its members. | Member | Returns | | --- | --- | | `dims` | the names of the sets it is over | | `n_columns` | the number of columns it occupies | | `labels_at(positions)` | the labels of the columns at those positions, per set; a position is a column less `start` | | `domain()` | the members it has | | `terms()` | its coefficients over `(*dims, COLUMN)` | | `variable[sets]` | a one-term expression referencing it | ```python import numpy as np from nimopt import Model, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) open_plant = m.var("open_plant", (P,), lower=0.0, upper=1.0, integer=True) print(x.dims, x.n_columns) print(open_plant.n_columns) print(m.n_columns) print(m.integrality()) ``` Output: ```text ('P', 'W') 6 2 8 [0 0 0 0 0 0 1 1] ``` Each variable occupies the next range of the single column space of the model. `m.n_columns` counts every column declared so far. ## A variable over no dimension The bracket of a variable lists the dimensions it is declared over. A variable over no dimension takes no bracket and enters a row on its own. It is one column: a value-at-risk level, a budget slack, or a bound shared by every row of a family. `theta[()]` is the same term written out. A variable over one or more dimensions expresses no term until it is read. Using one without a bracket raises `TypeError` and reports the reading it requires. The same rule applies to a parameter, read as `cost[G, T]`, and as `k` over no dimension. Comparing a variable expresses a row, and `==` between two variables expresses a row as well. A list of variables therefore cannot be searched with `in` or `.index`. Both compare their items and raise `TypeError` and report the first variable they compare. A dict and a set match on identity. Store variables in one of them, or search them by `name`. ```python import numpy as np from nimopt import Model, Set, Sum S = Set("S", np.array(["s1", "s2"])) m = Model("cvar", sense="min") theta = m.var("theta", (), lower=-np.inf) p = m.var("p", (S,)) m.constraint("tail", theta - Sum(S, p[S]) >= 0.0) m.set_objective(theta) print(m.n_columns, m.n_rows) print(m.constraints["tail"].relation) ``` Output: ```text 3 1 theta - Sum(S, p[S]) >= 0 ``` ## `COLUMN` and `ROW` The dimension names `nimopt` reserves. `COLUMN` is `"__column__"` and `ROW` is `"__row__"`. Both are written so that no ordinary set name collides with them. A variable's terms are an array over `(*dims, COLUMN)`, and a constraint's block is one over `(ROW, COLUMN)`. That is the whole correspondence between a model and its matrix: the column space is a dimension, and the array is the matrix. ```python import numpy as np from nimopt import COLUMN, ROW, Model, Set P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) print(COLUMN, ROW) print(x.terms().dims) print(x.domain().dims) ``` Output: ```text __column__ __row__ ('P', 'W', '__column__') ('P', 'W') ``` A caller writes neither name. Both appear when a `nimblend` array from inside a model is inspected. --- # /nimblend/arrays # nimblend arrays `nimblend` is the layer below `nimopt`. Its vocabulary is dimensions, labels, entries and alignment, and it contains no optimization term. A `nimopt` caller uses these names when reading a solution or inspecting what a model built. **Import from `nimblend` itself, never from a submodule. Never read an array's `.index` or `.data` or a domain's `.codes`.** Those are the raw index matrix, the value buffer and the ravelled members of the layer below the array layer. Reading them bypasses the contract, and a test in this repository fails on any of them. Do not assemble one either. `SparseArray.from_canonical` below takes an index matrix built by the caller, and that is the work of the array layer. No module of `nimopt` calls it, and a test enforces that. A `Domain` returns the array over its own members instead. ## `Array` A labeled N-dimensional array. `Array` is the contract both implementations satisfy, and the interface a caller writes against. `absence` declares the meaning of a coordinate the array does not have: `"empty"` that it contributes nothing, `"unknown"` that it was not modeled. Operators and reductions follow from that declaration. The declaration is part of the definition of the array, not a hint. | Member | Returns | | --- | --- | | `dims`, `shape`, `nnz` | the dimensions, their sizes, and the number of entries | | `coords` | the coordinate of each dimension, resolving its labels | | `absence` | the meaning of a coordinate the array does not have | | `as_empty()`, `as_unknown()` | the array under the other absence declaration | | `values()`, `coordinates()` | the entries and their multi-indices | | `domain()` | the coordinates the array has | | `sum`, `min`, `max`, `mean` | reductions over named dimensions | | `sel`, `restrict` | a selection by label, and by domain | | `rename`, `transpose`, `expand`, `conform` | reshaping the dimensions | | `shift`, `roll` | a lag that drops, and one that wraps | | `group` | entries combined into a destination | | `to_dense(fill=None)` | the entries as an ndarray | | `+`, `-`, `*`, `/`, unary `-` | arithmetic over one frame, and with a scalar | The arithmetic is part of the contract, not an implementation's own: a caller writing `coefficient * columns` is writing against `Array`. Both implementations behave alike, including over frames that differ. One frame nested inside the other broadcasts over the wider; frames sharing some dimensions align on those and multiply out the rest. Frames sharing no dimension raise. A product mixing the two implementations returns a `SparseArray`: a product intersects presence and has at most the entries of the sparse operand. Absence and zero remain distinct. A stored `0.0` is a coordinate present with the value zero. A coordinate the array does not have is a different case. ```python import numpy as np import nimblend as nb labels = {"A": np.array(["a0", "a1"]), "B": np.array(["b0", "b1"])} array = nb.SparseArray.from_dense(np.array([[1.0, 0.0], [0.0, 2.0]]), labels) print(isinstance(array, nb.Array)) print(array.dims, array.shape, array.absence) print(array.nnz) print(array.sum("B").to_dense()) ``` Output: ```text True ('A', 'B') (2, 2) empty 4 [1. 2.] ``` `from_dense` stores every cell it was given, and this array has four entries and not two. The zeros are stored, and a stored value is present. ## `SparseArray` Entries in canonical order, under a coordinate per dimension. An entry that is not stored is absent. It is what an operation returns when the result is sparse, and what `Solution.primal` returns for a variable over a subset: a dense frame there would be the grid the variable was declared to avoid. | Member | Returns | | --- | --- | | `SparseArray.from_dense(values, labels)` | every cell of an ndarray | | `SparseArray.from_canonical(index, values, coords, dims)` | entries already in order | | `as_empty()`, `as_unknown()` | the same entries under the other declaration | | `to_csr()` | the entries as compressed rows | ```python import numpy as np import nimblend as nb labels = {"A": np.array(["a0", "a1"]), "B": np.array(["b0", "b1"])} array = nb.SparseArray.from_dense(np.array([[1.0, 0.0], [0.0, 2.0]]), labels) print(array.absence) print(array.as_unknown().absence) print(array.to_dense()) ``` Output: ```text empty unknown [[1. 0.] [0. 2.]] ``` ## `DenseArray` An ndarray over labeled dimensions, distinguishing absence from zero. The storage of presence follows the absence declaration: the two declarations require opposite behavior from an operator. An `"unknown"` array marks absence with NaN. NaN propagates through arithmetic at no cost and requires no storage beside the values. An `"empty"` array stores a boolean mask: absence is the additive identity there, and substituting it costs less than marking it. `Solution.primal` returns one of these for a variable over a full product. The solver returns a value at every cell of the frame in column order, and the values reshape with no index built. ```python import numpy as np import nimblend as nb coords = { "A": nb.StoredCoord(np.array(["a0", "a1"])), "B": nb.StoredCoord(np.array(["b0", "b1"])), } array = nb.DenseArray(np.array([[1.0, 0.0], [0.0, 2.0]]), coords, ("A", "B")) print(array.absence, array.nnz) print(array.present) print(array.as_unknown().absence) ``` Output: ```text empty 4 [[ True True] [ True True]] unknown ``` ## Building an array from columns `nimblend` itself provides the two constructors a caller uses when the data is not already an ndarray. | Constructor | Returns | | --- | --- | | `from_long(dims, coords, labels, values)` | one label column per dimension and one value column | | `from_dense(values, labels)` | every cell of an ndarray | | `is_canonical(index, shape)` | whether buffers are in the order `from_canonical` takes | `from_long` resolves each label through the coordinate that dimension already has. A caller with coordinates of its own therefore gives its entries as labels and resolves no position itself. The columns are read in parallel and must have equal length. ```python import numpy as np import nimblend as nb coords = { "t": nb.StoredCoord(np.array([2030, 2040])), "r": nb.StoredCoord(np.array(["DE", "FR"])), } arr = nb.from_long( ("t", "r"), coords, {"t": np.array([2030, 2040, 2040]), "r": np.array(["DE", "DE", "FR"])}, np.array([5.0, 6.0, 7.0]), ) print(arr.nnz, arr.to_dense()[1, 1]) ``` Output: ```text 3 7.0 ``` A column of a different length raises `ValueError`; the message gives the column and both lengths. ```python raises=ValueError import numpy as np import nimblend as nb nb.from_long( ("t",), {"t": nb.StoredCoord(np.array([2030, 2040]))}, {"t": np.array([2030, 2040])}, np.array([1.0]), ) ``` Raises ValueError: ```text ValueError: label columns have lengths {'t': 2} and the value column has length 1; pass columns of equal length ``` `is_canonical` reports whether buffers meet the precondition of `SparseArray.from_canonical`: entries sorted by ravel key, with no repeat. Verifying it inside `from_canonical` would cost the ravel that this path avoids. ```python import numpy as np import nimblend as nb ordered = np.array([[0, 0, 1], [0, 1, 0]], dtype=np.int32) print(nb.is_canonical(ordered, (2, 2))) print(nb.is_canonical(ordered[:, ::-1].copy(), (2, 2))) ``` Output: ```text True False ``` ## The frame of a binary result `combined_dims(left, right)` returns the dimensions a binary operator's result has, from the dimensions of the two operands alone. A caller reads it before materialising either operand. A combination therefore reports its frame while its data is unbound. | Operands | Result | | --- | --- | | equal frames | that frame, in its order | | one frame nested in the other | the wider | | frames that overlap | the left, then the dimensions only the right has | | frames sharing no dimension | raises | ```python import nimblend as nb print(nb.combined_dims(("P", "Q"), ("Q", "R"))) print(nb.combined_dims(("P",), ("P", "Q"))) ``` Output: ```text ('P', 'Q', 'R') ('P', 'Q') ``` Frames sharing no dimension have no common dimension to align on, and `combined_dims` raises. Their combination would be an outer product. ```python raises=ValueError import nimblend as nb nb.combined_dims(("P",), ("Q",)) ``` Raises ValueError: ```text ValueError: frames ('P',) and ('Q',) share no dimension; pass operands that share a dimension ``` ## Adding several arrays `sum_arrays(arrays)` adds several arrays over one frame in one merge. The result equals adding them in order with `+`. A chain of `+` merges the running sum once per operand, and `sum_arrays` sorts the entries of all operands once. Under `"empty"` the result has every coordinate of any operand. Under `"unknown"` it has only the coordinates of every operand. The operands have the same dimensions, labels and absence. ```python import numpy as np import nimblend as nb labels = {"A": np.array(["a0", "a1"]), "B": np.array(["b0", "b1"])} first = nb.SparseArray.from_dense(np.array([[1.0, 0.0], [0.0, 2.0]]), labels) second = nb.SparseArray.from_dense(np.array([[10.0, 20.0], [30.0, 40.0]]), labels) print(nb.sum_arrays([first, second, first]).to_dense()) ``` Output: ```text [[12. 20.] [30. 44.]] ``` ## Densifying an unknown array An array declaring `"unknown"` that does not have every coordinate of its frame raises on `to_dense()` without a fill. There is no value it can place at the rest, and choosing one silently would invent a value. ```python raises=ValueError import numpy as np import nimblend as nb labels = {"A": np.array(["a0", "a1"]), "B": np.array(["b0", "b1"])} partial = nb.SparseArray.from_canonical( np.array([[0], [0]], dtype=np.int32), np.array([1.0]), { "A": nb.StoredCoord(labels["A"]), "B": nb.StoredCoord(labels["B"]), }, ("A", "B"), absence="unknown", ) partial.to_dense() ``` Raises ValueError: ```text ValueError: absence is 'unknown' and the array has no value at 3 of 4 coordinates; pass fill= to to_dense() ``` With a fill value, the grid is returned. ```python import numpy as np import nimblend as nb labels = {"A": np.array(["a0", "a1"]), "B": np.array(["b0", "b1"])} partial = nb.SparseArray.from_canonical( np.array([[0], [0]], dtype=np.int32), np.array([1.0]), { "A": nb.StoredCoord(labels["A"]), "B": nb.StoredCoord(labels["B"]), }, ("A", "B"), absence="unknown", ) print(partial.to_dense(fill=np.nan)) ``` Output: ```text [[ 1. nan] [nan nan]] ``` --- # /nimblend/domains # nimblend domains ## `Domain` A sorted, unique set of multi-indices over named dimensions. A coordinate resolves labels for one dimension, and a domain resolves them for a tuple of dimensions. It reports the multi-indices it has, the position of each, and the multi-index at a given position. It records which coordinates are present, and it records no numbering origin. | Constructor | Returns | | --- | --- | | `Domain.full(dims, coords)` | every coordinate of the product `dims` spans | | `Domain.from_labels(dims, coords, labels)` | a domain from one label column per dimension | | `Domain.from_coordinates(dims, coords, index)` | a domain from an index matrix of one row per dimension | | Member | Returns | | --- | --- | | `size`, `dims`, `shape`, `coords` | the number of members, the dimensions, their extents and their coordinates | | `is_full` | whether every coordinate of the product is present | | `coordinates(positions=None)` | the multi-index of each member, or of the members at `positions`, as an int32 index matrix | | `labels(positions=None)` | each member's label, or the labels of the members at `positions`, per dimension | | `intersect(other)` | the members both have | | `union(other)` | the members either has | | `difference(other)` | the members this one has and `other` does not | | `positions_of(array)` | each entry of `array` as its position here, `-1` where absent | | `positions_of_coordinates(index)` | each column of an index matrix as its position here, `-1` where absent | | `expand(dims, coords)` | every member crossed with the full extent of the named dimensions | | `transpose(*dims)` | the same members, over the dimensions in the order given, or reversed when none are given | | `project(dims)` | the distinct coordinates of the members over `dims`, in the order given | | `as_coord()` | the domain as a coordinate, each member at its rank | | `array(values, absence="empty")` | the members with one value each, as a `SparseArray` | | `identity(into, coord, start=0)` | each member paired with its own position along `into`, valued 1.0 | That table is the whole surface. **A domain's `codes` are the raw ravelled members of the layer below it. An array's `.index` and `.data` are its raw buffers. Never read them.** `coordinates()` and `labels()` report which members are present, and `positions_of_coordinates` reports the position of one. With `positions`, `coordinates` and `labels` decode those members alone. `as_coord` numbers them, and `array` and `identity` return an array over them. No layer above builds an index matrix. A test in this repository fails on a read of any of the three. The label columns of `from_labels` are read in parallel: the k-th entry of each column belongs to the same member. A domain is a list of members, not a cross product. ```python import numpy as np import nimblend as nb coords = { "P": nb.StoredCoord(np.array(["lisbon", "porto"])), "W": nb.StoredCoord(np.array(["berlin", "paris", "rome"])), } full = nb.Domain.full(("P", "W"), coords) pairs = nb.Domain.from_labels( ("P", "W"), coords, {"P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris"])}, ) print(full.size, pairs.size) print(pairs.coordinates()) print(pairs.labels()) print(full.intersect(pairs).size, full.difference(pairs).size) ``` Output: ```text 6 2 [[0 1] [0 1]] {'P': array(['lisbon', 'porto'], dtype='= 0`. `as_coord` returns the domain as a coordinate: a member's position is its rank among the members present. This is how a dimension over a subset of a product is numbered. ```python import numpy as np import nimblend as nb coords = { "P": nb.StoredCoord(np.array(["lisbon", "porto"])), "W": nb.StoredCoord(np.array(["berlin", "paris", "rome"])), } pairs = nb.Domain.from_labels( ("P", "W"), coords, {"P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris"])}, ) asked = np.array([[0, 1, 1], [0, 1, 2]], dtype=np.int32) print(pairs.positions_of_coordinates(asked)) print(pairs.positions_of_coordinates(asked) >= 0) numbered = pairs.as_coord() print(numbered.to_position(asked[:, :2])) ``` Output: ```text [ 0 1 -1] [ True True False] [0 1] ``` `("porto", "rome")` is not a member, and the query returns `-1`. The other two are the first and second members of the domain, and `as_coord()` numbers them 0 and 1. ## Crossing a domain with further dimensions `expand` replicates every member across the full extent of the given dimensions. That enumerates the coordinates a term can have, before the test of which ones it does have. The new dimensions are appended, and `transpose` reads the result in another order. The code of a member is its own code scaled by the appended extent, plus each position within it. The cross product is therefore arithmetic on the members, and it builds no index matrix. ```python import numpy as np import nimblend as nb coords = { "P": nb.StoredCoord(np.array(["lisbon", "porto"])), "W": nb.StoredCoord(np.array(["berlin", "paris", "rome"])), "H": nb.StoredCoord(np.array([0, 1])), } pairs = nb.Domain.from_labels( ("P", "W"), coords, {"P": np.array(["lisbon", "porto"]), "W": np.array(["berlin", "paris"])}, ) hourly = pairs.expand(("H",), coords) print(hourly.dims, hourly.size) print(hourly.labels()) print(hourly.transpose("H", "P", "W").dims) ``` Output: ```text ('P', 'W', 'H') 4 {'P': array(['lisbon', 'lisbon', 'porto', 'porto'], dtype='=` or `==` produces a `Relation`: the expression, the sense, and the right-hand side. A relation is not yet part of the model. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) m = Model("transport") x = m.var("x", (P, W)) rule = Sum(W, x[P, W]) <= supply[P] print(type(rule).__name__, rule.sense) ``` Output: ```text Relation <= ``` ## Adding a constraint `m.constraint(name, relation)` adds the relation to the model under a name and returns the `Constraint`. The name identifies the constraint's rows in the matrix and its dual values in the solution. A constraint produces one row per member of its expression's frame. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) supply_rows = m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) demand_rows = m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) print(supply_rows.n_rows, supply_rows.nnz) print(demand_rows.n_rows, demand_rows.nnz) print(m.n_rows, m.nnz) ``` Output: ```text 2 6 3 6 5 12 ``` The supply expression has frame `(P,)` and produces two rows; the demand expression has frame `(W,)` and produces three. Each supply row has three nonzeros, one per route out of its plant, and each demand row two, one per route into its warehouse: twelve nonzeros in total. ## The right-hand side The right-hand side is a scalar, applied to every row, or a parameter read at exactly the frame of the constraint, giving each row its own value. The parameter is read at its sets here as it is anywhere else: `supply[P]`, not `supply`. A parameter without a bracket raises `TypeError` and reports the reading it requires. Supply is indexed over `P`, and so are the supply rows. A parameter read over any other index set raises `ValueError`, and the message gives both index sets. ```python raises=TypeError import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply) ``` Raises TypeError: ```text TypeError: parameter 'supply' is over ('P',) and expresses no coefficient until it is read; read it at its sets as supply[P] ``` ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= demand[W]) ``` Raises ValueError: ```text ValueError: constraint 'supply' has free dimensions ('P',); its right-hand side 'demand' is over ('W',) ``` ## One bound per constraint Python evaluates the chained comparison `0 <= expr <= 10` as `(0 <= expr) and (expr <= 10)` and discards the first relation. `nimopt` raises `TypeError` on the chained form and drops no bound. Each bound is written as its own constraint. ```python raises=TypeError import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) m = Model("transport") x = m.var("x", (P, W)) 0.0 <= Sum(W, x[P, W]) <= 10.0 ``` Raises TypeError: ```text TypeError: a relation has no truth value; write each bound in its own constraint ``` Next: [Solving](/tutorial/solving). --- # /tutorial/solving # Solving The remaining pieces are the objective function and the solver call. ## The objective `m.set_objective(expression)` takes an expression with an empty frame. The model's `sense`, `"min"` by default, sets the direction. `m.solve()` assembles the matrix, calls HiGHS and returns a `Solution`. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(solution.status) print(solution.objective) ``` Output: ```text optimal 135.0 ``` The status is `optimal` and the objective value is 135. ## Status `status` reports the outcome of the solve. `status` and `feasible` are readable after any solve. `objective` and `primal` raise `ValueError` where `feasible` is False. They raise at status `unbounded` and `unbounded_or_infeasible` whatever `feasible` reports. `bound` and `gap` are `None` at those two statuses. `dual` raises `ValueError` where `status` is not `optimal`. Raising the demand of Berlin to 40 makes total demand 70 against total supply 55. The model is infeasible. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([40.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(solution.status) solution.objective ``` Raises ValueError: ```text infeasible ValueError: status is 'infeasible' and the solver reports no feasible point; read `status` before reading values ``` For an infeasible model whose cause is not evident, `m.session()` keeps the solver instance open and `diagnose()` returns the conflicting rows. See [Solvers](/reference/solvers). ## A solve stopped at a limit An option in `options()` stops the solver early. A solver stopped at a limit reports the best point it found, and `feasible` is True for it. `objective` and `primal` then return that point. `bound` returns what the solver proved about the optimum, an upper bound under sense `max` and a lower bound under sense `min`. `gap` returns the relative distance from the objective to that bound. A search stopped after one node returns the point the solver found there. The thread count of the solver determines that point. The example below reports the properties every such point has. ```python import numpy as np from nimopt import Model, Param, Set, Sum rng = np.random.default_rng(1) ITEM = Set("item", np.array([f"i{t}" for t in range(40)])) BIN = Set("bin", np.array([f"b{t}" for t in range(5)])) weight = Param.from_dense("weight", (BIN, ITEM), rng.uniform(1, 50, (5, 40))) value = Param.from_dense("value", (ITEM,), rng.uniform(1, 100, 40)) capacity = Param.from_dense("capacity", (BIN,), np.full(5, 306.0)) m = Model("knapsack", sense="max") x = m.var("x", (ITEM,), integer=True, upper=1.0) m.constraint("capacity", Sum(ITEM, weight[BIN, ITEM] * x[ITEM]) <= capacity[BIN]) m.set_objective(Sum(ITEM, value[ITEM] * x[ITEM])) solution = m.solve(options={"node_limit": 1}) print(f"status: {solution.status}") print(f"feasible: {solution.feasible}") print(f"the bound is above the objective: {solution.bound > solution.objective}") print(f"the gap is positive: {solution.gap > 0.0}") ``` Output: ```text status: solution_limit feasible: True the bound is above the objective: True the gap is positive: True ``` The status identifies the limit the solver stopped at. HiGHS reports a stop at `node_limit` as `solution_limit`. `gap` is `None` where the solver proved no bound, and `feasible` is False where it found no point. ## The matrix `m.assemble()` builds the coefficient matrix without a solver call and returns it in CSR form. `to_dense()` renders it for a model of this size, and `row_of(name)` gives the row range of a named constraint. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) assembled = m.assemble() print(m.n_rows, m.n_columns, m.nnz) print(assembled.row_of("demand")) print(assembled.to_dense()) ``` Output: ```text 5 6 12 slice(2, 5, None) [[1. 1. 1. 0. 0. 0.] [0. 0. 0. 1. 1. 1.] [1. 0. 0. 1. 0. 0.] [0. 1. 0. 0. 1. 0.] [0. 0. 1. 0. 0. 1.]] ``` The six columns are the routes, Lisbon's three followed by Porto's. Rows 0 and 1 are the supply rows, each with a 1 under its plant's three routes. Rows 2 to 4 are the demand rows, each with a 1 under the two routes into its warehouse. `row_of("demand")` returns that range. Next: [Reading the solution](/tutorial/reading-the-answer). --- # /tutorial/reading-the-answer # Reading the solution A solver returns primal and dual values as flat vectors. `nimopt` returns them as arrays over the sets each variable and constraint was declared on. ## Primals and duals `solution.primal(name)` returns a variable's values over its sets. `solution.dual(name)` returns a constraint's dual values over its frame: for the demand constraint, one value per warehouse. `dual` also takes a variable, and returns its reduced costs: the cost of the variable less the duals of the rows it appears in. A route with no shipment at the optimum reports the amount by which its cost exceeds those duals. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() shipped = solution.primal("x") print(shipped.dims) print(shipped.to_dense()) print(solution.dual("demand").dims) print(solution.dual("demand").to_dense()) ``` Output: ```text ('P', 'W') [[20. 0. 10.] [ 0. 15. 5.]] ('W',) [3. 1. 6.] ``` Rows are plants and columns are warehouses: Lisbon ships 20 to Berlin and 10 to Rome, Porto ships 15 to Paris and 5 to Rome. The dual of a demand row is the increase in total cost per additional unit of demand at that warehouse: 3 in Berlin, 1 in Paris, 6 in Rome. The supply constraint of Lisbon binds, and each marginal unit is served from Porto at the cost of the route from Porto. ## A variable over a subset If Porto cannot ship to Rome, the variable is declared over the five existing routes with `subset=`. The sixth route has no column, and the solution has values at the five members only. ```python import numpy as np from nimopt import Model, Param, Set, Sum, subset P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) routes = subset( (P, W), { "P": np.array(["lisbon", "lisbon", "lisbon", "porto", "porto"]), "W": np.array(["berlin", "paris", "rome", "berlin", "paris"]), }, ) m = Model("transport") x = m.var("x", (P, W), subset=routes) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() print(x.n_columns) print(solution.objective) print(type(solution.primal("x")).__name__) ``` Output: ```text 5 135.0 SparseArray ``` Five columns instead of six. The objective is unchanged at 135. Rome is served from Lisbon in both solutions. The array type differs: a variable over a full product returns a `DenseArray`, a variable over a subset a `SparseArray` with an entry per member and nothing elsewhere. ## Absence is not zero The model contains no decision for the route Porto to Rome. Every array a solution returns declares `absence="unknown"`. `to_dense()` raises `ValueError`. It supplies no value for a coordinate the model does not produce. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set, Sum, subset P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) routes = subset( (P, W), { "P": np.array(["lisbon", "lisbon", "lisbon", "porto", "porto"]), "W": np.array(["berlin", "paris", "rome", "berlin", "paris"]), }, ) m = Model("transport") x = m.var("x", (P, W), subset=routes) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() solution.primal("x").to_dense() ``` Raises ValueError: ```text ValueError: absence is 'unknown' and the array has no value at 1 of 6 coordinates; pass fill= to to_dense() ``` `to_dense(fill=...)` supplies the value for missing coordinates. `nan` distinguishes a missing route from a route with zero shipment. ```python import numpy as np from nimopt import Model, Param, Set, Sum, subset P = Set("P", np.array(["lisbon", "porto"])) W = Set("W", np.array(["berlin", "paris", "rome"])) cost = Param.from_dense("cost", (P, W), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) supply = Param.from_dense("supply", (P,), np.array([30.0, 25.0])) demand = Param.from_dense("demand", (W,), np.array([20.0, 15.0, 15.0])) routes = subset( (P, W), { "P": np.array(["lisbon", "lisbon", "lisbon", "porto", "porto"]), "W": np.array(["berlin", "paris", "rome", "berlin", "paris"]), }, ) m = Model("transport") x = m.var("x", (P, W), subset=routes) m.constraint("supply", Sum(W, x[P, W]) <= supply[P]) m.constraint("demand", Sum(P, x[P, W]) >= demand[W]) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) solution = m.solve() values = solution.primal("x") print(values.absence, values.nnz) print(values.to_dense(fill=np.nan)) ``` Output: ```text unknown 5 [[15. 0. 15.] [ 5. 15. nan]] ``` Lisbon serves Rome alone, and the Porto to Rome cell reads `nan`. A stored `0.0` would denote a route that exists and ships nothing. ## Summary The tutorial declared index sets and parameters, a decision variable, expressions, two constraint families and an objective. It solved the model and read the solution back over its sets. The [guides](/guides/subsets) cover variables over sparse networks, conditions on rows, lags, bounds from data, and models with millions of rows. The [worked models](/models) present ten complete formulations. --- # /guides/at-scale # Declaring a model at scale Models with millions of columns are declared the same way as small ones. A declaration computes shapes, not blocks. A constraint computes its row and nonzero counts when it is added, and no matrix exists until `assemble()` or `solve()` is called. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array([f"p{i}" for i in range(200)])) W = Set("W", np.array([f"w{i}" for i in range(100)])) cost = Param.from_dense("cost", (P, W), np.ones((200, 100))) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= 1.0) m.constraint("demand", Sum(P, x[P, W]) >= 1.0) m.set_objective(Sum(P, W, cost[P, W] * x[P, W])) print(m.n_columns) print(m.n_rows) print(m.nnz) ``` Output: ```text 20000 300 40000 ``` Twenty thousand columns, three hundred rows and forty thousand nonzeros are declared, and the model contains no matrix. The parameter is the data of the caller. The model has added two constraints and an objective, each a symbolic expression recording the variable, the coefficient and the sets that are summed. ## Building the matrix once `assemble()` builds the matrix into one buffer and returns it in CSR form. Each constraint writes its rows into its own slice. The matrix therefore exists once, and one expression at a time is materialised. ```python import numpy as np from nimopt import Model, Set, Sum P = Set("P", np.array([f"p{i}" for i in range(200)])) W = Set("W", np.array([f"w{i}" for i in range(100)])) m = Model("transport") x = m.var("x", (P, W)) m.constraint("supply", Sum(W, x[P, W]) <= 1.0) m.constraint("demand", Sum(P, x[P, W]) >= 1.0) assembled = m.assemble() print(assembled.n_rows, assembled.n_cols) print(assembled.values.shape) print(assembled.row_of("demand")) ``` Output: ```text 300 20000 (40000,) slice(200, 300, None) ``` `indices` and `values` are views of that buffer; only `indptr` is built. `row_of(name)` gives the rows of a constraint as a slice. A caller locates a named row block in a matrix too large to print. ## Reading a large solution `to_dense()` is for a model small enough to print. At scale, read the CSR arrays, or read the solution through `primal()` and `dual()`. Both return labeled arrays over the sets, not offsets into a vector. A variable over a subset keeps its values sparse. Reading a solution over a sparse network therefore builds no grid. ## What a variable costs A member's column is computed from its multi-index and is not stored. A variable over millions of columns costs its members, not its columns. `n_columns` above is twenty thousand, and the variable stores a few numbers: the start of its block and the sizes of its sets. ## Progress A large model takes seconds to build. `progress=True` reports the build in a terminal. ```python skip="the report draws to a terminal, and a page is not one" from nimopt.models import storage model = storage.definition().build(storage.data(60), progress=True) solution = model.solve(options={"time_limit": 300.0, "log": True}, progress=True) ``` A solve stopped at `time_limit` reports `feasible` True where the solver found a point. `objective`, `primal`, `bound` and `gap` then read that point. See [Solving](/tutorial/solving). --- # /guides/bounds-from-parameters # Bounds from a parameter A capacity per generator, an energy limit per battery, a flow limit per line: bounds usually come from data and differ by member. `lower=` and `upper=` accept a number, applied to every column, or a `Param`, giving each member its own value. ```python import numpy as np from nimopt import Model, Param, Set G = Set("G", np.array(["g0", "g1"])) cap = Param.from_dense("cap", (G,), np.array([5.0, 7.0])) m = Model("schedule") m.var("x", (G,), upper=cap) lower, upper = m.column_bounds() print(lower) print(upper) ``` Output: ```text [0. 0.] [5. 7.] ``` ## Broadcasting a narrower parameter A parameter indexed over fewer dimensions than the variable is broadcast over the rest. A capacity per unit bounds every period of that unit. ```python import numpy as np from nimopt import Model, Param, Set G = Set("G", np.array(["g0", "g1"])) T = Set("T", np.array(["t0", "t1", "t2"])) cap = Param.from_dense("cap", (G,), np.array([5.0, 7.0])) m = Model("schedule") m.var("x", (G, T), upper=cap) print(m.column_bounds()[1]) ``` Output: ```text [5. 5. 5. 7. 7. 7.] ``` Six columns. The three periods of each unit take the capacity of that unit. ## A bound has a value at every column A bound with no value at some member of the variable raises `ValueError`, and the message gives the member. A dense parameter covers its product by construction. A long-form parameter can omit a member. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set S = Set("S", np.array(["a", "b", "c"])) cap = Param.from_long("cap", (S,), {"S": np.array(["a", "c"])}, np.array([1.0, 3.0])) m = Model("bounds") m.var("x", (S,), upper=cap) m.column_bounds() ``` Raises ValueError: ```text ValueError: the upper bound 'cap' has no value at member ('b',) of variable 'x'; give the bound a value at every member of the variable ``` The check runs when the bound vectors are built. `column_bounds()`, `assemble()` and `solve()` all build them. ## A bound over a dimension the variable lacks A parameter indexed over a dimension the variable is not declared on raises `ValueError`. A bound is per column, and a dimension the variable does not have selects no column. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set G = Set("G", np.array(["g0", "g1"])) W = Set("W", np.array(["w0", "w1"])) cap = Param.from_dense("cap", (W,), np.array([5.0, 7.0])) m = Model("schedule") m.var("x", (G,), upper=cap) ``` Raises ValueError: ```text ValueError: variable 'x' is declared over ('G',) and is not over ['W']; its upper bound 'cap' is declared over ('W',) ``` --- # /guides/coefficient-arithmetic # Coefficient arithmetic A coefficient is often derived from several parameters: fuel price divided by efficiency, a cost scaled by a factor. In `nimopt` a coefficient is a parameter read at its sets or an arithmetic combination of such readings. `+`, `-`, `*`, `/` and a power by a number combine them. The combination is symbolic: it contains references, derives its dimensions from its operands, and is evaluated once, when the term it multiplies is materialised. A derived coefficient can therefore appear in a definition before any data exists. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["base", "peak"])) T = Set("T", np.arange(3)) price = Param.from_dense("fuel_price", (G, T), np.full((2, 3), 30.0)) eta = Param.from_dense("efficiency", (G, T), np.array([[0.5] * 3, [0.4] * 3])) load = Param.from_dense("load", (T,), np.full(3, 100.0)) cap = Param.from_dense("capacity", (G, T), np.full((2, 3), 80.0)) unit_cost = price[G, T] / eta[G, T] m = Model("dispatch", sense="min") gen = m.var("gen", (G, T), lower=0.0, upper=cap) m.constraint("balance", Sum(G, gen[G, T]) == load[T]) m.set_objective(Sum(G, T, unit_cost[G, T] * gen[G, T])) print(unit_cost.name, unit_cost.dims) print(m.solve().objective) ``` Output: ```text (fuel_price / efficiency) ('G', 'T') 18900.0 ``` ## Reading a derived coefficient `unit_cost[G, T]` reads a combination the same way `price[G, T]` reads a parameter, and the sets given are checked against the combination's dimensions. A transposed or incomplete index raises `ValueError`. ```python raises=ValueError import numpy as np from nimopt import Param, Set G = Set("G", np.array(["base", "peak"])) T = Set("T", np.arange(3)) price = Param.from_dense("fuel_price", (G, T), np.full((2, 3), 30.0)) eta = Param.from_dense("efficiency", (G, T), np.array([[0.5] * 3, [0.4] * 3])) (price[G, T] / eta[G, T])[T, G] ``` Raises ValueError: ```text ValueError: coefficient (fuel_price / efficiency) is declared over ('G', 'T'); got ('T', 'G') ``` A bare parameter has no arithmetic: `price * 2.0` raises `TypeError`. Read the parameter at its sets first and combine the references. ## Alignment Two operands with the same dimensions align entry by entry. Operands whose dimensions nest or overlap align on the shared dimensions and broadcast over the rest, with the left operand's order first. Operands sharing no dimension raise `ValueError`. Their product would be an outer product, and a linear model does not require one. The same rule applies when a coefficient multiplies a variable. A number has no dimensions and scales. ```python raises=ValueError import numpy as np from nimopt import Param, Set G = Set("G", np.array(["base", "peak"])) T = Set("T", np.arange(3)) over_g = Param.from_dense("over_g", (G,), np.ones(2)) over_t = Param.from_dense("over_t", (T,), np.ones(3)) over_g[G] * over_t[T] ``` Raises ValueError: ```text ValueError: frames ('G',) and ('T',) share no dimension; pass operands that share a dimension ``` **A coefficient never adds a column the variable does not have.** The variable's members define the columns; a coefficient can only reduce which of them receive a nonzero. A coefficient with dimensions the variable lacks defines rows over those dimensions: this is how a term maps rows to columns. ## Division by zero A divisor that is zero raises `ZeroDivisionError`, and the message gives the coordinate. A quotient of infinity is not passed to a solver. The check covers a Python number, a NumPy scalar and a coefficient with a zero at any coordinate. A NumPy scalar divides to infinity where a Python number raises. ```python raises=ZeroDivisionError import numpy as np from nimopt import Param, Set G = Set("G", np.array(["base", "peak"])) T = Set("T", np.arange(3)) price = Param.from_dense("fuel_price", (G, T), np.full((2, 3), 30.0)) holed = Param.from_dense("holed", (G, T), np.array([[0.5, 0.0, 0.5], [0.4] * 3])) price[G, T] / holed[G, T] ``` Raises ZeroDivisionError: ```text ZeroDivisionError: divisor holed is zero at 1 coordinate(s), first at {'G': 'base', 'T': 1}; remove the zeros or divide by another parameter ``` ## Constants An expression consists of terms and a constant. `x + 1 <= 5` produces the row `x <= 4`: the constant moves to the right-hand side when the constraint is built, and into the reported objective value as a fixed cost. A constant adds no column. A constant alone is neither an objective nor a constraint, and raises. ```python import numpy as np from nimopt import Model, Param, Set, Sum P = Set("P", np.array(["a"])) one = Param.from_dense("one", (P,), np.ones(1)) m = Model("m", sense="max") x = m.var("x", (P,), upper=100.0) m.constraint("cap", one[P] * x[P] + 1.0 <= 5.0) m.set_objective(Sum(P, one[P] * x[P]) + 7.0) print(m.assemble().n_cols, m.assemble().row_upper) print(m.solve().objective) ``` Output: ```text 1 [4.] 11.0 ``` ## Forms that raise A line of modeling arithmetic produces a linear term, or raises with a message that gives the form to write instead. There is no third outcome. | Written | Form to write instead | | --- | --- | | `x[P] * y[P]` | a linear term has one variable; a coefficient multiplies it | | `x[P] ** 2` | a coefficient takes the power, and a variable multiplies it | | `x[P] / y[P]` | a variable in a denominator is not linear | | `2.0 / x[P]` | the same: write the reciprocal as a coefficient | | `x[P] / 0.0` | a divisor of zero is handled before it is passed to an expression | | `x[P] < 1.0` | `<=` and `>=`; an LP has no row for a strict inequality | | `x[P] > 1.0` | the same | | `x[P] != 1.0` | one bound per constraint | | `0.0 <= x[P] <= 1.0` | each bound as its own constraint | | `x[P] + c[P]` | a coefficient has no row until a variable multiplies it | | `c[P] * (x[P] + 1.0)` | a coefficient times a constant is one value per row; the constant goes to the right-hand side | | `abs(x[P])`, `min` and `max` | reduce with `Sum`, or bound the expression with two rows | | `np.sum(x[P])` | `Sum` and its sets | | `np.array([...]) * x[P]` | `Param.from_dense`, read at its sets | | `Sum(P, Sum(P, x[P]))` | a dimension is reduced once | | `Sum(P - 1, x[P])` | `Sum(P, x[P - 1])`: the lag belongs on the reference | | `T - 1.7` | a lag is a whole number of members | Each raises with a `nimopt` message, not with a bare Python error. The test suite executes the table and the replacement forms the messages give. ```python raises=TypeError import numpy as np from nimopt import Model, Set P = Set("P", np.array(["a", "b"])) m = Model("m") x = m.var("x", (P,)) np.array([1.0, 2.0]) * x[P] ``` Raises TypeError: ```text TypeError: a coefficient is a parameter; build one with `Param.from_dense` or `Param.from_long` and read it at its sets. A product of two expressions is not linear. ``` A NumPy array multiplying a term would let NumPy apply the operator and return an array of expressions. The expression types reject the ufunc, and the message describes how a coefficient is built. --- # /guides/conditions # Conditions on a sum and on a constraint Two cases require a condition. A constraint may sum over part of the members of a variable, such as the arcs of a network where the variable is indexed over the full product. A constraint may also apply to some members of its frame only, such as a capacity limit on one plant. `where=` covers both cases. `over=` declares the rows of a constraint explicitly. ## Restricting a sum `Sum(..., where=domain)` restricts each term to the members of `domain` before summing. ```python import numpy as np from nimopt import Model, Set, Sum, subset P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) arcs = subset( (P, W), {"P": np.array(["p0", "p0", "p1"]), "W": np.array(["w0", "w1", "w2"])}, ) m = Model("network") x = m.var("x", (P, W)) rows = m.constraint("capacity", Sum(W, x[P, W], where=arcs) <= 10.0) print(x.n_columns) print(rows.n_rows, rows.nnz) print(m.assemble().to_dense()) ``` Output: ```text 6 2 3 [[1. 1. 0. 0. 0. 0.] [0. 0. 0. 0. 0. 1.]] ``` The variable is over the full product and has six columns. The condition gives three of them a coefficient. A variable over a subset would have three columns from the start. Use a condition where the variable is over the product and one constraint reads part of it. Declare a subset where the model never uses the other members. ## Restricting the rows `m.constraint(..., where=domain)` takes a domain over the constraint's frame and keeps the rows in it. A row outside the condition is not produced. ```python import numpy as np from nimopt import Model, Set, Sum, subset P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) m = Model("network") x = m.var("x", (P, W)) only_p0 = subset((P,), {"P": np.array(["p0"])}) rows = m.constraint("capacity", Sum(W, x[P, W]) <= 10.0, where=only_p0) print(rows.n_rows) print(m.assemble().to_dense()) ``` Output: ```text 1 [[1. 1. 1. 0. 0. 0.]] ``` One row, for `p0`. `p1` has no capacity row. A condition over dimensions other than the frame of the constraint raises `ValueError`, and the message gives both index sets. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum, subset P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) m = Model("network") x = m.var("x", (P, W)) by_warehouse = subset((W,), {"W": np.array(["w0"])}) m.constraint("capacity", Sum(W, x[P, W]) <= 10.0, where=by_warehouse) ``` Raises ValueError: ```text ValueError: constraint 'capacity' has free dimensions ('P',); its condition is over ('W',) ``` ## Declaring the rows explicitly By default the rows of a constraint are derived from its terms: a row exists where every term has a value and the right-hand side has a value. A term with no value along a frame dimension removes the row. A row missing one of its terms would express a constraint that was not written. `over=domain` declares the rows instead of deriving them. A term with values at some of the rows contributes where it has them, and every row in the domain is produced. ```python import numpy as np from nimopt import Model, Set, Sum, product P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) m = Model("network") x = m.var("x", (P, W)) rows = m.constraint("capacity", Sum(W, x[P, W]) <= 10.0, over=product((P,))) print(rows.n_rows) ``` Output: ```text 2 ``` ## `over` or `where`, not both `over=` declares the rows and `where=` restricts them. Passing both raises `ValueError`. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum, product, subset P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) m = Model("network") x = m.var("x", (P, W)) m.constraint( "capacity", Sum(W, x[P, W]) <= 10.0, where=subset((P,), {"P": np.array(["p0"])}), over=product((P,)), ) ``` Raises ValueError: ```text ValueError: constraint 'capacity' is given over= and where= together; pass one of them ``` A condition on a sum and a condition on the constraint compose. The first restricts what is summed, and the second restricts which rows exist. --- # /guides/fixed-members # A member fixed at a label Initial conditions, terminal conditions and boundary rows reference one member of a set, such as the stored energy at the first period or the level at the last. A label in place of a set in a reference fixes that dimension at one member and removes it from the frame. ```python import numpy as np from nimopt import Model, Set G = Set("G", np.array(["g0", "g1"])) T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (G, T)) print(x[G, T].frame) print(x[G, "t0"].frame) ``` Output: ```text ('G', 'T') ('G',) ``` `T` is fixed at `t0`, and the reference is indexed over `G` alone. An initial condition is one row per unit, referencing the column of that unit at the first period. ```python import numpy as np from nimopt import Model, Set G = Set("G", np.array(["g0", "g1"])) T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (G, T)) rows = m.constraint("start", x[G, "t0"] <= 1.0) print(rows.n_rows, rows.nnz) print(m.assemble().to_dense()) ``` Output: ```text 2 2 [[1. 0. 0. 0. 0. 0.] [0. 0. 0. 1. 0. 0.]] ``` Two rows over six columns, each with one nonzero: the `t0` column of its own unit. ## A coefficient at a member A parameter takes a label the same way and yields the coefficients at that member. ```python import numpy as np from nimopt import Model, Param, Set G = Set("G", np.array(["g0", "g1"])) T = Set("T", np.array(["t0", "t1", "t2"])) rate = Param.from_dense("rate", (G, T), np.array([[2.0, 4.0, 5.0], [3.0, 1.0, 6.0]])) m = Model("schedule") x = m.var("x", (G, T)) rows = m.constraint("start", rate[G, "t0"] * x[G, "t0"] <= 1.0) print(rate[G, "t0"].dims) print(m.assemble().to_dense()) ``` Output: ```text ('G',) [[2. 0. 0. 0. 0. 0.] [0. 0. 0. 3. 0. 0.]] ``` The coefficients are the `t0` column of `rate`: 2.0 and 3.0. ## The label must be a member A label that is not a member of the set raises `ValueError`; the message gives the label and the dimension. ```python raises=ValueError import numpy as np from nimopt import Model, Set G = Set("G", np.array(["g0", "g1"])) T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (G, T)) x[G, "t9"] ``` Raises ValueError: ```text ValueError: variable 'x' is read at member 't9' of dimension 'T'; read it at a member that set contains ``` --- # /guides/highs-methods # Interior point and first-order methods `method=` selects the algorithm HiGHS runs on the assembled matrix. Four values select one algorithm each, and three further options configure the interior point and first-order methods. Every setting on this page runs in the test suite against the bundled models, and each solve returns the optimum of the model. Which method is fastest or smallest depends on the model. | Option | Choices | Applies to | | --- | --- | --- | | `method` | `choose`, `simplex`, `barrier`, `hipo`, `pdlp` | every solve | | `newton_system` | `choose`, `augmented`, `normaleq` | `hipo` | | `crossover` | `choose`, `off`, `on` | `barrier` and `hipo` | | `pdlp_tol` | a relative tolerance | `pdlp` | `barrier` runs IPX, the HiGHS interior point method built on a preconditioned conjugate gradient. `hipo` runs HiPO, an interior point method built on a direct factorization of the Newton system, parallel across the elimination tree. `newton_system` selects the augmented system or the normal equations, and `choose` leaves the selection to the solver. HiPO stores that factorization in memory and requires more memory than IPX on the same model. `threads` applies to its iterations and not to the crossover that may follow them. `pdlp` runs cuPDLP-C, a primal-dual hybrid gradient method that reads the matrix only through matrix-vector products. `crossover` controls whether an interior point method passes its solution to the simplex method to obtain a vertex. ```python skip="needs a HiGHS built with HiPO" solution = model.solve( options={ "method": "hipo", "newton_system": "augmented", "crossover": "off", "threads": 8, } ) ``` Gurobi and Mosek support `crossover` and `method` up to `barrier`. `newton_system`, `pdlp_tol`, `hipo` and `pdlp` are specific to HiGHS. Passing one of them to another solver raises and reports the option. The [solvers reference](/reference/solvers) records the same rule. ## Crossover An interior point method stops at a point inside the feasible region, within tolerance of the optimum on every constraint. Crossover moves that point to a vertex with the simplex method. A basis, an exact active set and duals at a vertex require that move. Crossover runs serially, and on a large model it can cost more than the interior point iterations before it. `crossover="off"` returns the interior point as the solution. Primals and duals are read the same way. A constraint that is tight at the optimum may be a tolerance away from equality, and a variable at a bound may be a tolerance inside it. The interior point is sufficient where the solution is read as quantities and prices and not as a basis. ## PDLP on a GPU A first-order method computes no factorization. Its memory is the matrix in two orientations, one for each product, plus working vectors of the row and column dimensions. That memory grows linearly with the problem, and a problem too large for a factorization still fits. **Where the memory goes.** HiGHS presolves on the CPU, in host memory, before any method runs: the original LP and its reduced form are both stored there while presolve runs, and that is the peak on the host. The reduced problem is transferred to the GPU. It is smaller than the problem as declared: presolve removes the rows that bound a single column and the columns it can fix. The GPU therefore needs memory for the reduced matrix twice, the working vectors, and the buffers of the sparse kernels, and nothing else. The host keeps the original and the reduced LP throughout, and postsolve maps the solution back through them. A HiGHS built without CUDA runs the same method on the CPU. The solver log, requested with `log=True`, reports the device the method runs on. **The tolerance.** `pdlp_tol` is the relative tolerance at which PDLP stops: the duality gap and the primal and dual residuals, each relative to the scale of the problem, must all fall below it. A looser tolerance stops sooner and returns a point farther from the optimum and from feasibility. A tighter tolerance costs more iterations, and each iteration is a pass over the matrix. The iterations are cheap and numerous, and the tolerance controls both the solve time and the accuracy of the solution. HiGHS checks the point PDLP returns against its own `feasibility_tol` and `optimality_tol` after postsolve. That check is stricter than the criterion of PDLP. PDLP measures its residuals relative to the scaled problem it iterates on. Postsolve maps the point back onto the original rows, where a residual that passed can exceed the tolerance. A point that meets `pdlp_tol` and fails the HiGHS check reports the status HiGHS calls unknown. `nimopt` raises on that status and reports it, and it reads no values from that point. On a model of any size the default `pdlp_tol` is usually too loose for HiGHS to accept the point. The two agree at a tolerance one to two orders tighter. Tighten `pdlp_tol` until the point passes, or loosen the feasibility and optimality tolerances to the accuracy the model requires. PDLP reports two iterates in its log, the running average marked `[A]` and the last marked `[L]`, and stops on whichever meets the tolerance first. Progress is not monotone: the gap can close and open again as the method restarts. The method returns the last row of the log, not the best row. ```python skip="needs a HiGHS built with CUDA to run on a GPU" solution = model.solve(options={"method": "pdlp", "pdlp_tol": 1e-8, "log": True}) ``` ## Installing a HiGHS with HiPO and GPU support HiGHS 1.15.1 stores the HiPO orderings and its BLAS in a separate library, `libhighs_extras`. `libhighs` loads that library at run time by name. The `highspy` wheel on PyPI and the conda-forge package ship without that library and without CUDA. Given `hipo`, such a HiGHS logs an error and runs simplex. `nimopt` raises instead and reports the missing library. Given `pdlp`, it runs the method on the CPU. Both come from a source build of the same version, with three pieces installed separately: 1. **The extras library**, from the `extern` directory of the HiGHS repository. That directory is a CMake project of its own. It builds METIS, AMD and RCM from the tree and links a BLAS. OpenBLAS from a conda environment is sufficient. ```bash cmake -S extern -B build-extras -G Ninja -DCMAKE_BUILD_TYPE=Release \ -DHIPO=ON -DBLA_VENDOR=OpenBLAS cmake --build build-extras ``` 2. **The wheel**, from the repository root, where the HiGHS build configuration enables HiPO. Pass the prefix of the conda environment through the `CMAKE_PREFIX_PATH` environment variable, not through `CMAKE_ARGS`. `CMAKE_ARGS` replaces the path the Python build adds for pybind11. ```bash CMAKE_PREFIX_PATH=$CONDA_PREFIX CMAKE_ARGS="-DBLA_VENDOR=OpenBLAS" \ uv build --wheel --python /bin/python -o dist . ``` For the GPU, the same command with a CUDA toolkit in the environment and the card's compute capability: ```bash CMAKE_PREFIX_PATH=$CONDA_PREFIX CUDACXX=$CONDA_PREFIX/bin/nvcc CUDAToolkit_ROOT=$CONDA_PREFIX \ CMAKE_ARGS="-DCUPDLP_GPU=ON -DCMAKE_CUDA_ARCHITECTURES=75 -DBLA_VENDOR=OpenBLAS" \ uv build --wheel --python /bin/python -o dist . ``` 3. **The placement.** The `libhighs` of the wheel searches its own directory for the extras library. Copy the built `libhighs_extras.so` into the `highspy` package directory of the environment the wheel is installed in. The extras library links the BLAS it was built against, and the GPU wheel links the CUDA runtime, cuBLAS and cuSPARSE. The conda environment that provides them must remain installed. The `libhighs` and `libcudalin` of the GPU wheel require the `lib` directory of that environment in their own run path. `patchelf --set-rpath` sets it on the unpacked wheel before the wheel is packed again. The run path of the extension module alone does not resolve the CUDA libraries. ```bash uv pip install --python /bin/python dist/highspy-1.15.1-*.whl cp build-extras/libhighs_extras.so /lib/python3.13/site-packages/highspy/ ``` A solve with `method="hipo"` verifies the installation. A HiGHS without the extras library raises before anything is solved. A HiGHS with the library returns the optimum. --- # /guides/lags # Lags Time-coupled constraints reference the previous period. A storage balance relates the stored energy at `t` to the stored energy at `t-1`. A ramp limit bounds the change in output between consecutive periods. `T - 1` is the set `T` lagged by one member, and `x[T - 1]` references the variable at the previous member. ## A lag that drops the boundary row ```python import numpy as np from nimopt import Model, Set T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (T,)) rows = m.constraint("carry", x[T] - x[T - 1] <= 0.0) print(rows.n_rows, rows.nnz) print(m.assemble().to_dense()) ``` Output: ```text 2 4 [[-1. 1. 0.] [ 0. -1. 1.]] ``` Each row references its own column and the previous one. The first member has no predecessor, and its row is not produced. Three members give two rows. `T + 1` references the following member by the same rule. ## A lag that wraps `T.cyclic` lags with wrap-around: the member before the first is the last. No row is dropped. ```python import numpy as np from nimopt import Model, Set T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (T,)) rows = m.constraint("carry", x[T] - x[T.cyclic - 1] <= 0.0) print(rows.n_rows, rows.nnz) print(m.assemble().to_dense()) ``` Output: ```text 3 6 [[ 1. 0. -1.] [-1. 1. 0.] [ 0. -1. 1.]] ``` Three rows, and the first references the last column: the `-1` in row 0 is in the final position. A storage balance over a repeating horizon is written this way, and the level at the end of the horizon enters the row of the first period. ## A lag applies to a reference, not to a sum `Sum` runs over the members of a set and takes the set itself. Passing a lagged set raises `ValueError`. ```python raises=ValueError import numpy as np from nimopt import Model, Set, Sum T = Set("T", np.array(["t0", "t1", "t2"])) m = Model("schedule") x = m.var("x", (T,)) Sum(T - 1, x[T]) ``` Raises ValueError: ```text ValueError: a sum is over the members of ['T'] and takes the set, not a lag of it; write the lag at the variable's reference ``` The lag belongs on the variable reference: `Sum(T, x[T - 1])`. ## A lag applies to a variable, not to a parameter Reading a parameter at a lag raises `ValueError`. A coefficient is indexed by the row it appears in, and a lag selects which column a row references. `rate[T] * x[T - 1]` applies the rate at `t` to the variable at `t-1`. ```python raises=ValueError import numpy as np from nimopt import Param, Set T = Set("T", np.array(["t0", "t1", "t2"])) rate = Param.from_dense("rate", (T,), np.array([1.0, 2.0, 3.0])) rate[T - 1] ``` Raises ValueError: ```text ValueError: parameter 'rate' is read at a lag ['T']; write the lag at the variable's reference ``` A lag is an integer number of members. A fractional lag raises `ValueError` and is not truncated to a different lag. ```python raises=ValueError import numpy as np from nimopt import Set T = Set("T", np.arange(3)) T - 1.7 ``` Raises ValueError: ```text ValueError: a lag is a whole number of members; got 1.7 ``` --- # /guides/piecewise # Piecewise-linear curves A fuel cost that rises in steps, an efficiency that changes with load, and a revenue that saturates are each a curve through a list of points. `piecewise` relates an expression `y` to an expression `x` through such points. `x` is on the curve. `sign` compares `y` with the curve: `==` sets `y` to the curve, `>=` bounds `y` below by it and `<=` bounds `y` above by it. The points are two parameters, `x_points` and `y_points`. Each is over some or all of the sets of `x` and over one breakpoint set. The declaration generates the variables and the constraints of its method. Each generated name begins with the declaration's name. ## A cost curve that is not convex `method="incremental"` is exact for breakpoints that are strictly increasing or strictly decreasing. It adds one continuous variable and one integer variable per segment. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["a", "b"])) B = Set("B", np.array(["b0", "b1", "b2", "b3"])) power = Param.from_dense("power", (G, B), [[0, 10, 20, 30], [0, 10, 20, 30]]) cost = Param.from_dense("cost", (G, B), [[0, 5, 30, 35], [0, 20, 30, 40]]) m = Model("curve") p = m.var("p", (G,)) c = m.var("c", (G,)) m.constraint("demand", Sum(G, p[G]) == 25.0) m.piecewise( "fuel", x=p[G], x_points=power[G, B], y=c[G], y_points=cost[G, B], sign=">=", method="incremental", ) m.set_objective(Sum(G, c[G])) s = m.solve() print(s.objective) print(s.primal("p").to_dense()) print(list(m.variables)) print(list(m.constraints)) ``` Output: ```text 30.0 [10. 15.] ['p', 'c', 'fuel_fill', 'fuel_order'] ['demand', 'fuel_x', 'fuel_y', 'fuel_order_bound', 'fuel_fill_order', 'fuel_order_link'] ``` The generated variables and constraints are declarations of the model. `primal`, `dual`, `row` and `explain` read them by their names. `x` and `y` are expressions. `p[G] + 5.0` is on the curve five units above `p`, and the incremental method reads that constant. The tangent method multiplies `x` by the slope of each segment, and a constant in `x` raises `ValueError`. Subtract it from `x_points` instead. ## A convex curve `method="tangent"` adds no variable. It adds one row per segment, and two rows that keep `x` between the first and the last breakpoint. The rows describe the curve exactly when the points are convex under `>=`, or concave under `<=`. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["a"])) B = Set("B", np.array(["b0", "b1", "b2"])) power = Param.from_dense("power", (G, B), [[0, 10, 20]]) cost = Param.from_dense("cost", (G, B), [[0, 10, 30]]) m = Model("convex") p = m.var("p", (G,)) c = m.var("c", (G,)) m.constraint("demand", Sum(G, p[G]) == 15.0) m.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "tangent") m.set_objective(Sum(G, c[G])) s = m.solve() print(s.objective) print(m.constraints["fuel_tangent"]) ``` Output: ```text 20.0 Constraint('fuel_tangent', ('G', 'fuel_segment'), 2 rows, 4 coefficients) ``` Points that are not convex under `>=` raise `ValueError` and identify the first entity at fault. The tangent rows of such points describe a different curve. ```python raises=ValueError import numpy as np from nimopt import Model, Param, Set G = Set("G", np.array(["a"])) B = Set("B", np.array(["b0", "b1", "b2"])) power = Param.from_dense("power", (G, B), [[0, 10, 20]]) cost = Param.from_dense("cost", (G, B), [[0, 20, 30]]) m = Model("concave") p = m.var("p", (G,)) c = m.var("c", (G,)) m.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "tangent") ``` Raises ValueError: ```text ValueError: piecewise 'fuel' has points that are not convex, required by sign '>=' at {'G': 'a'}; use method 'incremental' ``` ## The method from the data `method="auto"` generates the tangent rows where they describe the curve exactly: the sign is not `==`, no `active` is given, `x` has no constant and the points of every entity are convex under `>=` or concave under `<=`. Otherwise it generates the incremental declarations. `formulation` reports the method the model generated. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["a"])) B = Set("B", np.array(["b0", "b1", "b2"])) power = Param.from_dense("power", (G, B), [[0, 10, 20]]) for values in ([0, 10, 30], [0, 20, 30]): cost = Param.from_dense("cost", (G, B), [values]) m = Model("auto") p = m.var("p", (G,)) c = m.var("c", (G,)) m.constraint("demand", Sum(G, p[G]) == 15.0) fuel = m.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "auto") m.set_objective(Sum(G, c[G])) print(fuel.formulation, m.solve().objective) ``` Output: ```text tangent 20.0 incremental 25.0 ``` ## One curve for every entity `x_points` and `y_points` are over the breakpoint set and over as many of the sets of `x` as the curves differ along. Points over the breakpoint set alone give every entity the same curve. ## Entities with fewer breakpoints A parameter built from a table has entries only where the table lists them. An entity with fewer breakpoints lists its first ones, and the rest are absent. An entity with no breakpoint has no generated rows and no generated columns. Its `x` and `y` are not related by the declaration. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["a", "b", "c"])) B = Set("B", np.array(["b0", "b1", "b2"])) labels = { "G": np.array(["a", "a", "a", "b", "b"]), "B": np.array(["b0", "b1", "b2", "b0", "b1"]), } power = Param.from_long("power", (G, B), labels, [0.0, 10.0, 20.0, 0.0, 15.0]) cost = Param.from_long("cost", (G, B), labels, [0.0, 10.0, 30.0, 0.0, 30.0]) m = Model("short") p = m.var("p", (G,)) c = m.var("c", (G,)) m.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "incremental") print(m.variables["fuel_fill"].n_columns) print(m.constraints["fuel_x"].n_rows) ``` Output: ```text 3 2 ``` A breakpoint absent before a present one raises `ValueError`. A table with an entity of one breakpoint raises `ValueError`. ## A curve that a binary variable switches off `active=` takes a binary variable over the sets of `x`, or a sum of them. Where it is 1, `x` is on the curve. Where it is 0, `x` is 0 and `y` is compared with 0. The curve then starts at its first breakpoint, and a first breakpoint above 0 is a minimum output. `active=` is supported by `method="incremental"` only. A variable with bounds outside 0 and 1 and a scaled variable raise `ValueError`. A continuous variable raises `ValueError` under the default. A value between 0 and 1 scales every breakpoint, so the curve is met at a fraction of its first breakpoint and at a fraction of its cost. `relaxed=True` accepts that variable and declares the scaled curve, which is the linear relaxation of the switch. `relaxed=True` with no `active=` raises `ValueError`. ```python import numpy as np from nimopt import Model, Param, Set, Sum G = Set("G", np.array(["a"])) T = Set("T", np.array(["t0", "t1"])) B = Set("B", np.array(["b0", "b1", "b2"])) power = Param.from_dense("power", (G, B), [[10, 20, 30]]) cost = Param.from_dense("cost", (G, B), [[10, 15, 30]]) demand = Param.from_dense("demand", (T,), [5.0, 25.0]) m = Model("committed") p = m.var("p", (G, T)) c = m.var("c", (G, T)) on = m.var("on", (G, T), upper=1.0, integer=True) spot = m.var("spot", (T,)) m.constraint("balance", Sum(G, p[G, T]) + spot[T] == demand[T]) m.piecewise( "fuel", x=p[G, T], x_points=power[G, B], y=c[G, T], y_points=cost[G, B], sign=">=", method="incremental", active=on[G, T], ) m.set_objective(Sum(G, T, c[G, T]) + 1.2 * Sum(T, spot[T]) + 2.0 * Sum(G, T, on[G, T])) s = m.solve() print(s.objective) print(s.primal("on").to_dense()) print(s.primal("p").to_dense()) ``` Output: ```text 29.0 [[0. 1.]] [[ 0. 20.]] ``` In `t0` the demand of 5 is below the minimum output of 10, and the unit is off. ## A declaration over some of the entities `where` restricts a declaration to some entities of its breakpoints. It is a parameter, a tuple of sets or a domain over the sets of `x_points` other than the breakpoint set, the form of `where` on `Model.constraint`. The breakpoint checks, the generated columns and the generated rows cover the entities at its coordinates. `x`, `y` and `active` are compared at those coordinates only. Two declarations can share the breakpoints and split the entities. In the example, `a` has a status variable and points that are not convex, `b` has convex points, and `c` has no curve. `on` declares the curve of `a` with the incremental method and `active`. `free` declares the curve of `b` with the tangent method. A linear term prices the output of `c`. The variable `fuel` is declared over the entities with a curve, and `p` over every entity. ```python import numpy as np from nimopt import Model, Param, Set, Sum, subset G = Set("G", np.array(["a", "b", "c"])) T = Set("T", np.array(["t0", "t1"])) B = Set("B", np.array(["b0", "b1", "b2"])) points = { "G": np.array(["a", "a", "a", "b", "b", "b"]), "B": np.array(["b0", "b1", "b2", "b0", "b1", "b2"]), } power = Param.from_long("power", (G, B), points, [10.0, 20.0, 30.0, 0.0, 10.0, 20.0]) cost = Param.from_long("cost", (G, B), points, [10.0, 25.0, 30.0, 0.0, 10.0, 30.0]) committed = Param.from_long("committed", (G,), {"G": np.array(["a"])}, [1.0]) flexible = Param.from_long("flexible", (G,), {"G": np.array(["b"])}, [1.0]) price = Param.from_long("price", (G,), {"G": np.array(["c"])}, [3.0]) demand = Param.from_dense("demand", (T,), [5.0, 25.0]) curves = subset( (G, T), {"G": np.array(["a", "a", "b", "b"]), "T": np.array(["t0", "t1", "t0", "t1"])}, ) status = subset((G, T), {"G": np.array(["a", "a"]), "T": np.array(["t0", "t1"])}) m = Model("split") p = m.var("p", (G, T)) fuel = m.var("fuel", (G, T), subset=curves) u = m.var("u", (G, T), subset=status, upper=1.0, integer=True) m.constraint("balance", Sum(G, p[G, T]) == demand[T]) m.piecewise( "on", p[G, T], power[G, B], fuel[G, T], cost[G, B], ">=", "incremental", active=u[G, T], where=committed, ) m.piecewise( "free", p[G, T], power[G, B], fuel[G, T], cost[G, B], ">=", "tangent", where=flexible, ) m.set_objective( Sum(G, T, fuel[G, T]) + Sum(G, T, price[G] * p[G, T]) + 2.0 * Sum(G, T, u[G, T]) ) s = m.solve() print(s.objective) print(m.variables["on_fill"].n_columns, m.constraints["free_tangent"].n_rows) ``` Output: ```text 34.5 4 4 ``` `on_fill` has one column per segment of `a` and timestep. `free_tangent` has one row per segment of `b` and timestep. ## Saving a piecewise declaration A model file of version 4 contains the declaration under `piecewise`. The generated variables, constraints and parameters are not written. Loading the file generates them again. ```python from nimopt import Definition, Sum, dumps d = Definition("curve", sense="min") G, B = d.set("G"), d.set("B") power, cost = d.param("power", (G, B)), d.param("cost", (G, B)) p = d.var("p", (G,)) c = d.var("c", (G,)) d.constraint("demand", Sum(G, p[G]) == 25.0) d.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "incremental") d.set_objective(Sum(G, c[G])) print(dumps(d)) ``` Output: ```text version: 4 name: curve sense: min sets: [G, B] parameters: power: [G, B] cost: [G, B] variables: p: sets: [G] c: sets: [G] constraints: demand: relation: Sum(G, p[G]) == 25 piecewise: fuel: x: p[G] x_points: power[G, B] y: c[G] y_points: cost[G, B] sign: '>=' method: incremental objective: Sum(G, c[G]) ``` `version=3` writes a model in the format a reader of version 3 accepts. The generated declarations are written in place of the declaration, and the breakpoints are not written. A model loaded from that file contains the same rows and no piecewise declaration. A definition has no data to generate the rows from, and writing it as version 3 raises `ValueError`. ```python import numpy as np from nimopt import Model, Param, Set, Sum, dumps G = Set("G", np.array(["a"])) B = Set("B", np.array(["b0", "b1"])) power = Param.from_dense("power", (G, B), [[0, 10]]) cost = Param.from_dense("cost", (G, B), [[0, 20]]) m = Model("line") p = m.var("p", (G,)) c = m.var("c", (G,)) m.piecewise("fuel", p[G], power[G, B], c[G], cost[G, B], ">=", "tangent") print(dumps(m, version=3)) ``` Output: ```text version: 3 name: line sense: min sets: [G, fuel_segment] parameters: fuel_slope: [G, fuel_segment] fuel_intercept: [G, fuel_segment] fuel_x_low: [G] fuel_x_high: [G] variables: p: sets: [G] c: sets: [G] constraints: fuel_tangent: relation: c[G] - fuel_slope[G, fuel_segment] * p[G] >= fuel_intercept[G, fuel_segment] fuel_x_min: relation: p[G] >= fuel_x_low[G] fuel_x_max: relation: p[G] <= fuel_x_high[G] ``` --- # /guides/saving-and-loading # Saving and loading a model A definition writes a YAML file whose expressions are written as they are typed in Python. The file is the canonical form of the model: every derived coefficient written out, every term with its own sum and sign, and the constant last. ```python from nimopt import Definition, Sum, dumps d = Definition("dispatch", sense="min") G, T = d.set("G"), d.set("T") price, eta = d.param("price", (G, T)), d.param("eta", (G, T)) cap, load = d.param("cap", (G, T)), d.param("load", (T,)) gen = d.var("gen", (G, T), upper=cap) d.constraint("balance", Sum(G, gen[G, T]) == load[T]) d.set_objective(Sum(G, T, 2 * (price[G, T] / eta[G, T]) * gen[G, T])) print(dumps(d)) ``` Output: ```text version: 4 name: dispatch sense: min sets: [G, T] parameters: price: [G, T] eta: [G, T] cap: [G, T] load: [T] variables: gen: sets: [G, T] upper: cap constraints: balance: relation: Sum(G, gen[G, T]) == load[T] objective: Sum(G, T, ((price[G, T] / eta[G, T]) * 2) * gen[G, T]) ``` The structure section is the schema of the data: every set, and every parameter with its dimensions. `build` raises `ValueError` for a mapping that omits one of them. ## Reading a file back `loads` reads text and `load` reads a path. A file without data returns a `Definition`, whose next step is `build(data)`. ```python import numpy as np from nimopt import loads text = """ version: 3 name: dispatch sense: min sets: [G, T] parameters: cost: [G] load: [T] variables: gen: sets: [G, T] upper: 10.0 constraints: balance: relation: Sum(G, gen[G, T]) == load[T] objective: Sum(G, T, cost[G] * gen[G, T]) """ d = loads(text) m = d.build( { "G": np.array(["a", "b"]), "T": np.arange(2), "cost": np.array([1.0, 3.0]), "load": np.array([12.0, 15.0]), } ) print(m.solve().objective) ``` Output: ```text 41.0 ``` ## Editing by hand The text is read through the same operators a Python model is built from. An edit is accepted where Python accepts it, and is normalized the same way. Adding a scalar, reordering terms and reversing a comparison all parse. The file written back is the canonical form. ```python from nimopt import dumps, loads edited = """ version: 3 name: dispatch sense: min sets: [G, T] parameters: cost: [G] load: [T] variables: gen: sets: [G, T] constraints: balance: relation: load[T] == Sum(G, gen[G, T]) * 2 + 1 - 1 objective: Sum(G, T, cost[G] * gen[G, T]) """ print(dumps(loads(edited)).splitlines()[-2]) ``` Output: ```text relation: 2 * Sum(G, gen[G, T]) == load[T] ``` An edit that raises in Python raises here with the same message. A construct outside the expression syntax raises and reports it. ```python raises=ValueError from nimopt import loads loads( """ version: 3 name: dispatch sense: min sets: [G, T] variables: gen: sets: [G, T] constraints: peak: relation: max(gen[G, T]) <= 10 """ ) ``` Raises ValueError: ```text ValueError: 'max(gen[G, T]) <= 10': the syntax supports one call; write Sum ``` ## Data inline, for a model small enough to read A built model writes its data into the file with `inline=True`. A set is a list, and a dense parameter is nested lists. A parameter with values at some coordinates of its product is a table of the dimensions then `value`. A file containing data loads to a built `Model`. ```python from nimopt import dumps, loads from nimopt.models import transport m = transport.definition().build(transport.data()) text = dumps(m, inline=True) print(text[text.index("data:") :]) print(loads(text).solve().objective) ``` Output: ```text data: P: [p0, p1, p2, p3] W: [w0, w1, w2, w3, w4, w5] cost: columns: [P, W, value] rows: - [p0, w0, 1.1322210842282328] - [p0, w1, 7.506161913602179] - [p0, w4, 1.3277881914895575] - [p1, w0, 6.835972487871987] - [p1, w3, 8.302044618221775] - [p1, w4, 5.853086206137439] - [p2, w2, 5.348999931723383] - [p2, w3, 8.480579390302147] - [p2, w4, 7.526828432972257] - [p3, w1, 1.0219080013611848] - [p3, w2, 7.859234212700555] - [p3, w3, 1.2686846024437148] supply: [60.0, 60.0, 60.0, 60.0] demand: [10.0, 10.0, 10.0, 10.0, 10.0, 10.0] 100.99601811246072 ``` ## Datetime members A set whose members are `datetime64` or `timedelta64` round trips through both data sources, in every unit. Inline, such a set is written as a mapping of `dtype` and `members`: a `datetime64` member as its ISO 8601 string, and a `timedelta64` member as its integer count of the unit in the `dtype`. A member fixed in a relation is written as quoted text. ```python import numpy as np from nimopt import Model, Param, Set, Sum, dumps, loads T = Set("T", np.array(["2030-01-01T00", "2030-01-01T01"], dtype="datetime64[h]")) cost = Param.from_dense("cost", (T,), np.array([2.0, 5.0])) m = Model("dispatch", sense="min") gen = m.var("gen", (T,), upper=10.0) m.constraint("start", gen["2030-01-01T00"] == 4.0) m.constraint("total", Sum(T, gen[T]) >= 6.0) m.set_objective(Sum(T, cost[T] * gen[T])) text = dumps(m, inline=True) print(text[text.index("constraints:") :]) print(loads(text).solve().objective) ``` Output: ```text constraints: start: relation: gen['2030-01-01T00'] == 4 total: relation: Sum(T, gen[T]) >= 6 objective: Sum(T, cost[T] * gen[T]) data: T: dtype: datetime64[h] members: [2030-01-01T00, 2030-01-01T01] cost: [2.0, 5.0] 18.0 ``` The fixed member is written as a string and is read back to the same member. A string is parsed as ISO 8601, and a `datetime.datetime`, a `datetime.date` and a `datetime64` of another unit are converted. An integer against a `timedelta64` set is a count of that set's own unit. A conversion that is not exact raises `ValueError` instead of truncating. ```python raises=ValueError import numpy as np from nimopt import Model, Set T = Set("T", np.array(["2030-01-01T00", "2030-01-01T01"], dtype="datetime64[h]")) m = Model("dispatch") gen = m.var("gen", (T,)) gen["2030-01-01T00:30"] ``` Raises ValueError: ```text ValueError: member '2030-01-01T00:30' does not convert exactly to datetime64[h] at dimension 'T' of variable 'gen'; write a member in the unit of that dimension ``` `Model.row` takes the same text for a datetime dimension, and a `Row` displays each datetime coordinate in it. ## Data beside the file, for a large model `save` writes a model's file and its data as an `.npz` beside it, under the file's stem. The file records the sidecar name. `load` reads both. The file records no path and no machine name. ```python import tempfile from pathlib import Path from nimopt import load, save from nimopt.models import storage m = storage.definition().build(storage.data()) with tempfile.TemporaryDirectory() as held: path = Path(held) / "storage.yaml" save(m, path) print(sorted(p.name for p in Path(held).iterdir())) print(path.read_text().splitlines()[-1]) print(load(path).solve().status) ``` Output: ```text ['storage.npz', 'storage.yaml'] data: storage.npz optimal ``` A definition's file takes data from the caller instead: `load(path, data=...)` with the mapping `build` takes or the path of an `.npz`. A file that contains data and a `data=` together raises `ValueError`: one model takes one data source. ## A file that describes its own format `instructions=True` on `save` and `dumps` writes a comment block at the top of the file. The block is the same in every file. It describes the format, not the model: the keys and their order, the defaults, the rules that determine which rows a constraint has, and the expression syntax. A reader with one file interprets it without this package. The block is a YAML comment. A file with it loads to the same model as one without it, and writing the loaded model with the flag gives the same text. ```python from nimopt import Definition, dumps, loads d = Definition("dispatch", sense="min") T = d.set("T") load, gen = d.param("load", (T,)), d.var("gen", (T,)) d.constraint("balance", gen[T] == load[T]) text = dumps(d, instructions=True) print("\n".join(text.splitlines()[:5])) print(dumps(loads(text)) == dumps(d)) print(dumps(loads(text), instructions=True) == text) ``` Output: ```text # --- Reading this file -------------------------------------------------- # A nimopt model file, format version 4. The keys are written in this # order, and no other key is accepted: version, name, sense, sets, # aliases, parameters, variables, constraints, piecewise, objective, # data. Only version, name and sense are required. True True ``` --- # /guides/subsets # A variable over a subset In a network model, most pairs of a set product are not connected: a plant serves some warehouses, a line joins two of many buses. A variable declared over the full product has a column for every pair, including the ones that do not exist. `subset=` restricts the variable to the members that do, and the others have no column at all. ## Declaring the subset `subset(sets, columns)` lists members of a set product by label, one column per set, read in parallel. `m.var(..., subset=arcs)` declares the variable over those members. ```python import numpy as np from nimopt import Model, Set, Sum, product, subset P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) arcs = subset( (P, W), {"P": np.array(["p0", "p0", "p1"]), "W": np.array(["w0", "w1", "w2"])}, ) m = Model("network") x = m.var("x", (P, W), subset=arcs) print(product((P, W)).size) print(x.n_columns) print(m.assemble().to_dense()) ``` Output: ```text 6 3 [] ``` The product has six members and the subset three, and `x` has three columns. The assembled matrix is empty: no constraint has been added. A model over a sparse network costs its arcs, not the grid that contains them. ## By label or by position `subset` takes labels and `subset_of` takes integer positions. Both read their columns in parallel: the k-th entry of each column belongs to the same member. A subset is a list of members, not a cross product of its columns. ```python import numpy as np from nimopt import Set, subset, subset_of P = Set("P", np.array(["p0", "p1"])) W = Set("W", np.array(["w0", "w1", "w2"])) by_label = subset( (P, W), {"P": np.array(["p0", "p1"]), "W": np.array(["w0", "w2"])}, ) by_index = subset_of((P, W), np.array([[0, 1], [0, 2]])) print(by_label.size, by_index.size) print(by_label.labels()) ``` Output: ```text 2 2 {'P': array(['p0', 'p1'], dtype='= 0 constraint balance (T) Sum(G, gen[T, G]) == load[T] objective min Sum(T, G, cost[G] * gen[T, G]) + Sum(T, G, no_load[G] * on[T, G]) ``` Snapshots are uncoupled. `reference` enumerates every on-off subset per snapshot and takes the cheapest feasible one. The optimum is computed without a solver. ```python from nimopt.models import commitment inputs = commitment.data() model = commitment.definition().build(inputs) solution = model.solve() print(int(model.integrality().sum()), "binary columns of", model.n_columns) print(solution.objective, commitment.reference(inputs)) ``` Output: ```text 12 binary columns of 24 13800.0 13800.0 ``` Both unit rows are produced for every generator and snapshot. ```python from nimopt.models import commitment model = commitment.definition().build(commitment.data()) print(model.absent("capacity")) ``` Output: ```text capacity 12 of 12 rows stated by terms ``` --- # /models/dispatch # Dispatch `nimopt.models.dispatch` is least-cost dispatch of a generator fleet against a load. One variable `p` is indexed over snapshots and generators, there is one balance row per snapshot, and each generator has a cost. Every other model in the corpus adds one axis to this one. ```text minimize Σ_{t,g} cost[g] · p[t,g] subject to Σ_g p[t,g] == load[t] for each snapshot t 0 ≤ p[t,g] ≤ p_max[g] ``` The balance row has no coefficient. A sum over a dimension requires none, and the corpus writes no coefficient a model does not require. ```python from nimopt.models import dispatch print(dispatch.definition().explain()) ``` Output: ```text dispatch min not built sets snapshot · generator parameters p_max (generator) · load (snapshot) · cost (generator) variables p (snapshot×generator) [0.0, p_max] constraint balance (snapshot) Sum(generator, p[snapshot, generator]) == load[snapshot] objective min Sum(snapshot, generator, cost[generator] * p[snapshot, generator]) ``` Snapshots are independent, and the optimum is the merit order per snapshot. `reference` computes it without a solver. ```python from nimopt.models import dispatch inputs = dispatch.data() solution = dispatch.definition().build(inputs).solve() print(solution.status) print(solution.objective, dispatch.reference(inputs)) ``` Output: ```text optimal 1920.0 1920.0 ``` The balance row is produced for every snapshot in the load, and `absent` reports no dropped row. ```python from nimopt.models import dispatch model = dispatch.definition().build(dispatch.data()) print(model.absent("balance")) ``` Output: ```text balance 6 of 6 rows stated by terms ``` --- # /models/expansion # Expansion `nimopt.models.expansion` is a two-stage stochastic program. The first stage builds capacity. The second stage dispatches it against a demand and a fuel price given by the scenario. The shape of the first-stage variable expresses the information available to it: `cap` is indexed by technology alone, and one capacity applies to every scenario. `p` and `shed` are indexed by scenario as well and may differ across it. ```text minimize Σ_g capital[g] · cap[g] + Σ_{s,g,t} weight[s] · cost[s,g] · p[s,g,t] + Σ_{s,t} weight[s] · voll[s] · shed[s,t] subject to p[s,g,t] ≤ cap[g] for each s, g, t Σ_g p[s,g,t] + shed[s,t] == demand[s,t] for each s, t cap, p, shed ≥ 0 ``` `weight` is the probability of a scenario. The second and third sums are therefore an expectation, and the objective is the committed capital plus the expected cost of the recourse. No row relates the capacity of one scenario to the capacity of another: the missing dimension expresses non-anticipativity. ```python from nimopt.models import expansion print(expansion.definition().explain()) ``` Output: ```text expansion min not built sets S · G · T parameters capital (G) · cost (S,G) · demand (S,T) · weight (S) · voll (S) variables cap (G) [0.0, inf] · p (S×G×T) [0.0, inf] · shed (S×T) [0.0, inf] constraint capacity (S,G,T) p[S, G, T] - cap[G] <= 0 constraint balance (S,T) Sum(G, p[S, G, T]) + shed[S, T] == demand[S, T] objective min Sum(G, capital[G] * cap[G]) + Sum(S, G, T, (weight[S] * cost[S, G]) * p[S, G, T]) + Sum(S, T, (weight[S] * voll[S]) * shed[S, T]) ``` `cap` has one column per technology, and `p` has one per scenario, technology and hour. The capacity row is declared over all three dimensions, and the variable it bounds has one. The rows of every scenario read the same column, and the capacity is therefore a shared decision. ```python from nimopt.models import expansion model = expansion.definition().build(expansion.data()) print(model.explain()) ``` Output: ```text expansion min 75 columns · 72 rows · 180 nonzeros sets G 3 · S 3 · T 6 parameters demand (S,T) 18 · capital (G) 3 · weight (S) 3 · cost (S,G) 9 · voll (S) 3 variables cap (G) 3 cols [0.0, inf] · p (S×G×T) 54 cols [0.0, inf] · shed (S×T) 18 cols [0.0, inf] constraint capacity (S,G,T) p[S, G, T] - cap[G] <= 0 54 rows 108 nz constraint balance (S,T) Sum(G, p[S, G, T]) + shed[S, T] == demand[S, T] 18 rows 72 nz objective min Sum(G, capital[G] * cap[G]) + Sum(S, G, T, (weight[S] * cost[S, G]) * p[S, G, T]) + Sum(S, T, (weight[S] * voll[S]) * shed[S, T]) ``` A technology is available in full wherever it is built, and no row couples one hour to the next. The recourse in each scenario-hour is therefore the merit order of the built capacity against that demand, with the remainder unserved at `voll`. The capacity is written as bands: the band of the cheapest technology, then the next, and last the band between the most expensive technology and lost load. The total separates into one term per band. Each term is convex in the level of its band and changes slope only at a demand. `reference` minimizes the terms one at a time and adds them, and the result is arithmetic over the inputs, not a second solve. ```python from nimopt.models import expansion inputs = expansion.data() solution = expansion.definition().build(inputs).solve() print(solution.objective, expansion.reference(inputs)) ``` Output: ```text 29703.899999999994 29703.899999999994 ``` The data describes a cost frontier: capital falls as marginal cost rises, and lost load is more expensive than the most expensive technology. All three technologies are built. The capacity stacks to 180 and leaves the peak hour of 225 in the cold scenario short. Shedding 45 costs less than a fourth band used in one hour of one scenario. ```python import numpy as np from nimopt.models import expansion inputs = expansion.data() solution = expansion.definition().build(inputs).solve() print(np.cumsum(np.asarray(solution.primal("cap").values()))) print(np.asarray(solution.primal("shed").values()).reshape(3, 6)) ``` Output: ```text [119. 153. 180.] [[ 0. 0. 0. 0. 0. 0.] [ 0. 0. 0. 0. 0. 0.] [ 0. 0. 0. 45. 0. 0.]] ``` The separation applies only to data whose technologies form a frontier and whose bands stack. `reference` raises on any other data, and it returns no number that is not the optimum. ```python raises=ValueError from nimopt.models import expansion inputs = expansion.data() inputs["capital"] = inputs["capital"][::-1] print(expansion.reference(inputs)) ``` Raises ValueError: ```text ValueError: capital does not fall from base to what follows it; pass a capital cost that falls across the merit order ``` --- # /models/fleet # Fleet `nimopt.models.fleet` is the same problem as `dispatch`. It is declared as one variable per unit over the snapshots alone, and the terms of the units are added into one balance row. The optimum is the same merit order, and the cost of the declaration differs. ```text minimize Σ_t Σ_u cost_u[t] · u[t] subject to Σ_u u[t] == load[t] for each snapshot t 0 ≤ u[t] ≤ p_max_u[t] for each unit u ``` `definition` takes a scale here: the number of variables is a property of the declaration, not of the data. ```python from nimopt.models import fleet print(fleet.definition().explain()) ``` Output: ```text fleet min not built sets T parameters load (T) · p_max_g0_0 (T) · cost_g0_0 (T) · p_max_g1_0 (T) · cost_g1_0 (T) · p_max_g2_0 (T) · cost_g2_0 (T) variables g0_0 (T) [0.0, p_max_g0_0] · g1_0 (T) [0.0, p_max_g1_0] · g2_0 (T) [0.0, p_max_g2_0] constraint balance (T) g0_0[T] + g1_0[T] + g2_0[T] == load[T] objective min Sum(T, cost_g0_0[T] * g0_0[T]) + Sum(T, cost_g1_0[T] * g1_0[T]) + Sum(T, cost_g2_0[T] * g2_0[T]) ``` ```python from nimopt.models import fleet inputs = fleet.data() model = fleet.definition().build(inputs) solution = model.solve() print(len(model.variables), "variables,", model.n_columns, "columns") print(solution.objective, fleet.reference(inputs)) ``` Output: ```text 3 variables, 12 columns 8900.0 8900.0 ``` One balance row per hour, and every unit appears in every row. ```python from nimopt.models import fleet model = fleet.definition().build(fleet.data()) print(model.absent("balance")) ``` Output: ```text balance 4 of 4 rows stated by terms ``` --- # /models # Worked models `nimopt.models` contains ten models. Each is a module with three functions. | Name | Returns | | --- | --- | | `definition()` | a `Definition`: the formulation with no data bound | | `data(scale=1)` | the inputs, at a given size | | `reference(data)` | the optimal objective, computed by direct arithmetic | The three functions have different callers. This documentation calls `explain()`, a benchmark calls `build(data(100))`, and a test compares a solve against `reference(data(1))`. A reference is arithmetic over the inputs and calls nothing from `nimopt`. A formulation error is therefore not checked against a copy of itself. | Model | Exercises | | --- | --- | | [`dispatch`](/models/dispatch) | the baseline formulation | | [`transport`](/models/transport) | a sparse network over a subset of a product | | [`storage`](/models/storage) | temporal coupling and a cyclic lag | | [`nodal`](/models/nodal) | grouping through a lookup parameter | | [`commitment`](/models/commitment) | binary columns | | [`fleet`](/models/fleet) | many small declarations rather than one large one | | [`profiled`](/models/profiled) | a bound that varies by hour | | [`sector`](/models/sector) | mixed density: dense in one axis, sparse in another | | [`expansion`](/models/expansion) | a two-stage stochastic program: capacity before the scenario, dispatch after | | [`recourse`](/models/recourse) | a binary first stage taken before the scenario is known | ```python from nimopt.models import dispatch inputs = dispatch.data() solution = dispatch.definition().build(inputs).solve() print(solution.status, solution.objective, dispatch.reference(inputs)) ``` Output: ```text optimal 1920.0 1920.0 ``` --- # /models/nodal # Nodal `nimopt.models.nodal` groups generators into buses through a lookup parameter. `at[G, B]` has an entry where generator `g` is located at bus `b`. Multiplying the generation by it maps a row over generators to a row over buses. The coefficient introduces `B`, a dimension no variable has, and the balance is indexed over the dimensions the lookup defines. ```text minimize Σ_{t,g} cost[g] · gen[t,g] subject to Σ_g at[g,b] · gen[t,g] == demand[b,t] for each bus b and hour t 0 ≤ gen[t,g] ≤ p_max[g] ``` ```python from nimopt.models import nodal print(nodal.definition().explain()) ``` Output: ```text nodal min not built sets T · G · B parameters at (G,B) · p_max (G) · cost (G) · demand (B,T) variables gen (T×G) [0.0, p_max] constraint balance (B,T) Sum(G, at[G, B] * gen[T, G]) == demand[B, T] objective min Sum(T, G, cost[G] * gen[T, G]) ``` Each bus meets its own demand from the generators sited at it. The optimum is a merit order per bus and hour. ```python from nimopt.models import nodal inputs = nodal.data() solution = nodal.definition().build(inputs).solve() print(solution.objective, nodal.reference(inputs)) ``` Output: ```text 11850.0 11850.0 ``` Every bus-hour has a row. A generator sited elsewhere is a term the lookup removes from that row, not a row that is dropped. ```python from nimopt.models import nodal model = nodal.definition().build(nodal.data()) print(model.absent("balance")) ``` Output: ```text balance 6 of 6 rows stated by terms term absent G='g0_0', B='b1_0', T=0 gen absent-coefficient (at) term absent G='g0_0', B='b1_0', T=1 gen absent-coefficient (at) term absent G='g0_0', B='b1_0', T=2 gen absent-coefficient (at) term absent G='g1_0', B='b1_0', T=0 gen absent-coefficient (at) term absent G='g1_0', B='b1_0', T=1 gen absent-coefficient (at) term absent G='g1_0', B='b1_0', T=2 gen absent-coefficient (at) term absent G='g2_0', B='b0_0', T=0 gen absent-coefficient (at) term absent G='g2_0', B='b0_0', T=1 gen absent-coefficient (at) term absent G='g2_0', B='b0_0', T=2 gen absent-coefficient (at) term absent G='g3_0', B='b0_0', T=0 gen absent-coefficient (at) term absent G='g3_0', B='b0_0', T=1 gen absent-coefficient (at) term absent G='g3_0', B='b0_0', T=2 gen absent-coefficient (at) ``` --- # /models/profiled # Profiled `nimopt.models.profiled` is a dispatch whose capacity varies by hour. `dispatch` bounds a generator by a single number. Here `p_max` is a parameter over generators and snapshots. A solar unit is bounded by its hourly availability, and a thermal unit by its rating. ```text minimize Σ_{t,g} cost[g] · gen[t,g] subject to Σ_g gen[t,g] == load[t] for each snapshot t 0 ≤ gen[t,g] ≤ profile[g,t] ``` The profile is indexed `(G, T)` and the variable `(T, G)`. A bound is read in the dimension order of the variable it bounds, and both orderings select the same columns. ```python from nimopt.models import profiled print(profiled.definition().explain()) ``` Output: ```text profiled min not built sets T · G parameters profile (G,T) · cost (G) · load (T) variables gen (T×G) [0.0, profile] constraint balance (T) Sum(G, gen[T, G]) == load[T] objective min Sum(T, G, cost[G] * gen[T, G]) ``` Each snapshot is independent. The optimum is the merit order against the capacities of that hour. ```python from nimopt.models import profiled inputs = profiled.data() solution = profiled.definition().build(inputs).solve() print(solution.objective, profiled.reference(inputs)) ``` Output: ```text 28119.536003699297 28119.536003699293 ``` Every hour has a balance row. A generator whose profile is zero is a column bounded to zero, not a dropped row. ```python from nimopt.models import profiled model = profiled.definition().build(profiled.data()) print(model.absent("balance")) ``` Output: ```text balance 8 of 8 rows stated by terms ``` --- # /models/recourse # Recourse `nimopt.models.recourse` is `commitment` under uncertainty. The demand and the fuel price are a scenario, and the on-off decision is taken before either is given. `on` is indexed by hour and unit alone, and `p` and `shed` are indexed by the scenario as well. One commitment applies to every scenario, and the binary column is therefore taken under uncertainty. ```text minimize Σ_{t,g} no_load[g] · on[t,g] + Σ_{s,g,t} weight[s] · cost[s,g] · p[s,g,t] + Σ_{s,t} weight[s] · voll[s] · shed[s,t] subject to p[s,g,t] ≤ p_max[g] · on[t,g] for each s, g, t p[s,g,t] ≥ p_min[g] · on[t,g] for each s, g, t Σ_g p[s,g,t] + shed[s,t] == demand[s,t] for each s, t on ∈ {0,1}, p, shed ≥ 0 ``` The first stage costs the same in every scenario. The second and third sums are weighted by the probability of a scenario, and the objective is a commitment charge plus the expected cost of the recourse. ```python from nimopt.models import recourse print(recourse.definition().explain()) ``` Output: ```text recourse min not built sets S · G · T parameters p_max (G) · p_min (G) · no_load (G) · cost (S,G) · demand (S,T) · weight (S) · voll (S) variables on (T×G) [0.0, 1.0] integer · p (S×G×T) [0.0, inf] · shed (S×T) [0.0, inf] constraint capacity (S,G,T) p[S, G, T] - p_max[G] * on[T, G] <= 0 constraint minimum (S,G,T) p[S, G, T] - p_min[G] * on[T, G] >= 0 constraint balance (S,T) Sum(G, p[S, G, T]) + shed[S, T] == demand[S, T] objective min Sum(T, G, no_load[G] * on[T, G]) + Sum(S, G, T, (weight[S] * cost[S, G]) * p[S, G, T]) + Sum(S, T, (weight[S] * voll[S]) * shed[S, T]) ``` No row couples one hour to the next, and the commitment is chosen hour by hour. `reference` enumerates every on-off subset of the fleet and scores each by its expected recourse across the scenarios. That enumeration is exact and cheap: three units make eight subsets. ```python from nimopt.models import recourse inputs = recourse.data() solution = recourse.definition().build(inputs).solve() print(solution.objective, recourse.reference(inputs)) ``` Output: ```text 16744.8125 16744.8125 ``` A committed unit runs at least its minimum, and the model has no sink for unwanted energy. A unit whose minimum exceeds the demand of the mildest scenario is therefore not committed: its output would violate the balance row of that scenario. The first hour demands 38.25 in the mildest scenario. The minimum of `base` is 40.0 and its fuel is the cheapest of the three, and it is not committed in that hour. The third hour is the opposite case: the fleet supplies 260.0 against a coldest demand of 268.75, every unit is committed, and the remainder is shed. ```python import numpy as np from nimopt.models import recourse inputs = recourse.data() solution = recourse.definition().build(inputs).solve() print(inputs["demand"].round(2)) print(np.asarray(solution.primal("on").values()).reshape(4, 3)) print(np.asarray(solution.primal("shed").values()).reshape(3, 4)) ``` Output: ```text [[ 38.25 119. 182.75 80.75] [ 45. 140. 215. 95. ] [ 56.25 175. 268.75 118.75]] [[0. 1. 0.] [1. 1. 0.] [1. 1. 1.] [1. 0. 0.]] [[0. 0. 0. 0. ] [0. 0. 0. 0. ] [0. 0. 8.75 0. ]] ``` The mean demand in the first hour is 43.875. `base` serves that demand at the lowest cost. A model given that one number commits `base`, and that commitment violates the balance row of the mildest scenario, whose probability is 0.5. The scenario dimension on `demand` excludes the commitment, and the missing scenario dimension on `on` applies that exclusion to every scenario. --- # /models/sector # Sector `nimopt.models.sector` has mixed density. The region-technology map is sparse: a technology exists in some regions and not in others. Every sited pair runs in every hour. The generation variable takes its members from a parameter over the sited pairs crossed with the whole horizon. It is therefore sparse in one axis and dense in the other. ```text minimize Σ_{(r,k) sited, t} cost[r,k] · gen[r,k,t] subject to Σ_k gen[r,k,t] == demand[r,t] for each region r and hour t 0 ≤ gen[r,k,t] ≤ capacity[r,k] for each sited (r,k) and hour t ``` ```python from nimopt.models import sector print(sector.definition().explain()) ``` Output: ```text sector min not built sets R · K · T parameters sited (R,K,T) · capacity (R,K) · cost (R,K) · demand (R,T) variables gen (R×K×T) over sited [0.0, capacity] constraint balance (R,T) Sum(K, gen[R, K, T]) == demand[R, T] objective min Sum(R, K, T, cost[R, K] * gen[R, K, T]) ``` Each region meets its own demand from the technologies sited in it. The optimum is a merit order per region and hour. ```python from nimopt.models import sector inputs = sector.data() model = sector.definition().build(inputs) solution = model.solve() print( model.n_columns, "columns of a possible", len(inputs["R"]) * len(inputs["K"]) * len(inputs["T"]), ) print(solution.objective, sector.reference(inputs)) ``` Output: ```text 16 columns of a possible 24 18350.0 18350.0 ``` Every region-hour has a balance row. The capacity bound applies to the sited pairs alone, and no balance row is dropped. ```python from nimopt.models import sector model = sector.definition().build(sector.data()) print(model.absent("balance")) ``` Output: ```text balance 8 of 8 rows stated by terms ``` --- # /models/storage # Storage `nimopt.models.storage` dispatches a generator fleet and a set of batteries against an hourly load. The `state_of_charge` row references the previous hour through `T.cyclic - 1`. The row at the first hour therefore references the last hour, and every hour has a row. The ramp row references `T - 1`. The first hour has no predecessor, and that row is not produced. The model exercises both lag rules. ```text minimize Σ_{g,t} cost[g,t] · gen[g,t] subject to Σ_g gen[g,t] + Σ_s discharge[s,t] − Σ_s charge[s,t] == load[t] soc[s,t] − soc[s,t−1] − charge_eta[s,t] · charge[s,t] + discharge_eta[s,t] · discharge[s,t] == 0 (t−1 wraps) gen[g,t] − gen[g,t−1] ≤ ramp_limit[g,t] (t=0 dropped) gen[g,t] ≤ capacity[g,t] charge[s,t] ≤ power[s,t] discharge[s,t] ≤ power[s,t] soc[s,t] ≤ energy[s,t] ``` Every limit is a constraint and not a bound, and each therefore has a dual value. ```python from nimopt.models import storage print(storage.definition().explain()) ``` Output: ```text storage min not built sets T · G · S parameters cost (G,T) · capacity (G,T) · ramp_limit (G,T) · load (T) · power (S,T) · energy (S,T) · charge_eta (S,T) · discharge_eta (S,T) variables gen (G×T) [0.0, inf] · charge (S×T) [0.0, inf] · discharge (S×T) [0.0, inf] · soc (S×T) [0.0, inf] constraint balance (T) Sum(G, gen[G, T]) + Sum(S, discharge[S, T]) - Sum(S, charge[S, T]) == load[T] constraint state_of_charge (S,T) soc[S, T] - soc[S, T.cyclic - 1] - charge_eta[S, T] * charge[S, T] + discharge_eta[S, T] * discharge[S, T] == 0 constraint generation_limit (G,T) gen[G, T] <= capacity[G, T] constraint ramp (G,T) gen[G, T] - gen[G, T - 1] <= ramp_limit[G, T] constraint charge_limit (S,T) charge[S, T] <= power[S, T] constraint discharge_limit (S,T) discharge[S, T] <= power[S, T] constraint energy_limit (S,T) soc[S, T] <= energy[S, T] objective min Sum(G, T, cost[G, T] * gen[G, T]) ``` The batteries are lossy, with a round-trip efficiency of `0.95 · 0.93`, and the costs of the fleet span 50.0 to 55.0. Shifting energy through the store costs more than it saves. The store is therefore idle, and the optimum is the hourly merit order. A store that cycles requires data with a wider cost spread, and the benchmarks supply it. ```python from nimopt.models import storage inputs = storage.data() solution = storage.definition().build(inputs).solve() print(solution.objective, storage.reference(inputs)) print("store moved:", abs(solution.primal("charge").to_dense()).max()) ``` Output: ```text 145449.3739227024 145449.37392270242 store moved: 0.0 ``` Every hour has a `state_of_charge` row: that lag wraps. The ramp row at the first hour is not produced: that lag does not wrap. ```python from nimopt.models import storage model = storage.definition().build(storage.data()) print(model.absent("state_of_charge")) print() print(model.absent("ramp")) ``` Output: ```text state_of_charge 24 of 24 rows stated by terms ramp 69 of 72 rows stated by terms row absent G='base0', T=0 term-does-not-reach (gen) row absent G='mid0', T=0 term-does-not-reach (gen) row absent G='peak0', T=0 term-does-not-reach (gen) ``` --- # /models/transport # Transport `nimopt.models.transport` ships from plants to warehouses over an incomplete network: a plant serves a band of nearby warehouses and not all of them. The cost parameter has one entry per arc, and the flow variable takes its members from that parameter. The model therefore has one column per arc, not one per cell of the plant-warehouse product. ```text minimize Σ_{(p,w) ∈ arcs} cost[p,w] · flow[p,w] subject to Σ_w flow[p,w] ≤ supply[p] for each plant p Σ_p flow[p,w] ≥ demand[w] for each warehouse w flow[p,w] ≥ 0 for each arc (p,w) ``` ```python from nimopt.models import transport print(transport.definition().explain()) ``` Output: ```text transport min not built sets P · W parameters cost (P,W) · supply (P) · demand (W) variables flow (P×W) over cost [0.0, inf] constraint supply (P) Sum(W, flow[P, W]) <= supply[P] constraint demand (W) Sum(P, flow[P, W]) >= demand[W] objective min Sum(P, W, cost[P, W] * flow[P, W]) ``` Supply is twice the total demand of the band of a plant. No supply row binds therefore, and each warehouse buys from the cheapest plant connected to it. `reference` computes that sum. ```python from nimopt.models import transport inputs = transport.data() model = transport.definition().build(inputs) solution = model.solve() print(model.n_columns, "columns for", len(inputs["cost"][1]), "arcs") print(solution.objective, transport.reference(inputs)) ``` Output: ```text 12 columns for 12 arcs 100.99601811246072 100.99601811246073 ``` Arcs are drawn from every warehouse but the last. No plant serves the last warehouse, and it has no demand row. `absent` reports the row and the rule that dropped it. ```python from nimopt.models import transport model = transport.definition().build(transport.data()) print(model.absent("demand")) ``` Output: ```text demand 5 of 6 rows stated by terms row absent W='w5' term-does-not-reach (flow) ```