OpenGoddard

Search symbols

Type part of a name or a dotted path

Pseudospectral optimal control with NumPy and SciPy.

Classes

Name	Description
Condition	Accumulate dimensionless equality or inequality residuals.
Dynamics	Scaled state derivatives at a section's collocation nodes.
DynamicsModel	Pointwise physical dynamics function(x, u, t, obj, section, xp).
Guess	Construct physical initial guesses on a supplied time grid.
LinearCallback	Declare a cost/boundary/path callback affine in the scaled vector p.
Problem	Define a physical optimal-control problem and solve its LGL discretization.

Modules

• derivativesStructured collocation derivatives and optional NumPy automatic differentiation.
• hpResidual-guided h/hp adaptation inside unchanged logical phases.
• meshPhase-preserving mesh refinement and independent, between-node defect checks.
• optimizeOpenGoddard.optimize - Optimal Trajectories module with PseudoSpectral Method

Conditionclass#

Re-exported from OpenGoddard._helpersView source

class Condition(length=0, *, backend=np)

Accumulate dimensionless equality or inequality residuals.

Equality residuals must be zero; inequality residuals must be nonnegative. Scalar and one-dimensional inputs follow NumPy broadcasting rules. A positive reference unit scales each residual. Prefer problem.condition() inside callbacks so automatic differentiation selects the correct backend.

Parameters

Name	Type	Default	Description
length	int	0	Number of initial zero residuals, normally zero.
backend	module	np	NumPy-compatible module, supplied by the callback view.

Examples

>>> result = Condition()
>>> result.equal(5, 5)
>>> result()
array([0.])

For physical bounds 0 <= x <= 100 in an inequality callback, call result.lower_bound(x, 0) and result.upper_bound(x, 100).

addmethod#

View source

def add(arg, unit=1.0)

Append scalar or one-dimensional residuals using a positive scale.

change_valuemethod#

View source

def change_value(index, value)

Replace an existing dimensionless residual by its flattened index.

equalmethod#

View source

def equal(arg1, arg2, unit=1.0)

Append (arg1 - arg2)/unit; enforces arg1 = arg2 in an equality callback.

Parameters

Name	Type	Default	Description
arg1	float or array_like	required	Left side of the equation.
arg2	float or array_like	required	Right side of the equation.
unit	float	1.0	Positive reference unit for the residual.

lower_boundmethod#

View source

def lower_bound(arg1, arg2, unit=1.0)

Append (arg1 - arg2)/unit; enforces arg1 >= arg2 in inequalities.

Parameters

Name	Type	Default	Description
arg1	float or array_like	required	Physical values to constrain.
arg2	float or array_like	required	Physical lower bounds.
unit	float	1.0	Positive reference unit for the residual.

upper_boundmethod#

View source

def upper_bound(arg1, arg2, unit=1.0)

Append (arg2 - arg1)/unit; enforces arg1 <= arg2 in inequalities.

Parameters

Name	Type	Default	Description
arg1	float or array_like	required	Physical values to constrain.
arg2	float or array_like	required	Physical upper bounds.
unit	float	1.0	Positive reference unit for the residual.

Dynamicsclass#

Re-exported from OpenGoddard._helpersView source

class Dynamics(prob, section=0)

Scaled state derivatives at a section’s collocation nodes.

Assign physical dx/dt values with dx[state] = values. Calling dx() converts them to d(x/unit_states)/d(t/unit_time), flattened in state order for the LGL residual. Scalars are broadcast to every node. Unassigned states have zero derivative. For pointwise models that also support independent-grid error checks and structured derivatives, use DynamicsModel.

Examples

Dynamics class must be used in dynamics function. It is an example of the equation of motion of thrust and free fall. Thrust is controllable.

$$\dot{x} = v, \qquad \dot{v} = T/m - g$$

>>> def dynamics(prob, obj, section):
>>>     x = prob.states(0, section)
>>>     v = prob.states(1, section)
>>>     T = prob.controls(0, section)
>>>     g = 9.8
>>>     m = 1.0
>>>     dx = Dynamics(prob, section)
>>>     dx[0] = v
>>>     dx[1] = T / m - g
>>>     return dx()

