General concepts » Matrix storage and finite-element component ordering

Matrix storage, tensor storage and finite-element coefficient numbering describe different representations. Storage order determines the position of an entry in an object's memory. Component numbering determines the association between a coefficient and a basis function. Mathematical indexing remains independent of both representations.

Matrix storage

Matrix and FixedSizeMatrix use Eigen's default column-major storage. For a dense matrix with $m$ rows, entry $A_{ij}$ has offset $i+mj$ , with zero-based indices. Rows vary fastest within each column. SparseMatrix explicitly uses column-major compressed sparse storage; its entries are not a contiguous dense array.

SpatialMatrix also uses column-major storage, in a fixed three-by-three backing matrix. Its runtime dimensions select an active upper-left block. Consequently, its backing-storage offset is $i+3j$ , including when the active row count is smaller than three. The active entries of a rectangular spatial matrix must not be treated as a packed array with leading dimension equal to its runtime row count. Access through operator()(i,j) preserves the mathematical indexing.

Explicit row-major Eigen buffers occur in local integration kernels. These are internal storage choices and do not change the default storage of public matrix types or the numbering of finite-element coefficients.

Matrix finite-element coefficients

The matrix specializations of P0, P0g, P1 and H1 currently interleave row-major components within each scalar basis index. For an $m\times n$ range, scalar basis index $a$ and matrix component $(i,j)$ , the coefficient index is

\[ q=a(mn)+c,\qquad c=in+j. \]

The corresponding basis function is $\phi_a E_{ij}$ , where $E_{ij}$ is the matrix unit with a single nonzero entry. For a two-by-three range, the component sequence at each scalar basis index is

\[ (0,0),\ (0,1),\ (0,2),\ (1,0),\ (1,1),\ (1,2). \]

This numbering is distinct from column-major matrix storage. No Voigt scaling, symmetry assumption or component reduction is applied. A column-major component convention would instead use $c=i+mj$ ; that is not the current matrix-space coefficient convention.

Within these space families, local and global component expansion follows this rule relative to their respective scalar DOF maps. Generic algorithms must obtain local-to-global indices from the finite-element space rather than reconstructing them from this formula. Other conforming families, including future $H(\mathrm{div})$ and $H(\mathrm{curl})$ spaces, may have different basis functions, mappings and DOF distributions.

Spatial tensor storage and index meaning

SpatialTensor uses a fixed-capacity stack array with runtime extents. Active entries occupy a contiguous prefix, with the last index varying fastest. For rank $R$ , extents $d_0,\ldots,d_{R-1}$ and indices $i_0,\ldots,i_{R-1}$ , the storage offset is

\[ p=\sum_{r=0}^{R-1} i_r\prod_{s=r+1}^{R-1}d_s. \]

The empty product equals one. An extent of zero makes the active tensor empty. This is the tensor analogue of row-major storage and differs from the column-major backing storage of SpatialMatrix.

Index roles are specified by each mathematical operation. In particular, the gradient of a matrix field has entries $G_{ijk}=\partial_k A_{ij}$ , and a constitutive tensor acts through $(\mathcal C A)_{ij}=\sum_{kl}\mathcal C_{ijkl}A_{kl}$ . These definitions do not depend on storage order. Tensor storage does not define a tensor-valued finite-element DOF convention; the current matrix solution API does not introduce higher-rank tensor solution spaces.

Consequences for assembly

A change of component numbering permutes coefficient vectors and the rows and columns of assembled operators. Such a change requires coordinated updates to basis evaluation, component access, DOF maps, persistence and interoperability. Changing memory order alone does not redefine components.

For uncoupled matrix mass and stiffness forms in the componentwise spaces, interleaved numbering gives $A_h=A_{\mathrm{scalar}}\otimes I_{mn}$ . With $N$ scalar local basis functions, full dense traversal examines $(mn)^2N^2$ entries although only $mnN^2$ component pairs can contribute. Switching between row-major and column-major component numbering does not remove this structural work. Assembly must preserve the form's component coupling to avoid traversing absent pairs; coupled constitutive forms require their actual coupling pattern.