Blocks of right-hand sides

A matrix-like right-hand side holds one system per column. Every column is solved against the one factorization.

Factorize, then substitute a block

One factorization, three right-hand sides in one sweep.

\[\begin{split} A X = B, \quad A = \left(\begin{array}{ccc|c} 4 & 1 & 0 & 2 \\ 1 & 5 & 1 & 0 \\ 0 & 1 & 6 & 1 \\ \hline 2 & 0 & 1 & 7 \end{array}\right), \quad B = \begin{pmatrix} 14 & 7 & 8 \\ 14 & 7 & 3 \\ 24 & 8 & 6 \\ 33 & 10 & 5 \end{pmatrix}, \quad X = \begin{pmatrix} 1 & 1 & 2 \\ 2 & 1 & 0 \\ 3 & 1 & 1 \\ 4 & 1 & 0 \end{pmatrix} \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
const tfel::math::tmatrix<4, 3, double> B{14, 7, 8, 14, 7, 3, 24, 8, 6, 33, 10, 5};
tfel::math::tmatrix<4, 3, double> X;
tfel::math::fsarray<4, int> piv;

if (!tdls::factorize(A, piv)) return 1;

// a matrix-like right-hand side: one system per column, solved together
tdls::substitute(A, piv, B, X);

Factorize, then substitute a block in place

The same block, one buffer.

\[\begin{split} A X = B, \quad A = \left(\begin{array}{ccc|c} 4 & 1 & 0 & 2 \\ 1 & 5 & 1 & 0 \\ 0 & 1 & 6 & 1 \\ \hline 2 & 0 & 1 & 7 \end{array}\right), \quad B = \begin{pmatrix} 14 & 7 & 8 \\ 14 & 7 & 3 \\ 24 & 8 & 6 \\ 33 & 10 & 5 \end{pmatrix}, \quad X = \begin{pmatrix} 1 & 1 & 2 \\ 2 & 1 & 0 \\ 3 & 1 & 1 \\ 4 & 1 & 0 \end{pmatrix} \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
tfel::math::tmatrix<4, 3, double> Y{14, 7, 8, 14, 7, 3, 24, 8, 6, 33, 10, 5};
tfel::math::fsarray<4, int> piv;

if (!tdls::factorize(A, piv)) return 1;

// Y holds the three right-hand sides on entry and the three solutions on exit
tdls::substitute_inplace(A, piv, Y);

Solve a block

Factorization and block substitution in one call, two buffers.

\[\begin{split} A X = B, \quad A = \left(\begin{array}{ccc|c} 4 & 1 & 0 & 2 \\ 1 & 5 & 1 & 0 \\ 0 & 1 & 6 & 1 \\ \hline 2 & 0 & 1 & 7 \end{array}\right), \quad B = \begin{pmatrix} 14 & 7 & 8 \\ 14 & 7 & 3 \\ 24 & 8 & 6 \\ 33 & 10 & 5 \end{pmatrix}, \quad X = \begin{pmatrix} 1 & 1 & 2 \\ 2 & 1 & 0 \\ 3 & 1 & 1 \\ 4 & 1 & 0 \end{pmatrix} \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
const tfel::math::tmatrix<4, 3, double> B{14, 7, 8, 14, 7, 3, 24, 8, 6, 33, 10, 5};
tfel::math::tmatrix<4, 3, double> X;
tfel::math::fsarray<4, int> piv;

const bool ok = tdls::solve(A, piv, B, X);

Solve a block in place

Factorization and block substitution in one call, one buffer.

\[\begin{split} A X = B, \quad A = \left(\begin{array}{ccc|c} 4 & 1 & 0 & 2 \\ 1 & 5 & 1 & 0 \\ 0 & 1 & 6 & 1 \\ \hline 2 & 0 & 1 & 7 \end{array}\right), \quad B = \begin{pmatrix} 14 & 7 & 8 \\ 14 & 7 & 3 \\ 24 & 8 & 6 \\ 33 & 10 & 5 \end{pmatrix}, \quad X = \begin{pmatrix} 1 & 1 & 2 \\ 2 & 1 & 0 \\ 3 & 1 & 1 \\ 4 & 1 & 0 \end{pmatrix} \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
tfel::math::tmatrix<4, 3, double> Y{14, 7, 8, 14, 7, 3, 24, 8, 6, 33, 10, 5};
tfel::math::fsarray<4, int> piv;

// Y holds the three right-hand sides on entry and the three solutions on exit
const bool ok = tdls::solve_inplace(A, piv, Y);

A block of canonical columns

A matrix-like x: its columns receive the solutions of consecutive canonical columns, the tangent-operator pattern.

\[\begin{split} A X = \begin{pmatrix} e_0 & e_1 & e_2 \end{pmatrix}, \quad A = \left(\begin{array}{ccc|c} 4 & 1 & 0 & 2 \\ 1 & 5 & 1 & 0 \\ 0 & 1 & 6 & 1 \\ \hline 2 & 0 & 1 & 7 \end{array}\right) \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
tfel::math::tmatrix<4, 3, double> X;
tfel::math::fsarray<4, int> piv;

if (!tdls::factorize(A, piv)) return 1;

// a matrix-like x: its columns receive the solutions of e_0, e_1 and e_2
tdls::substitute_canonical(A, piv, 0, X);

Pass cutting

A wide block cut into passes of two columns, through the template argument of the entry point.

\[\begin{split} A X = B, \quad B = \left(\begin{array}{cc|cc|c} 14 & 7 & 8 & 3 & 4 \\ 14 & 7 & 3 & 5 & 1 \\ 24 & 8 & 6 & 2 & 0 \\ 33 & 10 & 5 & 7 & 2 \end{array}\right), \quad X = \left(\begin{array}{cc|cc|c} 1 & 1 & 2 & 0 & 1 \\ 2 & 1 & 0 & 1 & 0 \\ 3 & 1 & 1 & 0 & 0 \\ 4 & 1 & 0 & 1 & 0 \end{array}\right) \end{split}\]

tfel::math::tmatrix<4, 4, double> A{4, 1, 0, 2, 1, 5, 1, 0, 0, 1, 6, 1, 2, 0, 1, 7};
const tfel::math::tmatrix<4, 5, double> B{14, 7, 8, 3, 4, 14, 7,  3, 5, 1,
                                          24, 8, 6, 2, 0, 33, 10, 5, 7, 2};
tfel::math::tmatrix<4, 5, double> X;
tfel::math::fsarray<4, int> piv;

// five columns in passes of two: the pass width is the template argument
const bool ok = tdls::solve<2>(A, piv, B, X);

A runtime-sized block

A runtime-sized matrix as right-hand side of a runtime-sized system: the column count is read from the block.

\[\begin{split} A X = B, \quad A = \left(\begin{array}{ccc|cc} 5 & 1 & 0 & 0 & 1 \\ 1 & 6 & 1 & 0 & 0 \\ 0 & 1 & 7 & 1 & 0 \\ \hline 0 & 0 & 1 & 8 & 1 \\ 1 & 0 & 0 & 1 & 9 \end{array}\right), \quad B = \begin{pmatrix} 12 & 7 & 6 \\ 16 & 8 & 2 \\ 27 & 9 & 7 \\ 40 & 10 & 2 \\ 50 & 11 & 10 \end{pmatrix}, \quad X = \begin{pmatrix} 1 & 1 & 1 \\ 2 & 1 & 0 \\ 3 & 1 & 1 \\ 4 & 1 & 0 \\ 5 & 1 & 1 \end{pmatrix} \end{split}\]

tfel::math::matrix<double> A = {
    {5, 1, 0, 0, 1}, {1, 6, 1, 0, 0}, {0, 1, 7, 1, 0}, {0, 0, 1, 8, 1}, {1, 0, 0, 1, 9}};
const tfel::math::matrix<double> B = {
    {12, 7, 6}, {16, 8, 2}, {27, 9, 7}, {40, 10, 2}, {50, 11, 10}};
tfel::math::matrix<double> X(5, 3);
tfel::math::vector<int> piv(5);

// the column count is read from the block at run time
const bool ok = tdls::solve(A, piv, B, X);