number_of_stateattributeinstance attribute#

View source

number_of_state = prob.number_of_states[self.section]

sectionattributeinstance attribute#

View source

section = index(section, prob.number_of_section, 'section')

unit_statesattributeinstance attribute#

View source

unit_states = prob.unit_states

unit_timeattributeinstance attribute#

View source

unit_time = prob.unit_time

DynamicsModelclass#

Re-exported from OpenGoddard._modelsView source

class DynamicsModel(function)

Pointwise physical dynamics function(x, u, t, obj, section, xp).

x has shape (states, nodes), u (controls, nodes), and t (nodes,). Return an array with shape (states, nodes). Each output column must depend ONLY on the corresponding input column. Use xp for array operations; it is NumPy for values and autograd.numpy for automatic differentiation. Reductions or interpolation across nodes violate this contract.

Instances are callable with the legacy (problem, obj, section) signature. Local derivatives use simultaneous central differences, or reverse-mode AD of each state’s sum over independent nodes. No full dense dynamics Jacobian is formed. Non-smooth expressions require the user’s choice of formulation.

functionattributeinstance attribute#

View source

function = function

inputsmethodstaticmethod#

View source

def inputs(problem, section)

linearizemethod#

View source

def linearize(problem, obj, section, mode)

physicalmethod#

View source

def physical(x, u, time, obj, section)

Guessclass#

Re-exported from OpenGoddard._helpersView source

class Guess

Construct physical initial guesses on a supplied time grid.

constantmethodclassmethod#

View source

def constant(time, const: float) -> NDArray[np.float64]

Return one constant physical value per time node.

Parameters

Name	Type	Description
time	array_like	Time nodes; only their count is used.
const	float	Constant physical value.

Returns

ndarrayNDArray[np.float64]

Values of shape (N,).

cubicmethodclassmethod#

View source

def cubic(time, y0: float, yprime0: float, yf: float, yprimef: float) -> NDArray[np.float64]

Interpolate endpoint values and physical-time slopes with a cubic.

Parameters

Name	Type	Description
time	array_like	Finite, nondecreasing physical times with a positive span. Repeated interface times are supported.
y0	float	Physical value at the first time.
yprime0	float	Physical dy/dt at the first time.
yf	float	Physical value at the last time.
yprimef	float	Physical dy/dt at the last time.

Returns

ndarrayNDArray[np.float64]

Values of shape (N,).

linearmethodclassmethod#

View source

def linear(time, y0: float, yf: float) -> NDArray[np.float64]

Interpolate physical endpoint values linearly in time.

Parameters

Name	Type	Description
time	array_like	Finite, nondecreasing physical times with a positive span. Repeated interface times are supported.
y0	float	Physical value at the first time.
yf	float	Physical value at the last time.

Returns

ndarrayNDArray[np.float64]

Values of shape (N,).

plotmethodclassmethod#

View source

def plot(x, y, title='', xlabel='', ylabel='')

Plot physical initial-guess values against time.

Parameters

Name	Type	Default	Description
x	array_like	required	array on the horizontal axis of the plot
y	array_like	required	array on the vertical axis of the plot
title	str	''	title
xlabel	str	''	xlabel
ylabel	str	''	ylabel

zerosmethodclassmethod#

View source

def zeros(time) -> NDArray[np.float64]

Return one zero per time node.

Parameters

Name	Type	Description
time	array_like	Time nodes; only their count is used.

Returns

ndarrayNDArray[np.float64]

Values of shape (N,).

LinearCallbackclass#

Re-exported from OpenGoddard._modelsView source

class LinearCallback(function)

Declare a cost/boundary/path callback affine in the scaled vector p.

Its Jacobian is computed once per solve and cached. This is an explicit mathematical contract, not a heuristic: nonlinear expressions must not be wrapped. Coefficients, units and obj must remain constant during a solve.

functionattributeinstance attribute#

View source

function = function

Problemclass#

Re-exported from OpenGoddard.optimizeView source

