Skip to content

Format #276

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -47,10 +47,10 @@ Let `M`, `D`, `F` be matrix-based, diagonal-matrix-based, and function-based
```@example operator_algebra
using SciMLOperators, LinearAlgebra
N = 4
function f(v, u, p, t)
function f(v, u, p, t)
u .* v
end
function f(w, v, u, p, t)
function f(w, v, u, p, t)
w .= u .* v
end

Expand Down
4 changes: 1 addition & 3 deletions docs/pages.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,4 @@ pages = [
"FFT Tutorial" => "tutorials/fftw.md"
],
"interface.md",
"Premade Operators" => "premade_operators.md",

]
"Premade Operators" => "premade_operators.md"]
4 changes: 2 additions & 2 deletions docs/src/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@ act like “normal” functions for equation solvers. For example, if `A(v,u,p,t
has the same operation as `update_coefficients(A, u, p, t); A * v`, then `A`
can be used in any place where a differential equation definition
`(u,p,t) -> A(u, u, p, t)` is used without requiring the user or solver to do any extra
work.
work.

Another example is state-dependent mass matrices. `M(u,p,t)*u' = f(u,p,t)`.
When solving such an equation, the solver must understand how to "update M"
Expand All @@ -66,7 +66,7 @@ an extended operator interface with all of these properties, hence the
`AbstractSciMLOperator` interface.

!!! warn

This means that LinearMaps.jl is fundamentally lacking and is incompatible
with many of the tools in the SciML ecosystem, except for the specific cases
where the matrix-free operator is a constant!
Expand Down
69 changes: 36 additions & 33 deletions docs/src/tutorials/getting_started.md
Original file line number Diff line number Diff line change
Expand Up @@ -11,16 +11,16 @@ and updating operators.
Before we get into the deeper operators, let's show the simplest SciMLOperator:
`MatrixOperator`. `MatrixOperator` just turns a matrix into an `AbstractSciMLOperator`,
so it's not really a matrix-free operator but it's a starting point that is good for
understanding the interface and testing. To create a `MatrixOperator`, simply call the
understanding the interface and testing. To create a `MatrixOperator`, simply call the
constructor on a matrix:

```@example getting_started
using SciMLOperators, LinearAlgebra
A = [-2.0 1 0 0 0
1 -2 1 0 0
0 1 -2 1 0
0 0 1 -2 1
0 0 0 1 -2]
A = [-2.0 1 0 0 0
1 -2 1 0 0
0 1 -2 1 0
0 0 1 -2 1
0 0 0 1 -2]

opA = MatrixOperator(A)
```
Expand All @@ -29,7 +29,7 @@ The operators can do [operations as defined in the operator interface](@ref oper
matrix multiplication as the core action:

```@example getting_started
v = [3.0,2.0,1.0,2.0,3.0]
v = [3.0, 2.0, 1.0, 2.0, 3.0]
opA*v
```

Expand All @@ -43,7 +43,8 @@ mul!(w, opA, v)
```

```@example getting_started
α = 1.0; β = 1.0
α = 1.0;
β = 1.0
mul!(w, opA, v, α, β) # α*opA*v + β*w
```

Expand All @@ -64,18 +65,20 @@ For example, let's make the operator `A .* u + dt*I` where `dt` is a parameter
and `u` is a state vector:

```@example getting_started
A = [-2.0 1 0 0 0
1 -2 1 0 0
0 1 -2 1 0
0 0 1 -2 1
0 0 0 1 -2]
A = [-2.0 1 0 0 0
1 -2 1 0 0
0 1 -2 1 0
0 0 1 -2 1
0 0 0 1 -2]

function update_function!(B, u, p, t)
dt = p
B .= A .* u + dt*I
end

u = Array(1:1.0:5); p = 0.1; t = 0.0
u = Array(1:1.0:5);
p = 0.1;
t = 0.0
opB = MatrixOperator(copy(A); update_func! = update_function!)
```

Expand Down Expand Up @@ -123,18 +126,18 @@ With `FunctionOperator`, we directly define the operator application function `o
which means `w = opA(u,p,t)*v`. For example we can do the following:

```@example getting_started
function Afunc!(w,v,u,p,t)
function Afunc!(w, v, u, p, t)
w[1] = -2v[1] + v[2]
for i in 2:4
w[i] = v[i-1] - 2v[i] + v[i+1]
w[i] = v[i - 1] - 2v[i] + v[i + 1]
end
w[5] = v[4] - 2v[5]
nothing
end

function Afunc!(v,u,p,t)
function Afunc!(v, u, p, t)
w = zeros(5)
Afunc!(w,v,u,p,t)
Afunc!(w, v, u, p, t)
w
end

Expand All @@ -148,29 +151,29 @@ mfopA*v - opA*v
```

```@example getting_started
mfopA(v,u,p,t) - opA(v,u,p,t)
mfopA(v, u, p, t) - opA(v, u, p, t)
```

We can also create the state-dependent operator as well:

```@example getting_started
function Bfunc!(w,v,u,p,t)
function Bfunc!(w, v, u, p, t)
dt = p
w[1] = -(2*u[1]-dt)*v[1] + v[2]*u[1]
for i in 2:4
w[i] = v[i-1]*u[i] - (2*u[i]-dt)*v[i] + v[i+1]*u[i]
w[i] = v[i - 1]*u[i] - (2*u[i]-dt)*v[i] + v[i + 1]*u[i]
end
w[5] = v[4]*u[5] - (2*u[5]-dt)*v[5]
nothing
end

function Bfunc!(v,u,p,t)
function Bfunc!(v, u, p, t)
w = zeros(5)
Bfunc!(w,v,u,p,t)
Bfunc!(w, v, u, p, t)
w
end

mfopB = FunctionOperator(Bfunc!, zeros(5), zeros(5); u, p, t, isconstant=false)
mfopB = FunctionOperator(Bfunc!, zeros(5), zeros(5); u, p, t, isconstant = false)
```

```@example getting_started
Expand All @@ -187,19 +190,19 @@ operator for `A .* u` (since right now there is not a built in operator for vect
but that would be a fantastic thing to add!):

```@example getting_started
function Cfunc!(w,v,u,p,t)
function Cfunc!(w, v, u, p, t)
w[1] = -2v[1] + v[2]
for i in 2:4
w[i] = v[i-1] - 2v[i] + v[i+1]
w[i] = v[i - 1] - 2v[i] + v[i + 1]
end
w[5] = v[4] - 2v[5]
w .= w .* u
nothing
end

function Cfunc!(v,u,p,t)
function Cfunc!(v, u, p, t)
w = zeros(5)
Cfunc!(w,v,u,p,t)
Cfunc!(w, v, u, p, t)
w
end

Expand Down Expand Up @@ -227,11 +230,11 @@ adjoints, inverses, and more. For more information, see the [operator algebras t
Great! You now know how to be state/parameter/time-dependent operators and make them matrix-free, along with
doing algebras on operators. What's next?

* Interested in more examples of building operators? See the example of [making a fast fourier transform linear operator](@ref fft)
* Interested in more operators ready to go? See the [Premade Operators page](@ref premade_operators) for all of the operators included with SciMLOperators. Note that there are also downstream packages that make new operators.
* Want to make your own SciMLOperator? See the [AbstractSciMLOperator interface page](@ref operator_interface) which describes the full interface.
- Interested in more examples of building operators? See the example of [making a fast fourier transform linear operator](@ref fft)
- Interested in more operators ready to go? See the [Premade Operators page](@ref premade_operators) for all of the operators included with SciMLOperators. Note that there are also downstream packages that make new operators.
- Want to make your own SciMLOperator? See the [AbstractSciMLOperator interface page](@ref operator_interface) which describes the full interface.

How do you use SciMLOperators? Check out the following downstream pages:

* [Using SciMLOperators in LinearSolve.jl for matrix-free Krylov methods](https://docs.sciml.ai/LinearSolve/stable/tutorials/linear/)
* [Using SciMLOperators in OrdinaryDiffEq.jl for semi-linear ODE solvers](https://docs.sciml.ai/DiffEqDocs/stable/solvers/nonautonomous_linear_ode/)
- [Using SciMLOperators in LinearSolve.jl for matrix-free Krylov methods](https://docs.sciml.ai/LinearSolve/stable/tutorials/linear/)
- [Using SciMLOperators in OrdinaryDiffEq.jl for semi-linear ODE solvers](https://docs.sciml.ai/DiffEqDocs/stable/solvers/nonautonomous_linear_ode/)
6 changes: 3 additions & 3 deletions docs/src/tutorials/operator_algebras.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,10 +7,10 @@ in order to build more complex objects and using their operations.
```@example operator_algebra
using SciMLOperators, LinearAlgebra
N = 4
function f(v, u, p, t)
function f(v, u, p, t)
u .* v
end
function f(w, v, u, p, t)
function f(w, v, u, p, t)
w .= u .* v
end

Expand Down Expand Up @@ -56,4 +56,4 @@ L4 = cache_operator(L4, u)
# allocation-free evaluation
L2(w, v, u, p, t) # == mul!(w, L2, v)
L4(w, v, u, p, t, α, β) # == mul!(w, L4, v, α, β)
```
```
13 changes: 6 additions & 7 deletions src/SciMLOperators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -35,20 +35,20 @@ the following type of equation:
w = L(u,p,t)[v]
```

where `L[v]` is the operator application of ``L`` on the vector ``v``.
where `L[v]` is the operator application of ``L`` on the vector ``v``.

## Interface

An `AbstractSciMLOperator` can be called like a function in the following ways:

- `L(v, u, p, t)` - Out-of-place application where `v` is the action vector and `u` is the update vector
- `L(w, v, u, p, t)` - In-place application where `w` is the destination, `v` is the action vector, and `u` is the update vector
- `L(w, v, u, p, t, α, β)` - In-place application with scaling: `w = α*(L*v) + β*w`
- `L(v, u, p, t)` - Out-of-place application where `v` is the action vector and `u` is the update vector
- `L(w, v, u, p, t)` - In-place application where `w` is the destination, `v` is the action vector, and `u` is the update vector
- `L(w, v, u, p, t, α, β)` - In-place application with scaling: `w = α*(L*v) + β*w`

Operator state can be updated separately from application:

- `update_coefficients!(L, u, p, t)` for in-place operator update
- `L = update_coefficients(L, u, p, t)` for out-of-place operator update
- `update_coefficients!(L, u, p, t)` for in-place operator update
- `L = update_coefficients(L, u, p, t)` for out-of-place operator update

SciMLOperators also overloads `Base.*`, `LinearAlgebra.mul!`,
`LinearAlgebra.ldiv!` for operator evaluation without updating operator state.
Expand Down Expand Up @@ -183,7 +183,6 @@ M = MatrixOperator(zero(N, N); update_func = mat_update_func,
M(v, u, p, t) == zeros(N) # true
M(v, u, p, t; scale = 1.0) != zero(N)
```

"""
abstract type AbstractSciMLOperator{T} end

Expand Down
Loading
Loading