diff --git a/docs/src/getting_started/find_root.md b/docs/src/getting_started/find_root.md index 343f7daf41d..ee8e23a145a 100644 --- a/docs/src/getting_started/find_root.md +++ b/docs/src/getting_started/find_root.md @@ -3,17 +3,25 @@ A nonlinear system $$f(u) = 0$$ is specified by defining a function `f(u,p)`, where `p` are the parameters of the system. Many problems can be written in such as way that solving a nonlinear rootfinding problem gives the solution. -For example: +For example: -* Do you want to know ``u`` such that ``4^u + 6^u = 7^u``? Then solve +* Do you want to know ``u`` such that ``4^u + 6^u = 7^u``? Then solve ``f(u) = 4^u + 6^u - 7^u = 0`` for `u`! * If you have an ODE ``u' = f(u)``, what is the point where the solution will be completely still, i.e. `u' = 0`? -* -All of these problems are solved by using a numerical rootfinder. Let's solve +All of these problems are solved by using a numerical rootfinder. Let's solve our first rootfind problem! +## Required Dependencies + +The following parts of the SciML Ecosystem will be used in this tutorial: + +| Module | Description | +| ----------- | ----------- | +| [ModelingToolkit.jl](https://docs.sciml.ai/ModelingToolkit/stable/) | The symbolic modeling environment | +| [NonlinearSolve.jl](https://docs.sciml.ai/NonlinearSolve/stable/) | The symbolic modeling environment | + ## Problem Setup For example, the following solves the vector equation: @@ -51,43 +59,121 @@ sol = solve(prob,NewtonRaphson()) @show sol.u, prob.f(sol.u,prob.p) ``` -## Step 1: Import the Packages +## Step by Step Solution + +### Step 1: Import the Packages + +To do this tutorial we will need a few components: + +* [ModelingToolkit.jl, our modeling environment](https://docs.sciml.ai/ModelingToolkit/stable/) +* [NonlinearSolve.jl, the nonlinear system solvers](https://docs.sciml.ai/NonlinearSolve/stable/) + +To start, let's add these packages [as demonstrated in the installation tutorial](@ref installation): + +```julia +]add ModelingToolkit NonlinearSolve +``` + +Now we're ready. Let's load in these packages: ```@example first_rootfind # Import the packages using ModelingToolkit, NonlinearSolve ``` -## Step 2: Define the Nonlinear System +### Step 2: Define the Nonlinear System + +Now let's define our nonlinear system. We use the `ModelingToolkit.@variabes` statement to +declare our 3 state variables: ```@example first_rootfind # Define the nonlinear system @variables x=1.0 y=0.0 z=0.0 +``` + +Notice that we are using the form `state = initial condition`. This is a nice shorthand +for coupling an initial condition to our states. We now must similarly define our parameters, +which we can associate default values via the form `parameter = default value`. This looks +like: + +```@example first_rootfind @parameters σ=10.0 ρ=26.0 β=8/3 +``` + +Now we create an array of equations to define our nonlinear system that must be satisfied. +This looks as follows: + +!!! note + + Note that in ModelingToolkit and Symbolics, `~` is used for equation equality. This is + separate from `=` which is the "assignment operator" in the Julia programming language. + For example, `x = x + 1` is a valid assignment in a programming language, and it is + invalid for that to represent "equality", which is the reason why a separate operator + is used! +```@example first_rootfind eqs = [0 ~ σ*(y-x), 0 ~ x*(ρ-z)-y, 0 ~ x*y - β*z] +``` + +Finally we bring these pieces together, the equation along with its states and parameters, +define our `NonlinearSystem`: + +```@example first_rootfind @named ns = NonlinearSystem(eqs, [x,y,z], [σ,ρ,β]) ``` -## Step 3: Convert the Symbolic Problem to a Numerical Problem +### Step 3: Convert the Symbolic Problem to a Numerical Problem + +Now that we have simplified our system, let's turn it into a numerical problem to +approximate. This is done with the `NonlinearProblem` constructor, that transforms it from +a symbolic `ModelingToolkit` representation to a numerical `NonlinearSolve` +representation. We need to tell it the numerical details for whether to override any of the +default values for the initial conditions and parameters. + +In this case, we will use the default values for all our variables, so we will pass a +blank override `[]`. This looks like: ```@example first_rootfind # Convert the symbolic system into a numerical system prob = NonlinearProblem(ns,[]) ``` -## Step 4: Solve the Numerical Problem +If for example we did want to change the initial condition of `x` +to `2.0` and the parameter `σ` to `4.0`, we would do `[x => 2.0, σ => 4.0]`. This looks +like: + +```@example first_rootfind +prob2 = NonlinearProblem(ns,[x => 2.0, σ => 4.0]) +``` + +### Step 4: Solve the Numerical Problem + +Now we solve the nonlinear system. For this we choose a solver from the +[NonlinearSolve.jl's solver options](https://docs.sciml.ai/NonlinearSolve/dev/solvers/NonlinearSystemSolvers/) +We will choose `NewtonRaphson` as follows: ```@example first_rootfind # Solve the numerical problem sol = solve(prob,NewtonRaphson()) ``` -## Step 5: Analyze the Solution +### Step 5: Analyze the Solution + +Now let's check out the solution. First of all, what kind of thing is the `sol`? We can +see that by asking for its type: + +```@example first_rootfind +typeof(sol) +``` + +From this we can see that it is an `NonlinearSolution`. We can see the documentation for +how to use the `NonlinearSolution` by checking the +[NonlinearSolve.jl solution type page](https://docs.sciml.ai/NonlinearSolve/dev/basics/NonlinearSolution/) Thus the solution is stored as `.u`. What is the solution to our +nonlinear system and what is the final residual value? We can check it as follows: ```@example first_rootfind # Analyze the solution -@show sol.u, prob.f(sol.u,prob.p) -``` \ No newline at end of file +@show sol.u, sol.resid +``` diff --git a/docs/src/getting_started/first_optimization.md b/docs/src/getting_started/first_optimization.md index 63afc1b9541..46de973bf94 100644 --- a/docs/src/getting_started/first_optimization.md +++ b/docs/src/getting_started/first_optimization.md @@ -1,17 +1,25 @@ # [Solve your first optimization problem](@id first_opt) **Numerical optimization** is the process of finding some numerical values that -minimize some equation. +minimize some equation. -* How much fuel should you put into an airplane to have the minimum weight that +* How much fuel should you put into an airplane to have the minimum weight that can go to its destination? * What parameters should I choose for my simulation so that it minimizes the distance of its predictions from my experimental data? -* -All of these are examples of problems solved by numerical optimization. +All of these are examples of problems solved by numerical optimization. Let's solve our first optimization problem! +## Required Dependencies + +The following parts of the SciML Ecosystem will be used in this tutorial: + +| Module | Description | +| ----------- | ----------- | +| [Optimization.jl](https://docs.sciml.ai/Optimization/stable/) | The numerical optimization package | +| [OptimizationNLopt.jl](https://docs.sciml.ai/Optimization/stable/optimization_packages/nlopt/) | The NLopt optimizers we will use | + ## Problem Setup First, what are we solving? Let's take a look at the Rosenbrock equation: @@ -20,10 +28,10 @@ First, what are we solving? Let's take a look at the Rosenbrock equation: L(u,p) = (p_1 - u_1)^2 + p_2 * (u_2 - u_1)^2 ``` -What we want to do is find the values of ``u_1`` and ``u_2`` such that ``L`` -achieves its minimum value possible. We will do this under a few constraints: -we want to find this optima within some bounded domain, i.e. ``u_i \in [-1,1]``. -This should be done with the parameter values ``p_1 = 1.0`` and `p_2 = 100.0``. +What we want to do is find the values of ``u_1`` and ``u_2`` such that ``L`` +achieves its minimum value possible. We will do this under a few constraints: +we want to find this optima within some bounded domain, i.e. ``u_i \in [-1,1]``. +This should be done with the parameter values ``p_1 = 1.0`` and `p_2 = 100.0``. What should ``u = [u_1,u_2]`` be to achieve this goal? Let's dive in! !!! note @@ -34,7 +42,7 @@ What should ``u = [u_1,u_2]`` be to achieve this goal? Let's dive in! ## Copy-Pastable Code ```@example -# Import the package +# Import the package using Optimization, OptimizationNLopt # Define the problem to optimize @@ -50,35 +58,103 @@ sol = solve(prob,NLopt.LD_LBFGS()) @show sol.u, L(sol.u,p) ``` -## Step 1: Import the packages +## Step by Step Solution + +### Step 1: Import the packages + +To do this tutorial we will need a few components: + +* [Optimization.jl](https://docs.sciml.ai/Optimization/stable/), the optimization interface. +* [OptimizationNLopt.jl](https://docs.sciml.ai/Optimization/stable/optimization_packages/nlopt/), the optimizers we will use. + +Note that Optimization.jl is an interface for optimizers, and thus we always have to choose +which optimizer we want to use. Here we choose to demonstrate `OptimizationNLopt` because +of its efficiency and versitility. But there are many other possible choices. Check out +the +[solver compatability chart](https://docs.sciml.ai/Optimization/stable/#Overview-of-the-Optimizers) +for a quick overview of what optimizer packages offer. + +To start, let's add these packages [as demonstrated in the installation tutorial](@ref installation): + +```julia +]add Optimization OptimizationNLopt +``` + +Now we're ready. Let's load in these packages: -```@example first_opt +```@example first_opt using Optimization, OptimizationNLopt ``` -## Step 2: Define the Optimization Problem +### Step 2: Define the Optimization Problem + +Now let's define our problem to optimize. We start by defining our loss function. In +Optimization.jl's `OptimizationProblem` interface, the states are given by an array +`u`. Thus we can designate `u[1]` to be `u_1` and `u[2]` to be `u_2`, similarly with our +parameters, and write out the loss function on a vector-defined state as follows: ```@example first_opt # Define the problem to optimize L(u,p) = (p[1] - u[1])^2 + p[2] * (u[2] - u[1]^2)^2 ``` +Now we need to define our `OptimizationProblem`. If you need help remembering how to define +the `OptimizationProblem`, you can always refer to the +[Optimization.jl problem definition page](https://docs.sciml.ai/Optimization/stable/API/optimization_problem/) + +Thus what we need to define is an initial condition `u0` and our parameter vector `p`. +We will make our initial condition have both values as zero, which is done by the Julia +shorthand `zeros(2)` that creates a vector `[0.0,0.0]`. We manually define the parameter +vector `p` to input our values. Then we set the lower bound and upper bound for the +optimization as follows: + ```@example first_opt u0 = zeros(2) p = [1.0,100.0] prob = OptimizationProblem(L, u0, p, lb = [-1.0,-1.0], ub = [1.0,1.0]) ``` -## Step 3: Solve the Optimization Problem +#### Note about defining uniform bounds + +Note that we can simplify the code a bit for the lower and upper bound definition by +using the Julia Base command `ones`, which returns a vector where each value is a one. +Thus for example, `ones(2)` is equivalent to `[1.0,1.0]`. Therefore `-1 * ones(2)` is +equivalent to `[-1.0, 1.0]`, meaning we could have written our problem as follows: + +```@example first_opt +prob = OptimizationProblem(L, u0, p, lb = -1 * ones(2), ub = ones(2)) +``` + +### Step 3: Solve the Optimization Problem + +Now we solve the `OptimizationProblem` that we have defined. This is done by passing +our `OptimizationProblem` `prob` along with a chosen solver to the `solve` command. At +the beginning we explained that we will use the `OptimizationNLopt` set of solvers, which +are +[documented in the OptimizationNLopt page](https://docs.sciml.ai/Optimization/stable/optimization_packages/nlopt/). +From here we are choosing the `NLopt.LD_LBFGS()` for its mixture of robustness and +performance. To perform this solve, we do the following: ```@example first_opt # Solve the optimization problem sol = solve(prob,NLopt.LD_LBFGS()) ``` -## Step 4: Analyze the Solution +### Step 4: Analyze the Solution + +Now let's check out the solution. First of all, what kind of thing is the `sol`? We can +see that by asking for its type: + +```@example first_opt +typeof(sol) +``` + +From this we can see that it is an `OptimizationSolution`. We can see the documentation for +how to use the `OptimizationSolution` by checking the +[Optimization.jl solution type page](https://docs.sciml.ai/Optimization/dev/API/optimization_solution/). Thus the solution is stored as `.u`. What is the solution to our +optimization and what is the final loss value? We can check it as follows: ```@example first_opt # Analyze the solution @show sol.u, L(sol.u,p) -``` \ No newline at end of file +``` diff --git a/docs/src/getting_started/first_simulation.md b/docs/src/getting_started/first_simulation.md index 16b35e3a50c..c7c3d16afeb 100644 --- a/docs/src/getting_started/first_simulation.md +++ b/docs/src/getting_started/first_simulation.md @@ -7,12 +7,23 @@ In this tutorial we will build and run our first simulation with SciML! This tutorial assumes that you have already installed Julia on your system. If you have not done so already, please [follow the installation tutorial first](@ref installation). -To build our simulation, we will use the [ModelingToolkit](https://mtk.sciml.ai/dev/) +To build our simulation, we will use the +[ModelingToolkit](https://docs.sciml.ai/ModelingToolkit/stable/) system for modeling and simulation. ModelingToolkit is a bit higher level than directly defining code for a differential equation system: it's a symbolic system that will automatically simplify our models, optimize our code, and generate compelling visualizations. Sounds neat? Let's dig in. +## Required Dependencies + +The following parts of the SciML Ecosystem will be used in this tutorial: + +| Module | Description | +| ----------- | ----------- | +| [ModelingToolkit.jl](https://docs.sciml.ai/ModelingToolkit/stable/) | The symbolic modeling environment | +| [DifferentialEquations.jl](https://docs.sciml.ai/DiffEqDocs/stable/) | The differential equation solvers | +| [Plots.jl](https://docs.juliaplots.org/stable/) | The plotting and visualization package | + ## Our Problem: Simulate the Lotka-Volterra Predator-Prey Dynamics The dynamics of our system are given by the @@ -82,8 +93,8 @@ plot(p1,p2,layout=(2,1)) To do this tutorial we will need a few components: -* [ModelingToolkit.jl, our modeling environment](https://mtk.sciml.ai/dev/) -* [DifferentialEquations.jl, the differential equation solvers](https://diffeq.sciml.ai/stable/) +* [ModelingToolkit.jl, our modeling environment](https://docs.sciml.ai/ModelingToolkit/stable/) +* [DifferentialEquations.jl, the differential equation solvers](https://docs.sciml.ai/DiffEqDocs/stable/) * [Plots.jl, our visualization tool](https://docs.juliaplots.org/stable/) To start, let's add these packages [as demonstrated in the installation tutorial](@ref installation):