class Problem(time_init, nodes, number_of_states, number_of_controls, maxIterator=100, method='LGL', *, mesh=None)

Define a physical optimal-control problem and solve its LGL discretization.

Parameters

Name	Type	Default	Description
time_init	sequence of float	required	Strictly increasing physical phase times.
nodes	sequence of int	required	LGL node counts, at least two per phase.
number_of_states	sequence of int	required	Positive state counts per phase.
number_of_controls	sequence of int	required	Nonnegative control counts per phase.
maxIterator	int	100	Maximum solver attempts, including the first attempt. Another attempt is made only after an iteration-limit status. solve(maxiter=...) limits iterations within each attempt.
method	str	'LGL'	Only "LGL" is supported for optimal control.
mesh	optional	None	Per-section lists of (left, right, nodes) elements partitioning the phase coordinate [-1, 1]. Logical section numbers do not change. Node counts include both sides of internal interfaces.

All state, control and time getters and setters use physical units. Only p and bounds use scaled units. Setters for units preserve physical guesses and bounds. Node, state, control and section indices follow Python’s negative indexing convention and raise IndexError outside their range.

Attributes

Name	Type	Description
nodes	list[int]	Node counts per phase, including duplicated interfaces.
number_of_states	list[int]	State counts per phase.
number_of_controls	list[int]	Control counts per phase.
number_of_section	int	Number of logical phases.
number_of_param	ndarray	State plus control counts per phase.
number_of_variables	int	Length of p, including phase end times.
div	list[list[int]]	Exclusive ends of state/control blocks in p.
tau	list[ndarray]	LGL phase coordinates in [-1, 1].
w	list[ndarray]	Quadrature weights in phase coordinates.
D	list[ndarray]	Differentiation matrices in phase coordinates.
time	list[ndarray]	Physical node times, derived from current p.
time_all_section	ndarray	Times concatenated in phase order.
unit_states	list[list[float]]	State reference units per phase.
unit_controls	list[list[float]]	Control reference units per phase.
unit_time	float	Time reference unit.
p	ndarray	Scaled states, controls and phase end times, shape (N,).
bounds	list[tuple]	Scaled lower/upper bounds for every entry of p.
dynamics	list[Callable]	One scaled-derivative callback per phase, called as dynamics(problem, obj, section). Use Dynamics to assemble physical derivatives, or DynamicsModel for a pointwise physical model supporting derivatives and error estimation.
knot_states_smooth	list[bool]	Enforce continuity at phase knots with equal state counts. Disable it when supplying jump conditions.
cost	Callable | None	Terminal cost, called as cost(problem, obj).
running_cost	Callable | None	Physical-time integrand, shape (sum(nodes),), called as running_cost(problem, obj).
cost_derivative	function	Complete objective gradient in p.
equality	Callable | None	One-dimensional residuals equal to zero.
inequality	Callable | None	One-dimensional residuals >= 0.
array_module		NumPy; differentiable callback views supply autograd.numpy.
metadata	dict	User metadata, shallow-copied on refine/remesh.
result	OptimizeResult | None	Latest solver result, or None before solve.
result_history	list[OptimizeResult]	Attempts of the latest solve.

Callback views share states/controls (including all-section getters), time_start/time_final, time_knots, time/time_all_section, dimensions, units, array_module and condition(). Use this reading interface for portable callbacks; mutations and direct p indices do not support path sampling. Use unit setters and refine/remesh rather than editing configuration lists. Construct a new Problem to change state/control counts. Refine/remesh create a new mesh and reset solver results.

Dattributeinstance attribute#

View source

D = [rule[2] for rule in rules]

array_moduleattributeclass attributeinstance attribute#

View source

array_module = np

boundsattributeinstance attribute#

View source

bounds = [(None, None)] * self.number_of_variables

costattributeinstance attribute#

View source

cost: Callable[[Problem, Any], float] | None = None

cost_derivativeattributeinstance attribute#

View source

cost_derivative: Callable[[Problem, Any], NDArray[np.float64]] | None = None

divattributeinstance attribute#

