Quick Guide to Solving A Matrix Problem with Autorej in Python

 

Assume you wish to solve A*x = b:

1.       import autorej as rej                  (or whatever name you want)

2.       Create the matrix, A, preferably as  numpy.array. The number of rows, m, and columns, n, is arbitrary.

3.       Create the right-hand-side column of data, b, as a a 1-D array.

4.       If the data is to be read from a file, you may use rej.read_problem(filename). The matrix must be one row per line, with the corresponding elements of b appended at the end. So the data will be laid out like this, for a 4 by 4 problem:

A  A  A  A  b

A  A  A  A  b

A  A  A  A  b

A  A  A  A  b

5.       Create an autorej “Problem”:

P=rej.Problem(A,b)

 

6.       To solve A *x = b:

x = P.solve()

 

7.       To solve A * x = b with the constraint that all elements of x be 0 or positive:

x = P.solve_nonneg()

 

8.  To solve with one equation obeyed exactly (i.e. a constraint equation) put the equation coefficients in a 1-D array, called Row here, and call, for example:

x = P.solve_constrained(Row,rhs_value)

 

9.       After the problem is solved you may access data about the problem:

 

10.   To get the Singular Value Decomposition of A:

U,S,V=P.get_svd()

11.   To get the singular values of A:

Sing=P.get_singular_values()

 

12.   To get the standard condition number of the matrix A:

cond=P.get_condition_number()

13.   To get the estimated “usable rank” determined by autorej’s algorithms:

rank=P.get_usable_rank()

14.   To get the autorej’s estimated RMS error in b:

sigma=P.get_sigma()

 

15.   To get the Tikhonov regularization parameter determined by autorej:

lamda=P.get_lambda()

 

16.   To get the effective singular values used by autorej’s solution.

RegSing = P.get_regularized_singular_values()

 

17.   To get the elements of UT*b which are samples of the noise in b.

N=P.get_noise_samples()

 

18.   to get the Picard Condition vector, S+ * UT * b

PCV=P.get_PCV()

 

19.   To get the pseudoinverse matrix:

Pinv=P.get_pseudoinverse()

 

20.   To get the generalized pseudoinverse matrix with singular values less than epsilon zeroed:

Ginv=P.get_generalized_inverse(epsilon)

 

21.   To get the inverse of the regularized version of  the original matrix, A:

Ainv=P.get_regularized_inverse()

 

Note that using this inverse matrix to compute solutions to other similar sized matrix problems is not recommended unless you are confident the new problem has similar noise characteristic in both the matrix A and the right hand side vector b.