The sparse linear systems module implements solution methods for systems of consistent linear equations. Systems of linear equations take the form:
\[ Ax = b \]
where \(x \in \mathbb{R}^{n}\), \(A \in \mathbb{R}^{m \times n} \) and \(b \in \mathbb{R}^{m}\). This module accepts sparse matrix input formats for \(A\) and \(b\). We assume that there are no rows of \(A\) where all elements are zero.
The algorithms implemented in this module can handle large sparse square linear systems. Currently, the algorithms implemented in this module solve the linear system using direct or iterative methods.
linear_solver_sparse( tbl_source_lhs, tbl_source_rhs, tbl_result, lhs_row_id, lhs_col_id, lhs_value, rhs_row_id, rhs_value, grouping_cols := NULL, optimizer := 'direct', optimizer_params := 'algorithm = llt' )
Arguments
The name of the table containing the left hand side matrix. For the LHS matrix, the input data is expected to be of the following form:
{TABLE|VIEW} sourceName ( ... row_id FLOAT8, col_id FLOAT8, value FLOAT8, ... )
Each row represents a single equation. The rhs columns refer to the right hand side of the equations and the lhs columns refer to the multipliers on the variables on the left hand side of the same equations.
TEXT. The name of the table containing the right hand side vector. For the RHS matrix, the input data is expected to be of the following form:
{TABLE|VIEW} <em>sourceName</em> ( ... <em>row_id</em> FLOAT8, <em>value</em> FLOAT8 ... )
Each row represents a single equation. The rhs columns refer to the right hand side of the equations while the lhs columns refers to the multipliers on the variables on the left hand side of the same equations.
TEXT. The name of the table where the output is saved. Output is stored in the tabled named by the tbl_result argument. The table contains the following columns. The output contains the following columns:
solution | FLOAT8[]. The solution is an array with the variables in the same order as that provided as input in the 'left_hand_side' column name of the 'source_table' |
---|---|
residual_norm | FLOAT8. Scaled residual norm, defined as \( \frac{|Ax - b|}{|b|} \). This value is an indication of the accuracy of the solution. |
iters | INTEGER. Number of iterations required by the algorithm (only applicable for iterative algorithms) . The output is NULL for 'direct' methods. |
TEXT. The name of the column (in tbl_source_lhs) storing the 'col id' of the equations.
TEXT. The name of the column (in tbl_source_lhs) storing the 'value' of the equations.
TEXT. The name of the column (in tbl_source_rhs) storing the 'col id' of the equations.
TEXT. The name of the column (in tbl_source_rhs) storing the 'value' of the equations.
INTEGER. The number of variables in the linear system equations.
TEXT, default: 'direct'. Type of optimizer.
For each optimizer, there are specific parameters that can be tuned for better performance.
There are several algorithms that can be classified as 'direct' methods of solving linear systems. Madlib functions provide various algorithmic options available for users.
The following table provides a guideline on the choice of algorithm based on conditions on the A matrix, speed of the algorithms and numerical stability.
Algorithm | Conditions on A | Speed | Memory ---------------------------------------------------------- llt | Sym. Pos Def | ++ | --- ldlt | Sym. Pos Def | ++ | --- For speed '++' is faster than '+', which is faster than '-'. For accuracy '+++' is better than '++'. For memory, '-' uses less memory than '--'. Note: ldlt is often preferred over llt
There are several algorithms that can be classified as 'iterative' methods of solving linear systems. Madlib functions provide various algorithmic options available for users.
The following table provides a guideline on the choice of algorithm based on conditions on the A matrix, speed of the algorithms and numerical stability.
Algorithm | Conditions on A | Speed | Memory | Convergence ---------------------------------------------------------------------- cg-mem | Sym. Pos Def | +++ | - | ++ bicgstab-mem | Square | ++ | - | + precond-cg-mem | Sym. Pos Def | ++ | - | +++ precond-bicgstab-mem | Square | + | - | ++ For memory, '-' uses less memory than '--'. For speed, '++' is faster than '+'.
Algorithm Details:
cg-mem | In memory conjugate gradient with diagonal preconditioners. |
---|---|
bicgstab-mem | Bi-conjugate gradient (equivalent to performing CG on the least squares formulation of Ax=b) with incomplete LU preconditioners. |
precond-cg-mem | In memory conjugate gradient with diagonal preconditioners. |
bicgstab-mem | Bi-conjugate gradient (equivalent to performing CG on the least squares formulation of Ax=b) with incomplete LU preconditioners. |
Termination tolerance (applicable only for iterative methods) which determines the stopping criterion (with respect to residual norm) for iterative methods.
SELECT madlib.linear_solver_sparse();
DROP TABLE IF EXISTS sparse_linear_systems_lhs; CREATE TABLE sparse_linear_systems_lhs ( rid INTEGER NOT NULL, cid INTEGER, val DOUBLE PRECISION ); DROP TABLE IF EXISTS sparse_linear_systems_rhs; CREATE TABLE sparse_linear_systems_rhs ( rid INTEGER NOT NULL, val DOUBLE PRECISION ); INSERT INTO sparse_linear_systems_lhs(rid, cid, val) VALUES (0, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1); INSERT INTO sparse_linear_systems_rhs(rid, val) VALUES (0, 10), (1, 20), (2, 30);
SELECT madlib.linear_solver_sparse( 'sparse_linear_systems_lhs', 'sparse_linear_systems_rhs', 'output_table', 'rid', 'cid', 'val', 'rid', 'val', 4 );
\x on SELECT * FROM output_table;Result:
--------------------+------------------------------------- solution | {10,20,30,0} residual_norm | 0 iters | NULL
DROP TABLE IF EXISTS output_table; SELECT madlib.linear_solver_sparse( 'sparse_linear_systems_lhs', 'sparse_linear_systems_rhs', 'output_table', 'rid', 'cid', 'val', 'rid', 'val', 4, NULL, 'direct', 'algorithm=llt' );
DROP TABLE IF EXISTS output_table; SELECT madlib.linear_solver_sparse( 'sparse_linear_systems_lhs', 'sparse_linear_systems_rhs', 'output_table', 'rid', 'cid', 'val', 'rid', 'val', 4, NULL, 'iterative', 'algorithm=cg-mem, toler=1e-5' );