View source

div = [list(phase) for phase in self._layout.divisions]

dynamicsattributeinstance attribute#

View source

dynamics: list[Callable[[Problem, Any, int], NDArray[np.float64]] | None] = [None] * self.number_of_section

equalityattributeinstance attribute#

View source

equality: Callable[[Problem, Any], NDArray[np.float64]] | None = None

inequalityattributeinstance attribute#

View source

inequality: Callable[[Problem, Any], NDArray[np.float64]] | None = None

iteratorattributeinstance attribute#

View source

iterator = 0

knot_states_smoothattributeinstance attribute#

View source

knot_states_smooth = [True] * self.number_of_section - 1

maxIteratorattributeinstance attribute#

View source

maxIterator = integer(maxIterator, 'maxIterator', 1)

meshattributeinstance attribute#

View source

mesh = normalize_mesh(self.nodes, mesh)

metadataattributeinstance attribute#

View source

metadata: dict[str, Any] = {}

methodattributeinstance attribute#

View source

method = method

nodesattributeinstance attribute#

View source

nodes = [integer(n, 'nodes', 2) for n in nodes]

number_of_controlsattributeinstance attribute#

View source

number_of_controls = [integer(n, 'number_of_controls') for n in number_of_controls]

number_of_paramattributeinstance attribute#

View source

number_of_param = np.array(self.number_of_states) + self.number_of_controls

number_of_sectionattributeinstance attribute#

View source

number_of_section = len(self.nodes)

number_of_statesattributeinstance attribute#

View source

number_of_states = [integer(n, 'number_of_states', 1) for n in number_of_states]

number_of_variablesattributeinstance attribute#

View source

number_of_variables = self._layout.size

pattributeinstance attribute#

View source

p = np.zeros(self.number_of_variables)

resultattributeinstance attribute#

View source

result: optimize.OptimizeResult | None = None

result_historyattributeinstance attribute#

View source

result_history: list[optimize.OptimizeResult] = []

running_costattributeinstance attribute#

View source

running_cost: Callable[[Problem, Any], NDArray[np.float64]] | None = None

t0attributeinstance attribute#

View source

t0 = float(times[0])

tauattributeinstance attribute#

View source

tau = [rule[0] for rule in rules]

time_initattributeinstance attribute#

View source

time_init = times.tolist()

unit_controlsattributeinstance attribute#

View source

unit_controls = [[1.0] * n for n in self.number_of_controls]

unit_statesattributeinstance attribute#

View source

unit_states = [[1.0] * n for n in self.number_of_states]

unit_timeattributeinstance attribute#

View source

unit_time = 1.0

wattributeinstance attribute#

View source

w = [rule[1] for rule in rules]

timeproperty#

View source

time: list[NDArray[np.float64]]

Physical node times per phase, derived from the current p.

Returned arrays are independent values. Change phase end times with set_time_final(); editing these arrays does not change the problem. Both sides of phase and element interfaces are included.

time_all_sectionproperty#

View source

time_all_section: NDArray[np.float64]

Physical node times concatenated in phase order, including interfaces.

conditionmethod#

View source

def condition(length: int = 0) -> Condition

Build equality or inequality residuals in this view’s array backend.

Parameters

Name	Type	Default	Description
length	int	0	Number of initial zero residuals, normally zero.

Returns

ConditionCondition

A residual builder. Use equal() for equality callbacks; lower_bound()/upper_bound() produce nonnegative feasible residuals.

controlsmethod#

View source

def controls(control: int, section: int) -> NDArray[np.float64]

Return one control’s physical values in a phase.

Parameters

Name	Type	Description
control	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.

Returns

ndarrayNDArray[np.float64]

Independent array of shape (nodes[section],). Both sides of internal element interfaces are present. Use set_controls() to edit.

controls_all_sectionmethod#

View source

def controls_all_section(control: int) -> NDArray[np.float64]

Return a control’s physical values concatenated in phase order.

Parameters

Name	Type	Description
control	int	Variable index, which must exist in every phase.

Returns

ndarrayNDArray[np.float64]

