diff --git a/CMakePresets.json b/CMakePresets.json index 8fc6d55ead..ecae983934 100644 --- a/CMakePresets.json +++ b/CMakePresets.json @@ -9,22 +9,11 @@ { "name": "core", "hidden": true, - "generator": "Ninja", "binaryDir": "${sourceDir}/out/build/${presetName}", "installDir": "${sourceDir}/out/install/${presetName}" }, { - "name": "common", - "hidden": true, - "cacheVariables": { - "WITH_IPOPT": "ON", - "WITH_OPENMP": "ON", - "WITH_QPOASES": "ON", - "WITH_LAPACK": "ON" - } - }, - { - "name": "devel", + "name": "debug", "hidden": true, "cacheVariables": { "CMAKE_BUILD_TYPE": "Debug", @@ -39,16 +28,10 @@ "CMAKE_BUILD_TYPE": "Release" } }, - { - "name": "python", - "hidden": true, - "cacheVariables": { - "WITH_PYTHON": "ON" - } - }, { "name": "mingw", "hidden": true, + "generator": "Ninja", "cacheVariables": { "WITH_THREAD_MINGW": "ON", "CMAKE_C_COMPILER": "gcc", @@ -58,6 +41,7 @@ { "name": "cl", "hidden": true, + "generator": "Ninja", "cacheVariables": { "CMAKE_C_COMPILER": "cl", "CMAKE_CXX_COMPILER": "cl" @@ -66,6 +50,7 @@ { "name": "clang-cl", "hidden": true, + "generator": "Ninja", "cacheVariables": { "CMAKE_C_COMPILER": "clang-cl", "CMAKE_CXX_COMPILER": "clang-cl" @@ -76,6 +61,38 @@ } } }, + { + "name": "python", + "hidden": true, + "cacheVariables": { + "WITH_PYTHON": "ON" + } + }, + { + "name": "common", + "hidden": true, + "cacheVariables": { + "WITH_IPOPT": "ON", + "WITH_OPENMP": "ON", + "WITH_QPOASES": "ON", + "WITH_LAPACK": "ON" + } + }, + { + "name": "devel", + "displayName": "Linux developer", + "inherits": [ + "debug", + "python", + "core", + "common" + ], + "cacheVariables": { + "PYTHON_PREFIX": "${sourceDir}/out/install/${presetName}/python", + "WITH_CLANG": "ON", + "WITH_BUILD_REQUIRED": "ON" + } + }, { "name": "windevel", "displayName": "Windows developer (MinGW)", @@ -83,7 +100,7 @@ "common", "mingw", "core", - "devel", + "debug", "python" ] }, @@ -93,7 +110,7 @@ "inherits": [ "cl", "core", - "devel" + "debug" ] } ] diff --git a/docs/examples/cplusplus/multiple_shooting_from_scratch.cpp b/docs/examples/cplusplus/multiple_shooting_from_scratch.cpp index e542de6446..8392693e42 100644 --- a/docs/examples/cplusplus/multiple_shooting_from_scratch.cpp +++ b/docs/examples/cplusplus/multiple_shooting_from_scratch.cpp @@ -79,7 +79,7 @@ int main(){ SXDict dae = {{"x", x}, {"p", u}, {"ode", ode}, {"quad", quad}}; // Create an integrator (CVodes) - Function F = integrator("integrator", "cvodes", dae, {{"t0", 0}, {"tf", tf/ns}}); + Function F = integrator("integrator", "cvodes", dae, 0, tf/ns); // Total number of NLP variables int NV = nx*(ns+1) + nu*ns; diff --git a/docs/examples/cplusplus/rocket_single_shooting.cpp b/docs/examples/cplusplus/rocket_single_shooting.cpp index 129f4db38a..986e3a873a 100644 --- a/docs/examples/cplusplus/rocket_single_shooting.cpp +++ b/docs/examples/cplusplus/rocket_single_shooting.cpp @@ -109,11 +109,9 @@ int main(){ opts["interpolation_order"] = 1; opts["number_of_finite_elements"] = 1000; } - opts["t0"] = t0; - opts["tf"] = tf; // Create integrator - Function F = integrator("integrator", plugin, dae, opts); + Function F = integrator("integrator", plugin, dae, t0, tf, opts); // control for all segments MX U = MX::sym("U",nu); diff --git a/docs/users_guide/source/daebuilder.rst b/docs/users_guide/source/daebuilder.rst index 107edcb565..87e8c8acf3 100644 --- a/docs/users_guide/source/daebuilder.rst +++ b/docs/users_guide/source/daebuilder.rst @@ -7,99 +7,95 @@ The |DaeBuilder| class The |DaeBuilder| class in |casadi| is an auxiliary class intended to facilitate the modeling complex dynamical systems for later use with optimal -control algorithms. This class can be seen as a low-level alternative to +control algorithms. This class has lower level of abstraction than a physical modeling language such as Modelica (cf. :numref:`sec-modelica`), while still being higher level than working directly with |casadi| symbolic -expressions. Another important usage it to provide an interface to -physical modeling languages and software and be a building blocks for -developing domain specific modeling environments. +expressions. In particular, it can be used to interface physical models +available in the Functional Mock-up Interface (FMI) format or be used +as a building block for constructing a domain specific modeling environments. -Using the |DaeBuilder| class consists of the following steps: - -* Step-by-step constructing a structured system of differential-algebraic equations (DAE) or, alternatively, importing an existing model from Modelica -* Symbolically reformulate the DAE -* Generate a chosen set of |casadi| functions to be used for e.g. optimal control or C code generation - -In the following sections, we describe the mathematical formulation of the class -and its intended usage. +There are in principle three different ways to use |DaeBuilder|: + * One can use the class to step-by-step construct a structured system of differential-algebraic equations (DAE), which can then be used to interface to simulation or optimization in CasADi. There is also experimental support for exporting the models in the FMI format for use in other tools. + * One can use the class to import existing models available in the FMI format. As of CasADi 3.6, we support import of standard FMUs, where the model equations are available only by calls to DLLs. The FMU import support has been tested for challenging FMUs, but is still under development as of this writing. + * One can use the class to import existing models in a symbolic XML format. This format was supported in the JModelica.org tool and has experimental support in OpenModelica. This format is not in active development due to the limited availability of tools to generate models in the format. .. _sec-daebuilder_io: -Mathematical formulation ------------------------- +.. Mathematical formulation +.. ------------------------ -The |DaeBuilder| class uses a relatively rich problem formulation that -consists of a set of input expressions and a set of output expressions, each -defined by a string identifier. The choice of expressions was inspired by the -*functional mock-up interface* (FMI) version 2.0 -[#f3]_ +.. The |DaeBuilder| class uses a relatively rich problem formulation that +.. consists of a set of input expressions and a set of output expressions, each +.. defined by a string identifier. The choice of expressions was inspired by the +.. *functional mock-up interface* (FMI) version 2.0 +.. [#f3]_ -Input expressions -^^^^^^^^^^^^^^^^^ +.. Input expressions +.. ^^^^^^^^^^^^^^^^^ -'t' - Time :math:`t` +.. 't' +.. Time :math:`t` -'c' - Named constants :math:`c` +.. 'c' +.. Named constants :math:`c` -'p' - Independent parameters :math:`p` +.. 'p' +.. Independent parameters :math:`p` -'v' - Dependent variables :math:`v`, depends on other variables and, acyclically, on other :math:`v` +.. 'v' +.. Dependent variables :math:`v`, depends on other variables and, acyclically, on other :math:`v` -'x' - Differential state :math:`x`, defined by an explicit ODE +.. 'x' +.. Differential state :math:`x`, defined by an explicit ODE -'s' - Differential state :math:`s`, defined by an implicit ODE +.. 's' +.. Differential state :math:`s`, defined by an implicit ODE -'sdot' - Time derivatives implicitly defined differential state :math:`\dot{s}` +.. 'sdot' +.. Time derivatives implicitly defined differential state :math:`\dot{s}` -'z' - Algebraic variable, defined by an algebraic equation +.. 'z' +.. Algebraic variable, defined by an algebraic equation -'q' - Quadrature state :math:`q`. A differential state that may not appear in - the right-hand-side and hence can be calculated by quadrature formulas. +.. 'q' +.. Quadrature state :math:`q`. A differential state that may not appear in +.. the right-hand-side and hence can be calculated by quadrature formulas. -'y' - Output variables :math:`y` +.. 'y' +.. Output variables :math:`y` -'u' - Input variables :math:`u` +.. 'u' +.. Input variables :math:`u` -Output expressions -^^^^^^^^^^^^^^^^^^ +.. Output expressions +.. ^^^^^^^^^^^^^^^^^^ -The above input expressions are used to define the following output expressions: +.. The above input expressions are used to define the following output expressions: -'vdef' - Explicit expression for calculating :math:`v` +.. 'vdef' +.. Explicit expression for calculating :math:`v` -'ode' - The explicit ODE right-hand-side: - :math:`\dot{x} = \text{ode}(t,v,x,s,z,u,p)` +.. 'ode' +.. The explicit ODE right-hand-side: +.. :math:`\dot{x} = \text{ode}(t,v,x,s,z,u,p)` -'dae' - The implicit ODE right-hand-side: - :math:`\text{dae}(t,v,x,s,z,u,p,\dot{s}) =0` +.. 'dae' +.. The implicit ODE right-hand-side: +.. :math:`\text{dae}(t,v,x,s,z,u,p,\dot{s}) =0` -'alg' - The algebraic equations: - :math:`\text{alg}(t,v,x,s,z,u,p) = 0` +.. 'alg' +.. The algebraic equations: +.. :math:`\text{alg}(t,v,x,s,z,u,p) = 0` -'quad' - The quadrature equations: - :math:`\dot{q} = \text{quad}(t,v,x,s,z,u,p)` +.. 'quad' +.. The quadrature equations: +.. :math:`\dot{q} = \text{quad}(t,v,x,s,z,u,p)` -'ydef' - Explicit expressions for calculating :math:`y` +.. 'ydef' +.. Explicit expressions for calculating :math:`y` .. _sec-daebuilder_syntax: @@ -188,9 +184,17 @@ list of functions, see the C++ API documentation for |DaeBuilder|. .. _sec-modelica: -Import of OCPs from Modelica +Symbolic import of OCPs from Modelica ---------------------------- +Note: JModelica.org is no longer offered by Modelon. Its closed-source successor code, OCT, +does retain CasADi interoperability however. Details of how to use OCT to generate CasADi +expressions is best described in OCT's user guide. The text in the following refers to the +legacy Modelica interoperability support. + +Legacy symbolic import from XML files +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + An alternative to model directly in |casadi|, as above, is to use an advanced physical modeling language such as Modelica to specify the model. For this, |casadi| offers interoperability with the open-source `JModelica.org `_ compiler, which @@ -237,35 +241,35 @@ using the ``parse_fmi`` command: ocp.parse_fmi('modelDescription.xml') -Symbolic reformulation ----------------------- +.. Symbolic reformulation +.. ---------------------- -One of the original purposes of the |DaeBuilder| class was to reformulate -a *fully-implicit DAE*, typically coming from Modelica, to a semi-explicit -DAE that can be used more readily in optimal control algorithms. +.. One of the original purposes of the |DaeBuilder| class was to reformulate +.. a *fully-implicit DAE*, typically coming from Modelica, to a semi-explicit +.. DAE that can be used more readily in optimal control algorithms. -This can be done by the \python{make_implicit} command: +.. This can be done by the \python{make_implicit} command: -.. side-by-side:: - .. code-block:: python +.. .. side-by-side:: +.. .. code-block:: python - ocp.make_explicit() +.. ocp.make_explicit() - && +.. && - .. code-block:: octave +.. .. code-block:: octave - ocp.make_explicit(); +.. ocp.make_explicit(); -Other useful commands available for an instance ``ocp`` of |DaeBuilder| include: +.. Other useful commands available for an instance ``ocp`` of |DaeBuilder| include: -* ``print ocp`` Print the optimal optimal control problem to screen -* ``ocp.scale_variables()`` Scale all variables using the *nominal* attribute for each variable -* ``ocp.eliminate_d()`` Eliminate all independent parameters from the symbolic expressions +.. * ``print ocp`` Print the optimal optimal control problem to screen +.. * ``ocp.scale_variables()`` Scale all variables using the *nominal* attribute for each variable +.. * ``ocp.eliminate_d()`` Eliminate all independent parameters from the symbolic expressions -For a more detailed description of this class and its functionalities, we again -refer to the API documentation. +.. For a more detailed description of this class and its functionalities, we again +.. refer to the API documentation. Function factory ---------------- diff --git a/docs/users_guide/source/function.rst b/docs/users_guide/source/function.rst index 49efa398a5..c4198eb20a 100644 --- a/docs/users_guide/source/function.rst +++ b/docs/users_guide/source/function.rst @@ -257,14 +257,14 @@ is a DAE of semi-explicit form with quadratures: .. math:: \begin{aligned} - \dot{x} &= f_{\text{ode}}(t,x,z,p), \qquad x(0) = x_0 \\ - 0 &= f_{\text{alg}}(t,x,z,p) \\ - \dot{q} &= f_{\text{quad}}(t,x,z,p), \qquad q(0) = 0 + \dot{x} &= f_{\text{ode}}(t,x,z,p,u), \qquad x(0) = x_0 \\ + 0 &= f_{\text{alg}}(t,x,z,p,u) \\ + \dot{q} &= f_{\text{quad}}(t,x,z,p,u), \qquad q(0) = 0 \end{aligned} For solvers of *ordinary* differential equations, the second equation and the algebraic variables :math:`z` must be absent. -An integrator in |casadi| is a function that takes the state at the initial time ``x0``, a set of parameters ``p``, and a guess for the algebraic variables (only for DAEs) ``z0`` and returns the state vector ``xf``, algebraic variables ``zf`` and the quadrature state ``qf``, all at the final time. +An integrator in |casadi| is a function that takes the state at the initial time ``x0``, a set of parameters ``p`` and controls ``u``, and a guess for the algebraic variables (only for DAEs) ``z0`` and returns the state vector ``xf``, algebraic variables ``zf`` and the quadrature state ``qf`` at a number of output times. The control vector ``u`` is assumed to be piecewise constant and has the same grid discretization as the output grid. The freely available `SUNDIALS suite `_ (distributed along with |casadi|) contains the two popular integrators CVodes and IDAS for ODEs and DAEs respectively. These integrators have support for forward and adjoint sensitivity analysis and when used via |casadi|'s Sundials interface, |casadi| will automatically formulate the Jacobian information, which is needed by the backward differentiation formula (BDF) that CVodes and IDAS use. Also automatically formulated will be the forward and adjoint sensitivity equations. @@ -298,9 +298,9 @@ An integrator, using the "idas" plugin, can be created using the syntax: F = integrator('F', 'idas', dae); disp(F) - -Integrating this DAE from 0 to 1 with :math:`x(0)=0`, :math:`p=0.1` and using the guess :math:`z(0)=0`, can -be done by evaluating the created function object: +This will result in an integration from :math:`t_0=0` until :math:`t_f=1`, i.e. a single output time. +We can evaluate the function object using the initial condition :math:`x(0)=0`, parameter :math:`p=0.1` +and the guess for the algebraic variable at the initial time :math:`z(0)=0` as follows: .. side-by-side:: .. exec-block:: python @@ -322,7 +322,7 @@ be done by evaluating the created function object: r = F('x0',0,'z0',0,'p',0.1); disp(r.xf) -The time horizon is assumed to be fixed [#f1]_ and can be changed from its default [0, 1] by setting the options "t0" and "tf". +Note that the time horizon is always fixed. It can be changed to its default values :math:`t_0=0` and :math:`t_f=1` by adding two additional argument to the constructor, after the DAE. :math:`t_f` can either be a single value or a vector of values. It may include the initial time. Sensitivity analysis ^^^^^^^^^^^^^^^^^^^^ @@ -810,6 +810,64 @@ The first input and output are used for accumulating, while the remainder inputs The ``map``/``mapaccum`` operation exhibits a graph size and initialization time that scales logarithmically with :math:`n`. +Conditional evaluation +^^^^ + +It is possible to include conditional evaluation of expressions in CasADi expression graphs by constucting ``conditional`` function instances. This function takes a number of existing ``Function`` instances, :math:`f_1`, :math:`f_2`, :math:`f_n` as well as a "default" function :math:`f_{def}`. All these functions must have the same input and output signatures, i.e. the same number of inputs and outputs with the same dimensions: + +.. side-by-side:: + .. exec-block:: python + + x = SX.sym("x") + f0 = Function("f0",[x],[sin(x)]) + f1 = Function("f1",[x],[cos(x)]) + f2 = Function("f2",[x],[tan(x)]) + f_cond = Function.conditional('f_cond', [f0, f1], f2) + print(cond_f) + + .. exec-block:: octave + + x = SX.sym('x'); + f0 = Function('f0',{x},{sin(x)}); + f1 = Function('f1',{x},{cos(x)}); + f2 = Function('f2',{x},{tan(x)}); + f_cond = Function.conditional('f_cond', {f0, f1}, f2); + disp(F); + +The result is a new function instance with the same input/output signature, but with one additional input corresponding to an index. Its evaluation corresponds to: + +.. math:: + f_{\text{cond}}(c, x_1, x_2, \ldots, x_m) = \left\{ + \begin{array}{ll} + f_0(x_1, x_2, \ldots, x_m) & \text{if $c = 0$,} \\ + f_1(x_1, x_2, \ldots, x_m) & \text{if $c = 1$,} \\ + f_2(x_1, x_2, \ldots, x_m) & \text{otherwise} + \end{array} + \right. + +A function above can be missing (i.e. be a null pointer ``Function()``), in which case all outputs will be evaluated to NaN. Note that the evaluation is "short-circuiting", i.e. only the relevant function is evaluated. This also applies to any derivative calculation. + +A common special case is when there is only a single case in addition to the default case. This is equivalent to an if-then-else statement, which can be written with the shorthand: + +.. side-by-side:: + .. exec-block:: python + + x = SX.sym("x") + f_true = Function("f_true",[x],[cos(x)]) + f_false = Function("f_false",[x],[sin(x)]) + f_cond = Function.if_else('f_cond', f_true, f_false) + print(f_cond) + + .. exec-block:: octave + + x = SX.sym('x'); + f_true = Function('f_true',{x},{cos(x)}); + f_false = Function('f_false',{x},{sin(x)}); + f_cond = Function.if_else('f_cond', f_true, f_false); + disp(f_cond); + +Note that conditional expressions such as these can result in non-smooth expressions that may not converge if used if used in a gradient-based optimization algorithm. + .. rubric:: Footnotes .. [#f1] for problems with free end time, you can always scale time by introducing an extra parameter and substitute :math:`t` for a dimensionless time variable that goes from 0 to 1