How to automatically produce this kind of matrix in Julia?

General Usage
1

Hello!

I am unsure of how to make a proper MWE, so I will explain what I would like to achieve, which is the automatic construction of this matrix:

image
image730×204 18.6 KB
[Source: https://dual.sphysics.org/7thworkshop/workshop/program/7DSW_mDBC.pdf]

My code needs to adapt, in the sense that if it is a 2D case, then the last row and column from the matrix A is neglected. Currently what I do is using StaticArrays, SMatrix functionality, and construct column per column the matrix as shown below.

I have had to use an if statement, I wonder if I could do this without an if statement in a smart way?

Kind regards

EDIT: Code below is just to explain how I am handling it right now

                # Filling the Aᵧ matrix is done in column-major order
                xⱼᵢ = -xᵢⱼ

                if Dimensions == 2
                    Aᵧ[i] += SMatrix{DimensionsPlus, DimensionsPlus, FloatType, DimensionsPlus*DimensionsPlus}(
                        # First row
                        VⱼWᵢⱼ, (Vⱼ * ∇ᵢWᵢⱼ)...,
                        # Second row
                        xⱼᵢ[1] * VⱼWᵢⱼ, (xⱼᵢ[1] * Vⱼ * ∇ᵢWᵢⱼ)...,
                        # Third row
                        xⱼᵢ[2] * VⱼWᵢⱼ, (xⱼᵢ[2] * Vⱼ * ∇ᵢWᵢⱼ)...,
                    )
                elseif Dimensions == 3
                    Aᵧ[i] += SMatrix{DimensionsPlus, DimensionsPlus, FloatType, DimensionsPlus*DimensionsPlus}(
                        # First row
                        VⱼWᵢⱼ, (Vⱼ * ∇ᵢWᵢⱼ)...,
                        # Second row
                        xⱼᵢ[1] * VⱼWᵢⱼ, (xⱼᵢ[1] * Vⱼ * ∇ᵢWᵢⱼ)...,
                        # Third row
                        xⱼᵢ[2] * VⱼWᵢⱼ, (xⱼᵢ[2] * Vⱼ * ∇ᵢWᵢⱼ)...,
                        # Fourth row
                        xⱼᵢ[3] * VⱼWᵢⱼ, (xⱼᵢ[3] * Vⱼ * ∇ᵢWᵢⱼ)...,
                    )
2

Isn’t this matrix the product A'B of two matrices A and B, where the different values j go on the rows of the matrices and A contains rows W_{gj}, \partial_xW_{gj}, \partial_yW_{gj},\partial_zW_{gj}, so it is the hcat of a 1 column matrix (W_{gj})_j and a derivative matrix (\partial_i W_{gj})_{j, i \in \{x, y, z\}}? This then works nicely as soon as you can get the derivative matrix automatically? Similarly, the matrix B has rows V_j, (x_j-x_g)V_j,\ldots. Would that help?

1 Like
3

Thanks a lot!

I managed to do it like this:

                # Filling the Aᵧ matrix is done in column-major order
                xⱼᵢ = -xᵢⱼ
                first_column = [VⱼWᵢⱼ; Vⱼ * ∇ᵢWᵢⱼ]
                Aᵧ[i] += SMatrix{DimensionsPlus, DimensionsPlus, FloatType, DimensionsPlus*DimensionsPlus}(
                    first_column...,
                    ((xⱼᵢ * first_column')')...
                )

Kind regards

2 Likes