Shape (sum(nodes),), including duplicated interface nodes.

estimate_errormethod#

View source

def estimate_error(obj: Any = None, check_factor: int = 2) -> dict[str, Any]

Estimate dimensionless residuals between the LGL nodes.

Parameters

Name	Type	Default	Description
obj	Any	None	User data passed to the model and inequality callback.
check_factor	int	2	Integer >= 2. Check each element at check_factor * element.nodes + 1 Gauss points.

Returns

dictdict[str, Any]

section_errors (per-phase ODE indicators), max_error (their maximum), path_violation (scaled inequality/box-bound violation), and elements (per-element error, smoothness, overshoot and split point). Indicators are residual estimates, not rigorous state-error bounds. Nonlinear path constraints are checked by sampling, not certified.

Raises

ValueError

A phase lacks DynamicsModel, has nonpositive duration, or uses varying nodal box bounds.

index_controlsmethod#

View source

def index_controls(control: int, section: int, index: int | None = None) -> int

Return the index of a control node in the scaled vector p.

Parameters

Name	Type	Default	Description
control	int	required	Variable index, supporting negative indices.
section	int	required	Phase index, supporting negative indices.
index	int | None	None	Local node index, supporting negative indices; None selects the first node.

Returns

intint

Index into p, whose values use reference-unit scaling.

index_statesmethod#

View source

def index_states(state: int, section: int, index: int | None = None) -> int

Return the index of a state node in the scaled vector p.

Parameters

Name	Type	Default	Description
state	int	required	Variable index, supporting negative indices.
section	int	required	Phase index, supporting negative indices.
index	int | None	None	Local node index, supporting negative indices; None selects the first node.

Returns

intint

Index into p, whose values use reference-unit scaling.

index_time_finalmethod#

View source

def index_time_final(section: int) -> int

Return the index of a phase end time in the scaled vector p.

method_LGmethod#

View source

def method_LG(n) -> tuple[NDArray[np.float64], ...]

Legendre-Gauss Pseudospectral method Gauss nodes are roots of $P_n(x)$.

Parameters

Name	Type	Description
n	int	number of nodes

Returns

tuple[NDArray[np.float64], ...]

tuple[ndarray, ndarray, ndarray]: Nodes (n,), weights (n,), and D (n, n).

method_LGLmethod#

View source

def method_LGL(n) -> tuple[NDArray[np.float64], ...]

Legendre-Gauss-Lobatto Pseudospectral method Gauss-Lobatto nodes are roots of $P'_{n-1}(x)$, plus -1 and 1.

Parameters

Name	Type	Description
n	int	number of nodes

Returns

tuple[NDArray[np.float64], ...]

tuple[ndarray, ndarray, ndarray]: Nodes (n,), weights (n,), and D (n, n).

References

Fariba Fahroo and I. Michael Ross. “Advances in Pseudospectral Methods for Optimal Control”, AIAA Guidance, Navigation and Control Conference and Exhibit, Guidance, Navigation, and Control and Co-located Conferences http://dx.doi.org/10.2514/6.2008-7309

method_LGRmethod#

View source

def method_LGR(n) -> tuple[NDArray[np.float64], ...]

Legendre-Gauss-Radau Pseudospectral method Gauss-Radau nodes are roots of $P_n(x) + P_{n-1}(x)$.

Parameters

Name	Type	Description
n	int	number of nodes

Returns

tuple[NDArray[np.float64], ...]

tuple[ndarray, ndarray, ndarray]: Nodes (n,), weights (n,), and D (n, n).

plotmethod#

View source

def plot(title_comment='')

plot inner variables that to be optimized

Parameters

Name	Type	Default	Description
title_comment	str	''	string for title

refinemethod#

View source

def refine(nodes) -> Problem

Return a new problem with interpolated guesses and new node counts.

Parameters

Name	Type	Description
nodes	sequence of int	Total node counts per phase, including duplicated interfaces. At least two nodes per element are required.

Returns

ProblemProblem

A new problem preserving element boundaries, physical units, constant box bounds and callbacks. Solver results are reset. Metadata containers are shallow-copied; nested objects remain shared.

