Fix tensor application to slides and fix README

Author: ChrisRackauckas

README.md

@@ -41,52 +41,61 @@ julia> Pkg.add("SciMLOperators")
 Let `M`, `D`, `F` be matrix-based, diagonal-matrix-based, and function-based
 `SciMLOperators` respectively.
 
-```julia
+Let `M`, `D`, `F` be matrix-based, diagonal-matrix-based, and function-based
+`SciMLOperators` respectively.
+
+```@example operator_algebra
+using SciMLOperators, LinearAlgebra
 N = 4
-f(u, p, t) = u .* u
-f(v, u, p, t) = v .= u .* u
+function f(v, u, p, t) 
+    u .* v
+end
+function f(w, v, u, p, t) 
+    w .= u .* v
+end
+
+u = rand(4)
+p = nothing # parameter struct
+t = 0.0     # time
 
 M = MatrixOperator(rand(N, N))
 D = DiagonalOperator(rand(N))
-F = FunctionOperator(f, zeros(N), zeros(N))
+F = FunctionOperator(f, zeros(N), zeros(N); u, p, t)
 ```
 
 Then, the following codes just work.
 
-```julia
+```@example operator_algebra
 L1 = 2M + 3F + LinearAlgebra.I + rand(N, N)
 L2 = D * F * M'
 L3 = kron(M, D, F)
-L4 = M \ D
+L4 = lu(M) \ D
 L5 = [M; D]' * [M F; F D] * [F; D]
 ```
 
 Each `L#` can be applied to `AbstractVector`s of appropriate sizes:
 
-```julia
-p = nothing # parameter struct
-t = 0.0     # time
-
-u = rand(N)
-v = L1(u, p, t) # == L1 * u
+```@example operator_algebra
+v = rand(N)
+w = L1(v, u, p, t) # == L1 * v
 
-u_kron = rand(N^3)
-v_kron = L3(u_kron, p, t) # == L3 * u_kron
+v_kron = rand(N^3)
+w_kron = L3(v_kron, u, p, t) # == L3 * v_kron
 ```
 
-For mutating operator evaluations, call `cache_operator` to generate
-in-place cache so the operation is nonallocating.
+For mutating operator evaluations, call `cache_operator` to generate an
+in-place cache, so the operation is nonallocating.
 
-```julia
+```@example operator_algebra
 α, β = rand(2)
 
 # allocate cache
 L2 = cache_operator(L2, u)
 L4 = cache_operator(L4, u)
 
 # allocation-free evaluation
-L2(v, u, p, t) # == mul!(v, L2, u)
-L4(v, u, p, t, α, β) # == mul!(v, L4, u, α, β)
+L2(w, v, u, p, t) # == mul!(w, L2, v)
+L4(w, v, u, p, t, α, β) # == mul!(w, L4, v, α, β)
 ```
 
 ## Roadmap

docs/src/tutorials/operator_algebras.md

@@ -1,7 +1,8 @@
 ## [Demonstration of Operator Algebras and Kron](@id operator_algebras)
 
 Let `M`, `D`, `F` be matrix-based, diagonal-matrix-based, and function-based
-`SciMLOperators` respectively.
+`SciMLOperators` respectively. Here are some examples of composing operators
+in order to build more complex objects and using their operations.
 
 ```@example operator_algebra
 using SciMLOperators, LinearAlgebra
@@ -55,15 +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, α, β)
-```
-
-The calling signature `L(v, u, p, t)`, for out-of-place evaluations, is
-equivalent to `L * v`, and the in-place evaluation `L(w, v, u, p, t, args...)`
-is equivalent to `LinearAlgebra.mul!(w, L, v, args...)`, where the arguments
-`u, p, t` are passed to `L` to update its state. More details are provided
-in the operator update section below.
-
-The `(v, u, p, t)` calling signature is standardized over the `SciML`
-ecosystem and is flexible enough to support use cases such as time-evolution
-in ODEs, as well as sensitivity computation with respect to the parameter
-object `p`.
+```

src/tensor.jl

@@ -172,7 +172,7 @@ function Base.:*(L::TensorProductOperator, v::AbstractVecOrMat)
     k = size(v, 2)
 
     U = reshape(v, (ni, no * k))
-    C = inner * U
+    C = stack([inner * _v for _v in eachcol(U)])
 
     V = outer_mul(L, v, C)