You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
This example introduces [Solution Manager](solution-manager.html), that will work as solution observer.
25
+
This example introduces [Solution Manager](solution-manager.html), that will work as a solution observer.
26
26
27
27
## Perovskite model
28
28
29
29
Perovskite model source file example: [perovskite.cpp](https://github.com/dae-cpp/dae-cpp/blob/master/examples/perovskite_model/perovskite.cpp)
30
30
{: .fs-5 .fw-400 }
31
31
32
-
A more sophisticated example, where we solve a big DAE system that describes potential distribution and ion concentration in a perovskite solar cell.
32
+
A more sophisticated example, where we solve a big DAE system that describes the potential distribution and ion concentration in a perovskite solar cell.
33
33
We show how to add parameters to the user-defined mass matrix, vector function, and Jacobian.
34
34
This example also demonstrates a custom [Solution Manager](solution-manager.html) implementation that can work as an observer and event function.
Copy file name to clipboardExpand all lines: docs/index.md
+2-2
Original file line number
Diff line number
Diff line change
@@ -23,13 +23,13 @@ to be solved in the interval $$t \in [0, t_\mathrm{end}]$$ with the given initia
23
23
24
24
### How does it work
25
25
26
-
The DAE solver uses implicit Backward Differentiation Formulae (BDF) of orders I-IV with adaptive time stepping. Every time step, the BDF integrator reduces the original DAE system to a system of nonlinear equations, which is solved using iterative [Quasi-Newton](https://en.wikipedia.org/wiki/Quasi-Newton_method) root-finding algorithm. The Quasi-Newton method reduces the problem further down to a system of linear equations, which is solved using [Eigen](https://eigen.tuxfamily.org/index.php?title=Main_Page), a versatile and fast C++ template library for linear algebra.
26
+
The DAE solver uses implicit Backward Differentiation Formulae (BDF) of orders I-IV with adaptive time stepping. Every time step, the BDF integrator converts the original DAE system to a system of nonlinear equations, which is solved using iterative [Quasi-Newton](https://en.wikipedia.org/wiki/Quasi-Newton_method) root-finding algorithm. The Quasi-Newton method reduces the problem further down to a system of linear equations, which is solved using [Eigen](https://eigen.tuxfamily.org/index.php?title=Main_Page), a versatile and fast C++ template library for linear algebra.
27
27
Eigen's sparse solver performs two steps: factorization (decomposition) of the Jacobian matrix and the linear system solving itself. This gives us the numerical solution of the entire DAE system at the current time step. Finally, depending on the convergence rate of the Quasi-Newton method, variability of the solution, and user-defined accuracy, the DAE solver adjusts the time step size and initiates a new iteration in time.
28
28
29
29
### The main features of the solver
30
30
31
31
- Header only, no pre-compilation required.
32
-
- Uses automatic (algorithmic, exact) differentiation ([autodiff](https://autodiff.github.io/)) to compute the Jacobian matrix, if it is not provided by the user.
32
+
- Uses [automatic](https://en.wikipedia.org/wiki/Automatic_differentiation) (algorithmic, exact) differentiation (using [autodiff](https://autodiff.github.io/)) to compute the Jacobian matrix, if it is not provided by the user.
33
33
- Fourth-order implicit BDF time integrator that preserves accuracy even when the time step rapidly changes.
34
34
- A very flexible and customizable variable time stepping algorithm based on the solution stability and variability.
35
35
- Mass matrix can be non-static (can depend on time) and it can be singular.
@@ -112,7 +112,7 @@ For more information about defining the matrix in sparse format, refer to the [S
112
112
The solver provides a helper class `daecpp::JacobianAutomatic` to compute the Jacobian matrix for the given vector function $$\mathbf{f}(\mathbf{x}, t)$$ algorithmically using [`autodiff`](https://autodiff.github.io/) package.
113
113
For relatively small systems, the user does not even need to define the Jacobian, not even automatic one. This will be handled by the solver itself. Behind the scenes, the solver will create an automatic Jacobian matrix object for the user, if analytic Jacobian is not provided.
114
114
115
-
However, in some cases, it can be difficult to derive the Jacobian matrix analytically. Or if the user has provided the Jacobian matrix, but the solution diverges, which means there might be a bug in the Jacobian matrix definition. In these cases, the user can try to feed the automatic Jacobian matrix to the solver explicitly, or print out (or save to a file) both user-defined and automatic Jacobians for comparison to find possible errors in the matrix definition.
115
+
For big systems, it is highly recommended to provide analytic (user-defined) Jacobian. However, in some cases, it can be difficult to derive the Jacobian matrix analytically. Or the user has provided the Jacobian matrix, but the solution diverges, which means there might be a bug in the Jacobian matrix definition. In these cases, the user can try to feed the automatic Jacobian matrix to the solver explicitly, or print out (or save to a file) both user-defined and automatic Jacobians for comparison to find possible errors in the matrix definition.
116
116
117
117
In order to construct an automatic Jacobian matrix object, the user needs to provide the vector function object `rhs` (see [Vector Function class](vector-function.html) section):
@@ -91,7 +91,7 @@ In the example above, we created an object `identityMatrix` which holds an ident
91
91
92
92
## Zero mass matrix
93
93
94
-
Similar ot the example above, a helper class `daecpp::MassMatrixZero` can be used to construct a zero matrix (a matrix filled by zeros only). This matrix can be useful to solve pure algebraic systems without time derivatives.
94
+
Similar ot the example above, a helper class `daecpp::MassMatrixZero` can be used to construct a zero matrix (a matrix filled with zeros only). This matrix can be useful to solve pure algebraic systems without time derivatives.
95
95
96
96
```cpp
97
97
daecpp::MassMatrixZero zeroMatrix; // Creates empty mass matrix for algebraic system
Copy file name to clipboardExpand all lines: docs/solution-manager.md
+2-2
Original file line number
Diff line number
Diff line change
@@ -78,7 +78,7 @@ public:
78
78
{: .warning }
79
79
Vectors `x_sol` and `t_sol` should be defined before calling the `UserDefinedSolutionManager` constructor and must live until the end of the computation. They shouldn't be passed to the constructor as temporary objects and must **not** go out of scope before finishing the system solving (otherwise, this will lead to dangling references `&m_x_sol` and `&m_t_sol` in the class definition above).
80
80
81
-
Similar to the mass matrix, vector function, Jacobian matrix definitions, inhereting the `daecpp::SolutionManager` class is a good practice (it serves as a blueprint), but it is not necessary. The user is allowed to define custom Solution Managers without inhereting `daecpp::SolutionManager`.
81
+
Similar to the mass matrix, vector function and Jacobian matrix definitions, inhereting the `daecpp::SolutionManager` class is a good practice (it serves as a blueprint), but it is not necessary. The user is allowed to define custom Solution Managers without inhereting `daecpp::SolutionManager`.
82
82
83
83
----
84
84
@@ -128,7 +128,7 @@ sol.print({0, 1, 42});
128
128
129
129
# Solution class
130
130
131
-
Solution class `daecpp::Solution` is derived from `daecpp::SolutionManager` and serves as a basic solution *observer* that writes solution vectors `x` and the corresponding times `t` every time step or every specific time from `t_output` vector (that can be optionally provided in the constructor) into `daecpp::SolutionHolder` class object. The Solution class is also used in the [`daecpp::System`](solve.html) class as a default Solution Manager.
131
+
Solution class `daecpp::Solution` is derived from `daecpp::SolutionManager` and serves as a basic solution *observer* that writes solution vectors `x` and the corresponding times `t` every time step or every specific time from `t_output` vector (that can be optionally provided in the constructor) into `daecpp::SolutionHolder` class object. The Solution class is also used in the [`daecpp::System`](solve.html#system-class) class as a default Solution Manager.
Copy file name to clipboardExpand all lines: docs/solve.md
+6-6
Original file line number
Diff line number
Diff line change
@@ -18,7 +18,7 @@ In this section, we consider both approaches to start the computation.
18
18
19
19
# System class
20
20
21
-
The `System` class serves as a wrapper for lower level [`daecpp::solve(...)`](#daecppsolve-function) function calls. Contains the [Solver Options](solver-options.html) and [Solution Holder](solution-manager.html#solution-holder-class) objects.
21
+
The `System` class serves as a wrapper for lower level [`daecpp::solve(...)`](#daecppsolve-function) function calls. The class allows the user to easily set up the DAE system, solve it, and get the result. The class already contains the [Solver Options](solver-options.html) and [Solution Holder](solution-manager.html#solution-holder-class) objects.
22
22
23
23
## System class public members
24
24
@@ -51,8 +51,8 @@ Here is the summary of the `solve(...)` method:
51
51
| Parameter | Type(s) | Description |
52
52
| --------- | ------- | ----------- |
53
53
| `x0` | `const daecpp::state_vector &` | The initial condition (initial state vector) |
54
-
| `t` | `const double` **or** `const std::vector<double> &` | Integration interval `[0, t]` **or** a vector of integration times (useful if the user need the output at particular times `t`) |
55
-
| `jacobian` | User-defined [Jacobian Matrix](jacobian-matrix.html) | **(Optional)** Jacobian matrix (matrix of the RHS derivatives) |
54
+
| `t` | `const double` **or** `const std::vector<double> &` | Integration interval `[0, t]` **or** a vector of integration times (useful if the user needs the output at particular times `t`) |
55
+
| `jacobian` | User-defined [Jacobian Matrix](jacobian-matrix.html) | **(Optional)** User-defined Jacobian matrix (matrix of the RHS derivatives) |
56
56
57
57
## Examples
58
58
@@ -82,7 +82,7 @@ int status = my_system.solve({0, 1}, 1.0, MyJacobian()); // Solves the system wi
82
82
83
83
if(!status)
84
84
{
85
-
my_system.sol.print(); // Prints solution on screen
85
+
my_system.sol.print(); // Prints solution on screen if solution is successfull
86
86
}
87
87
```
88
88
@@ -110,9 +110,9 @@ Here is the summary of the `solve(...)` function:
110
110
| `rhs` | User-defined [Vector Function](vector-function.html) | Vector function (the RHS) of the DAE system |
111
111
| `jac` | User-defined [Jacobian Matrix](jacobian-matrix.html) | **(Optional)** Jacobian matrix (matrix of the RHS derivatives) |
112
112
| `x0` | `const daecpp::state_vector &` | The initial condition (initial state vector) |
113
-
| `t` | `const double` **or** `const std::vector<double> &` | Integration interval `[0, t]` **or** a vector of integration times (useful if the user need the output at particular times `t`) |
113
+
| `t` | `const double` **or** `const std::vector<double> &` | Integration interval `[0, t]` **or** a vector of integration times (useful if the user needs the output at particular times `t`) |
Copy file name to clipboardExpand all lines: docs/solver-options.md
+2-2
Original file line number
Diff line number
Diff line change
@@ -6,7 +6,7 @@ nav_order: 7
6
6
7
7
# Solver Options class
8
8
9
-
The Solver Options class, `daecpp::SolverOptions`, allows the user to redefine and fine-tune the DAE solver options. To do so, the user should create an object of the Solver Options class, update parameters of interest by accessing public variables (see [example](#example) below), and then provide the updated object to the solver function `solve`.
9
+
The Solver Options class, `daecpp::SolverOptions`, allows the user to redefine and fine-tune the DAE solver options. To do so, the user should create an object of the Solver Options class, update parameters of interest by accessing public variables (see [example](#example) below), and then provide the updated object to the solver function [`solve`](solve.html).
10
10
11
11
The table below contains a complete list of all options accessible via the Solver Options class.
12
12
@@ -25,7 +25,7 @@ The table below contains a complete list of all options accessible via the Solve
25
25
|`Newton_scheme`|`unsigned`|`1`| Non-linear solver algorithm: <br> `0` - Classic Newton method (usually the most stable but the slowest), <br> `1` (default) - Quasi-Newton method I (balanced with focus on stability, updates Jacobian and performs factorization every 2nd iteration), <br> `2` - Quasi-Newton method II (balanced with focus on speed, updates Jacobian and performs factorization every 3rd iteration), <br> `3` - Quasi-Newton method III (can be the fastest but less stable, may require tweaking the time step increase/decrease thresholds, updates Jacobian and performs factorization every 4th iteration). |
26
26
|`is_mass_matrix_static`|`bool`|`false`| If `true`, the mass matrix will be updated only once (at the beginning of computation). Setting this option as `true` slightly speeds up the computation, but the mass matrix must be static (i.e., it must be independent on time). |
27
27
|`max_Jacobian_updates`|`int`|`8`| Maximum number of the Jacobian matrix updates and factorizations per time step. If the algorithm fails to converge after `max_Jacobian_updates` Jacobian matrix updates and factorizations per time step, the Newton iterations will be considered as diverged. The solver will try to roll back and decrease the time step. |
28
-
|`max_Newton_failed_attempts`|`unsigned`|`3`| If the Newton method fails to converge 'max_Newton_failed_attempts' times in a row, the solver will try to update Jacobian matrix every single iteration next time step (for Quasi-Newton methods). |
28
+
|`max_Newton_failed_attempts`|`unsigned`|`3`| If the Newton method fails to converge `max_Newton_failed_attempts` times in a row, the solver will try to update Jacobian matrix every single iteration next time step (for Quasi-Newton methods). |
29
29
|`dt_increase_threshold_delta`|`int`|`0`| Time step amplification threshold delta. Can be negative or positive integer number. The solver increases the time step if the number of successful Newton iterations per time step is less than the threshold. <br> **Example:** If the solver increases the time step too often, decrease the time step amplification threshold by setting `dt_increase_threshold_delta` to `-1`. |
30
30
|`dt_decrease_threshold_delta`|`int`|`0`| Time step reduction threshold delta. Can be negative or positive integer number. The solver decreases the time step if the number of successful Newton iterations per time step is greater than the threshold. <br> **Example:** If the solver struggles to converge but keeps doing many Newton iterations without reducing the time step, decrease the time step reduction threshold by setting `dt_decrease_threshold_delta` to `-1` (or `-2` and less). |
31
31
|`solution_variability_control`|`bool`|`true`| Turns ON/OFF solution variability control. Solution variability control tightens up the adaptive time stepping algorithm, and it is ON by default. Switching it OFF can lead to a significant speed boost for big systems, but it can also lead to instability. |
Copy file name to clipboardExpand all lines: docs/sparse-matrix.md
+3-3
Original file line number
Diff line number
Diff line change
@@ -24,7 +24,7 @@ std::cout << M.dense(3) << '\n'; // Prints matrix M on screen assuming it's a 3x
24
24
The three-array format defines a sparse matrix using two (unsigned) integer vectors of indices `i` and `j` and one floating-point vector `A` of the corresponding non-zero elements of the matrix. Only non-zero elements of the sparse matrix should be defined and stored. Theoretically, all three arrays can be empty. This will define a zero matrix where all elements are zeros.
25
25
26
26
{: .note }
27
-
All three arrays `i`, `j`, and `A` should be the same size. This can be (and will be) checked by calling the `void check()` method of the class `daecpp::sparse_matrix`.
27
+
All three arrays `i`, `j`, and `A` should be the same size. This can be (and will be) checked by calling the [`void check()`](#void-check-const) method of the class `daecpp::sparse_matrix`.
28
28
29
29
Here is an example of a sparse matrix:
30
30
@@ -64,7 +64,7 @@ The row and column numeration starts from `0`.
64
64
For `float_type` and `int_type` definition, refer to the [Prerequisites](prerequisites.html#dae-cpp-types) section.
65
65
66
66
{: .note }
67
-
All three vectors are initially empty an NOT preallocated
67
+
All three vectors are initially empty and **not** preallocated.
68
68
69
69
## Sparse Matrix class methods
70
70
@@ -189,4 +189,4 @@ Converts matrix from `dae-cpp` three-array format to `Eigen::SparseMatrix` forma
189
189
## `std::size_t N_elements() const`
190
190
191
191
<br>
192
-
Returns the number of non-zero elements in the matrix.
192
+
Returns the number of non-zero elements in the matrix (including duplicates, if any).
Copy file name to clipboardExpand all lines: docs/vector-function.md
+1-1
Original file line number
Diff line number
Diff line change
@@ -7,7 +7,7 @@ nav_order: 4
7
7
8
8
# Vector Function class
9
9
10
-
The Vector function class defines the (nonlinear) vector function (the RHS) $$\mathbf{f}(\mathbf{x}, t)$$ of the DAE system written in the matrix-vector form:
10
+
The Vector Function class defines the (nonlinear) vector function (the RHS) $$\mathbf{f}(\mathbf{x}, t)$$ of the DAE system written in the matrix-vector form:
0 commit comments