Raises

ValueError

Counts are invalid or box bounds vary between nodes.

remeshmethod#

View source

def remesh(mesh) -> Problem

Return a warm start on per-phase lists of (left, right, nodes) elements.

Element endpoints are phase coordinates in [-1, 1]. Logical section indices remain unchanged. States are continuous at internal interfaces; controls have independent left and right limits.

Parameters

Name	Type	Description
mesh	sequence	Per-phase lists of (left, right, nodes) elements partitioning [-1, 1]. Each element includes both endpoints.

Returns

ProblemProblem

A new problem with interpolated guesses, preserved physical units, constant bounds and callbacks, and reset solver results. Metadata containers are shallow-copied. Private extension caches are not copied.

Raises

ValueError

The mesh is invalid or box bounds vary between nodes.

set_controlsmethod#

View source

def set_controls(control: int, section: int, value: ArrayLike) -> None

Set one control’s initial guess in physical units.

Parameters

Name	Type	Description
control	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
value	ArrayLike	Finite one-dimensional values of shape (nodes[section],).

Raises

ValueError

Values have the wrong shape or are not finite.

set_controls_all_sectionmethod#

View source

def set_controls_all_section(control: int, value_all_section: ArrayLike) -> None

Set a control’s physical initial guess across all phases.

Parameters

Name	Type	Description
control	int	Variable index, which must exist in every phase.
value_all_section	ArrayLike	Finite values of shape (sum(nodes),), in phase order and including duplicated interface nodes.

set_controls_boundsmethod#

View source

def set_controls_bounds(control: int, section: int, lb: float | None, ub: float | None) -> None

Set constant physical box bounds at every node of a phase.

Parameters

Name	Type	Description
control	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
lb	float | None	Physical lower bound, or None for no lower bound.
ub	float | None	Physical upper bound, or None for no upper bound.

Notes

Bounds constrain mesh nodes. Use solve_adaptive() with a supported model to check for violations between nodes.

set_controls_bounds_all_sectionmethod#

View source

def set_controls_bounds_all_section(control: int, lb: float | None, ub: float | None) -> None

set value to bounds of control at all sections

Parameters

Name	Type	Description
control	int	control
lb	float or None	lower bound
ub	float or None	upper bound

set_statesmethod#

View source

def set_states(state: int, section: int, value: ArrayLike) -> None

Set one state’s initial guess in physical units.

Parameters

Name	Type	Description
state	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
value	ArrayLike	Finite one-dimensional values of shape (nodes[section],).

Raises

ValueError

Values have the wrong shape or are not finite.

set_states_all_sectionmethod#

View source

def set_states_all_section(state: int, value_all_section: ArrayLike) -> None

Set a state’s physical initial guess across all phases.

Parameters

Name	Type	Description
state	int	Variable index, which must exist in every phase.
value_all_section	ArrayLike	Finite values of shape (sum(nodes),), in phase order and including duplicated interface nodes.

set_states_boundsmethod#

View source

def set_states_bounds(state: int, section: int, lb: float | None, ub: float | None) -> None

Set constant physical box bounds at every node of a phase.

Parameters

Name	Type	Description
state	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
lb	float | None	Physical lower bound, or None for no lower bound.
ub	float | None	Physical upper bound, or None for no upper bound.

Notes

Bounds constrain mesh nodes. Use solve_adaptive() with a supported model to check for violations between nodes.

set_states_bounds_all_sectionmethod#

View source

def set_states_bounds_all_section(state: int, lb: float | None, ub: float | None) -> None

set value to bounds of state at all sections

Parameters

Name	Type	Description
state	int	state
lb	float or None	lower bound
ub	float or None	upper bound

set_time_finalmethod#

View source

def set_time_final(section: int, value: float) -> None

Set a phase’s end-time guess in physical units.

Parameters

Name	Type	Description
section	int	Phase index, supporting negative indices.
value	float	Finite physical end time. Times are immediately reflected in time and time_all_section. Forward ordering is enforced at solve.

Raises

ValueError

