I am new to Julia and I try to understand why a solve using the backslash operator allocates memory although the RHS vector is allocated, the solution vector is allocated and the factorisation is precomputed (here Cholesky factorisation). Is there a way to avoid such behaviour?
Below is a minimum example.
Cheers and thanks for the hints!
using LinearAlgebra, SparseArrays function Successive_solves!( x, Kc, b, nsolves ) for i=1:nsolves x .= Kc\b end return end function ThreePointStencil( N, dx ) # Assemble matrix (not the part of interest) cC = 2.0/dx^2 * ones(Float64, N); cC = 3.0/dx^2; cC[end] = 3.0/dx^2 cW = -1.0/dx^2 * ones(Float64, N); cW = 0.0 cE = -1.0/dx^2 * ones(Float64, N); cE[end] = 0.0 iC = 1:N iW = ones(Int64, N); iW[2:end-0] .= iC[1:end-1] iE = ones(Int64, N); iE[1:end-1] .= iC[2:end-0] I = [iC[:]; iC[:]; iC[:]] J = [iC[:]; iW[:]; iE[:]] V = [cC[:]; cW[:]; cE[:]] K = sparse( I, J, V ) droptol!(K, 1e-9) return K end # Input N = 1000000 dx = 1 # Sparse matrix K = ThreePointStencil( N, dx ) # Allocate arrays b = rand(N) # Right-hand-side array x = zeros(Float64, N) # Solution array # Factor Kc = cholesky(K) # One solve @time x .= Kc\b # why does this allocate memory although the factorisation done and the vectors (x, b) are allocated? # Successive solves nsolves = 10 @time Successive_solves!( x, Kc, b, nsolves ) # then, memory allocation depends on the number of solves...