Hello all, I am an Applied Mathematics PhD student at CU Boulder, and I want to implement fast matrix free finite difference methods using Julia. This appeared at first glance to be quite simple task, but when I compare the method to simply forming a sparse matrix, my method is slower. I believe t…

A couple of things going on here, but you need to carefully go line by line and eliminate all allocations in your core differentiation functions. Take a look at allocation tracking to figure out your worst culprits. You should know how much space you’ll need for each of your temporaries and output a…

[image] stillyslalom: Your computational kernel should specialize on the dimensionality and order of your operators, and that specialization will be made easier if you use static arrays for your operators, and pass those in as arguments as well. Alternatively, see this post for an example of …

[image] nicrummel: left = CartesianIndices(sz[1:dim-1]) right = CartesianIndices(sz[dim+1:N]) This pattern and the subsequent indexing will be type-unstable because the 1:dim-1 cannot be resolved at compile time. One option would be to change the dim::Int input to ::Val{dim} to lift i…

[image] stevengj: ghost cells What do you mean by that?

[image] RoyiAvital: What do you mean by [“ghost cells”]? See e.g. Arrays with periodic boundaries - #4 by stevengj

Follow Up Thanks for all the help everyone. I made some update to my code that made it work much better! My main takeaways are: Break functions into and outer general function and an inner specialized function. This allows the compiler to optimize the inner function based on the size/type of my …

[image] nicrummel: leftIx = CartesianIndices(ntuple(i->sz[i], Val(dim-1))) Note that this Val(dim-1) is type-unstable because dim is a runtime value. To make to make it compile-time, change the signature to function D(..., ::Val{dim}, ...) where {T<:Real,N,dim} and call it with Val(dim) ins…

Matrix Free Differentiation

General Usage

stevengj December 5, 2023, 8:42pm 6

See e.g. Arrays with periodic boundaries - #4 by stevengj

Topic		Replies	Views
Julia vs Matlab - building a Jacobian matrix New to Julia	30	4390	March 28, 2019
Efficient finite difference operators Performance diffeq	19	2685	August 8, 2018
Matrix construction operations - increasing speed New to Julia performance	13	1163	July 5, 2018
Making efficient some small subroutines Performance	22	1246	May 8, 2020
Help with this code (deallocation/pointers on sparse matrices) Performance array , sparse , matrices , sparsearrays	31	907	February 3, 2024

Matrix Free Differentiation

Related topics