The end time is not finite.

set_time_final_boundsmethod#

View source

def set_time_final_bounds(section: int, lb: float | None, ub: float | None) -> None

set value to bounds of time_final at specific section

Parameters

Name	Type	Description
section	int	section
lb	float or None	lower bound
ub	float or None	upper bound

set_unit_controlsmethod#

View source

def set_unit_controls(control: int, section: int, value: float) -> None

Change a control’s reference unit, preserving physical values and bounds.

Parameters

Name	Type	Description
control	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
value	float	Finite, strictly positive reference unit.

set_unit_controls_all_sectionmethod#

View source

def set_unit_controls_all_section(control: int, value: float) -> None

set a canonical unit value to the control at all sections

Parameters

Name	Type	Description
control	int	control
value	float	value

set_unit_statesmethod#

View source

def set_unit_states(state: int, section: int, value: float) -> None

Change a state’s reference unit, preserving physical values and bounds.

Parameters

Name	Type	Description
state	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.
value	float	Finite, strictly positive reference unit.

set_unit_states_all_sectionmethod#

View source

def set_unit_states_all_section(state: int, value: float) -> None

set a canonical unit value to the state at all sections

Parameters

Name	Type	Description
state	int	state
value	float	value

set_unit_timemethod#

View source

def set_unit_time(value: float) -> None

set a canonical unit value to the time

Parameters

Name	Type	Description
value	float	value

solvemethod#

View source

def solve(obj: Any = None, display_func: Callable[[], None] = _dummy_func, *, solver: str = 'SLSQP', derivatives: str = 'finite-difference', **options: Any) -> optimize.OptimizeResult

Solve the LGL transcription and return a SciPy OptimizeResult.

solver selects SLSQP (default) or sparse trust-constr. derivatives='finite-difference' preserves arbitrary legacy callbacks. 'structured' assembles exact spectral blocks and local central differences, and 'autograd' uses optional automatic differentiation. Both structured modes require DynamicsModel in every section. AD callbacks use view.array_module and view.condition(); they must not mutate arrays or coerce traced values to Python floats. With trust-constr, AD also supplies matrix-free Hessian-vector products.

dynamics must contain a callable for each section; equality and inequality may be None. running_cost returns one integrand value per node, concatenated in section order, in physical time. The solver applies each section’s (tf - ti) / 2 quadrature factor.

cost_derivative, when set, must differentiate the complete objective (including running cost) with respect to the scaled vector p. All SLSQP options, including ftol, maxiter, eps and disp, are forwarded. Concurrent workers are unsupported because callbacks share this mutable Problem instance.

Solver-specific options are forwarded: SLSQP accepts ftol/eps; trust-constr accepts gtol/xtol/barrier_tol. They are not interchangeable. cost_derivative is used in finite-difference and structured modes; autograd differentiates the complete objective itself.

maxIterator limits total attempts, including the first. A restart occurs only after an iteration limit (SLSQP status 9, trust-constr status 0). Other failures stop immediately and emit RuntimeWarning. Each call resets the restart counter. The final result is also saved in result and every attempt in result_history. display_func is called without arguments after each attempt, as in the legacy API. time and time_all_section always reflect the current parameters.

Parameters

Name	Type	Default	Description
obj	Any	None	User data passed unchanged to every model/cost/constraint callback.
display_func	Callable[[], None]	_dummy_func	Called without arguments after each solver attempt.
solver	str	'SLSQP'	“SLSQP” or “trust-constr”.
derivatives	str	'finite-difference'	“finite-difference”, “structured”, or “autograd”.
**options	Any	{}	SciPy solver options, e.g. maxiter, ftol (SLSQP) or gtol (trust-constr).

Returns

OptimizeResultoptimize.OptimizeResult

SciPy result with solver, derivatives and scaled equality_violation/inequality_violation diagnostics. The solution also replaces this problem’s p. Failure returns success=False and emits RuntimeWarning; it does not raise a convergence exception.

solve_adaptivemethod#

View source

def solve_adaptive(obj: Any = None, *, tolerance: float = 1e-05, max_refinements: int = 4, max_nodes: int = 100, growth: float = 1.5, strategy: str = 'p', max_elements: int = 24, min_element_nodes: int = 6, max_order: int = 14, path_check_factor: int = 4, **solve_options: Any) -> optimize.OptimizeResult

Solve and adapt the mesh, updating this problem to the final mesh.

Requires DynamicsModel in every phase. Callbacks must use physical getters and work with changing node counts. For h/hp path sampling, avoid direct p indices and mesh-dependent global reductions.

Parameters

Name	Type	Default	Description
obj	Any	None	User data passed to all callbacks.
tolerance	float	1e-05	Positive tolerance for ODE and scaled path residuals.
max_refinements	int	4	Maximum changes after the initial solve.
max_nodes	int	100	Maximum total nodes per phase, including interfaces.
growth	float	1.5	Node growth factor > 1 for p refinement.
strategy	str	'p'	“p” increases orders; “h” splits elements; “hp” chooses between increasing orders and splitting using error indicators.
max_elements	int	24	Maximum elements per phase for h/hp.
min_element_nodes	int	6	Minimum child node count for hp, at least four.
max_order	int	14	Node threshold for choosing h instead of p in hp.
path_check_factor	int	4	Sampling density for h/hp path constraints, at least two. Not used by the p strategy.
**solve_options	Any	{}	Options forwarded to solve(), including solver, derivatives and SciPy tolerances.

Returns

OptimizeResultoptimize.OptimizeResult

Final solver result. success reports optimizer convergence; mesh_converged separately reports meeting the mesh tolerance. mesh_history records every solve and residual estimate. h/hp also supply mesh_strategy and mesh_stop_reason. Stopping at a limit emits RuntimeWarning when the mesh has not converged.

statesmethod#

View source

def states(state: int, section: int) -> NDArray[np.float64]

Return one state’s physical values in a phase.

Parameters

Name	Type	Description
state	int	Variable index, supporting negative indices.
section	int	Phase index, supporting negative indices.

Returns

ndarrayNDArray[np.float64]

Independent array of shape (nodes[section],). Both sides of internal element interfaces are present. Use set_states() to edit.

states_all_sectionmethod#

View source

def states_all_section(state: int) -> NDArray[np.float64]

Return a state’s physical values concatenated in phase order.

Parameters

Name	Type	Description
state	int	Variable index, which must exist in every phase.

Returns

ndarrayNDArray[np.float64]

Shape (sum(nodes),), including duplicated interface nodes.

time_finalmethod#

View source

def time_final(section: int) -> float

Return a phase’s physical end time from the current scaled p.

Parameters

Name	Type	Description
section	int	Phase index; negative indices count from the last phase.

Returns

float

Physical time, equal to p[index_time_final(section)] * unit_time.

time_final_all_sectionmethod#

View source

def time_final_all_section() -> list[float]

Return physical phase end times in phase order.

time_knotsmethod#

View source

def time_knots() -> list[float]

Return physical phase boundaries [t0, tf0, tf1, …].

time_startmethod#

View source

def time_start(section: int) -> float

Return a phase’s physical start time.

Parameters

Name	Type	Description
section	int	Phase index; negative indices count from the last phase.

Returns

float

Physical time. The first phase starts at the fixed t0; later phases start at the preceding phase’s current end time.

time_to_taumethod#

View source

def time_to_tau(time)

Map a finite, nonconstant time vector to [-1, 1].

time_updatemethod#

View source

def time_update() -> NDArray[np.float64]

Return current physical node times concatenated in phase order.

Kept for compatibility. Times are derived from p on access, so explicit refreshes are unnecessary even after direct changes to p.

to_csvmethod#

View source

def to_csv(filename='OpenGoddard_output.csv', delimiter: str = ',') -> None

Write physical times, states and controls in phase/node order.

Parameters

Name	Type	Default	Description
filename	str or PathLike	'OpenGoddard_output.csv'	Output path. Missing components in phases with fewer states or controls are written as NaN.
delimiter	str	','	Column separator; a comma by default.