## Modelica_LinearSystems2.WorkInProgress.Math.Matrices

### Package Content

Name | Description
--- | ---
C_nullspace | Orthonormal nullspace of a complex matrix
C_QR | QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)
C_QR2 | QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)
C_rank | Rank of a complex matrix (computed with singular values)
C_singularValues | Compute singular values and left and right singular vectors of a complex matrix
C_solve2 | computes the solution to a complex system of linear equations A*X=B, using LU decomposition with partial pivoting and row interchanges
C_RQ | RQ decomposition of a rectangular complex matrix without column pivoting (A = R*Q)

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_nullspace

Orthonormal nullspace of a complex matrix

### Information

#### Syntax

>     Z = Matrices.nullspace(A);

>     (Z, nullity) = Matrices.nullspace(A);

#### Description

This function calculates a orthonormal basis Z=[z_1, z_2, ...] of the nullspace of a matrix A, i.e. A*z_i=0.

The nullspace is obtained by svd method. That is, matrix A is dcomposed into the matrices S, U, V:

>             H
>      A = USV

with the orthonormal matrices U and V and the matrix S with

>      S = [S1, 0]
>      S1 = [diag(s); 0]

with the singular values s={s1, s2, ..., sr} of A and r=rank(A). Note, that S has the same size as A. Since, U and V are orthonormal, we may write

>       H
>      U*A*V = [S1, 0].

Matrix S1 obviously has full column rank and therefore, the left n-r rows (n is the number of columns of A or S) of matrix V span a nullspace of A.

The nullity of matrix A is the dimension of the nullspace of A. In view of the above, it becomes clear that nullity holds

>      nullity = n - r

with

>      n = number of columns of matrix A
>      r = rank(A)

#### Example

>       A = [1, 2,  3, 1;
>            3, 4,  5, 2;
>           -1, 2, -3, 3];
>       (Z, nullity) = nullspace(A);

>       results in:

>       Z=[0.1715;
>         -0.686;
>          0.1715;
>          0.686]

>       nullity = 1

#### See also

Matrices.singularValues

Extends from Modelica.Icons.Function (Icon for a function).

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | input matrix

### Outputs

Type | Name | Description
--- | --- | ---
Complex | Z[size(A, 2), :] | orthonormal nullspace of matrix A
Integer | nullity | nullity, i.e. the dimension of the nullspace

### Modelica definition

    function C_nullspace "Orthonormal nullspace of a complex matrix"
      import Modelica_LinearSystems2;
      extends Modelica.Icons.Function;

      input Modelica_LinearSystems2.Math.Complex A[:,:] "input matrix";
      output Modelica_LinearSystems2.Math.Complex Z[size(A, 2),:]
        "orthonormal nullspace of matrix A";
      output Integer nullity "nullity, i.e. the dimension of the nullspace";

    protected
      Modelica_LinearSystems2.Math.Complex V[size(A, 2),size(A, 2)]
        "Right orthogonal matrix ";
      Real sigma[min(size(A, 1), size(A, 2))] "singular values";
      Integer rank "rank of matrix A";
      Real eps "tolerance for rank determination";
      Integer n=min(size(A, 1), size(A, 2));
      Integer i=n;

    algorithm
      (sigma,,V) := Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_singularValues(
                                                                           A);
      V := Modelica_LinearSystems2.Math.Complex.Internal.C_transpose(V);
      // rank computation
      eps := max(size(A, 1), size(A, 2))*max(sigma)*Modelica.Constants.eps;
      rank := 0;
      if n > 0 then
        while i > 0 loop
          if sigma[i] > eps then
            rank := i;
            i := 0;
          end if;
          i := i - 1;
        end while;
      end if;
      Z := V[:,rank + 1:size(A,2)];// nullspace computation
      nullity := size(A,2) - rank;// nullity

    end C_nullspace;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_QR

QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | Rectangular matrix with size(A,1) >= size(A,2)

### Outputs

Type | Name | Description
--- | --- | ---
Complex | Q[size(A, 1), size(A, 2)] | Rectangular matrix with orthonormal columns such that Q*R=A[:,p]
Complex | R[min(size(A, 1), size(A, 2)), size(A, 2)] | Square upper triangular matrix

### Modelica definition

    function C_QR
      "QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)"
      import Modelica_LinearSystems2.Math.Complex;

      input Complex A[:,:] "Rectangular matrix with size(A,1) >= size(A,2)";
      output Complex Q[size(A, 1),size(A, 2)]
        "Rectangular matrix with orthonormal columns such that Q*R=A[:,p]";
      output Complex R[min(size(A, 1), size(A, 2)),size(A, 2)]
        "Square upper triangular matrix";

    protected
      Integer nrow=size(A, 1);
      Integer ncol=size(A, 2);
      Integer minrowcol=min(nrow, ncol);

    algorithm
      assert(nrow >= ncol, "\nInput matrix A[" + String(nrow) + "," + String(ncol)
         +"] has more columns as rows. This is not allowed when calling Modelica.Matrices.C_QR(A).");

      if minrowcol > 0 then

        (Q,R) := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zgeqrf(
                                                                     A);

      else
        Q := fill(Complex(1), size(A, 1), 0);
        R := fill(Complex(0), 0,  0);
      end if;
    end C_QR;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_QR2

QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | Rectangular matrix with size(A,1) >= size(A,2)

### Outputs

Type | Name | Description
--- | --- | ---
Complex | Q[size(A, 1), size(A, 1)] | Rectangular matrix with orthonormal columns such that Q*R=A[:,p]
Complex | R[size(A, 1), size(A, 2)] | Square upper triangular matrix

### Modelica definition

    function C_QR2
      "QR decomposition of a rectangular complex matrix without column pivoting (A = Q*R)"
      import Modelica_LinearSystems2.Math.Complex;

      input Complex A[:,:] "Rectangular matrix with size(A,1) >= size(A,2)";
      output Complex Q[size(A, 1),size(A, 1)]
        "Rectangular matrix with orthonormal columns such that Q*R=A[:,p]";
      output Complex R[size(A,1),size(A, 2)] "Square upper triangular matrix";

    protected
      Complex A2[size(A, 1),size(A, 1)];
      Complex R2[size(A, 1),size(A, 1)];
      Integer nrow=size(A, 1);
      Integer ncol=size(A, 2);
      Integer minrowcol=min(nrow, ncol);

    algorithm
      assert(nrow >= ncol, "\nInput matrix A[" + String(nrow) + "," + String(ncol) + "] has more columns as rows.
    This is not allowed when calling Modelica.Matrices.C_QR(A).");

      if minrowcol > 0 then

        for l1 in 1:size(A, 1) loop
          for l2 in 1:size(A, 2) loop
            A2[l1, l2] := A[l1, l2];
          end for;
        end for;

        (Q,R2) := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zgeqrf(
                                                                      A2);
        for l1 in 1:size(R, 1) loop
          for l2 in l1:size(R, 2) loop
            R[l1, l2] := R2[l1, l2];
          end for;
        end for;

      else
        Q := fill(Complex(1), size(A, 1), 0);
        R := fill(Complex(0), 0, 0);
      end if;
    end C_QR2;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_rank

Rank of a complex matrix (computed with singular values)

### Information

Extends from Modelica.Icons.Function (Icon for a function).

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | Matrix
Real | eps | 0 | If eps > 0, the singular values are checked against eps; otherwise eps=max(size(A))*norm(A)*Modelica.Constants.eps is used

### Outputs

Type | Name | Description
--- | --- | ---
Integer | result | Rank of matrix A

### Modelica definition

    function C_rank
      "Rank of a complex matrix (computed with singular values)"
      import Modelica_LinearSystems2;
      extends Modelica.Icons.Function;
      input Modelica_LinearSystems2.Math.Complex A[:,:] "Matrix";
      input Real eps=0
        "If eps > 0, the singular values are checked against eps; otherwise eps=max(size(A))*norm(A)*Modelica.Constants.eps is used";
      output Integer result "Rank of matrix A";

    protected
      Integer n=min(size(A, 1), size(A, 2));
      Integer i=n;
      Real sigma[n];
      Real eps2;
    algorithm
      result := 0;
      if n > 0 then
        sigma := Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_singularValues(
                                                                        A);
        eps2 := if eps > 0 then eps else max(size(A))*sigma[1]*Modelica.Constants.eps;
        while i > 0 loop
          if sigma[i] > eps2 then
            result := i;
            i := 0;
          end if;
          i := i - 1;
        end while;
      end if;
    end C_rank;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_singularValues

Compute singular values and left and right singular vectors of a complex matrix

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | Square or rectangular matrix

### Outputs

Type | Name | Description
--- | --- | ---
Real | sigma[min(size(A, 1), size(A, 2))] | singular values
Complex | U[size(A, 1), size(A, 1)] | Left orthogonal matrix
Complex | VT[size(A, 2), size(A, 2)] | Transposed right orthogonal matrix

### Modelica definition

    function C_singularValues
      "Compute singular values and left and right singular vectors of a complex matrix"
      import Modelica_LinearSystems2.Math.Complex;

      input Complex A[:,:] "Square or rectangular matrix";
      output Real sigma[min(size(A, 1), size(A, 2))] "singular values";
      output Complex U[size(A, 1),size(A, 1)] "Left orthogonal matrix";
      output Complex VT[size(A, 2),size(A, 2)] "Transposed right orthogonal matrix";

    protected
      Integer info;
      Integer n=min(size(A, 1), size(A, 2)) "Number of singular values";

    algorithm
      if n > 0 then
        (sigma,U,VT,info) := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zgesvd(
                                                                                 A);
        assert(info == 0, "The numerical algorithm to compute the singular value decomposition did not converge");
      end if;
    end C_singularValues;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_solve2

computes the solution to a complex system of linear equations A*X=B, using LU decomposition with partial pivoting and row interchanges

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, size(A, 1)] |  |
Complex | B[size(A, 1), :] |  |

### Outputs

Type | Name | Description
--- | --- | ---
Complex | X[size(A, 1), size(B, 2)] |

### Modelica definition

    function C_solve2
      "computes the solution to a complex system of linear equations A*X=B, using LU decomposition with partial pivoting and row interchanges"
      import Modelica_LinearSystems2.Math.Complex;
      import Modelica_LinearSystems2.Math.Matrices;

      input Complex A[:,size(A, 1)];
      input Complex B[size(A, 1),:];
      output Complex X[size(A, 1),size(B, 2)];

    protected
      Integer info;
    algorithm
      (X,info) := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zgesv(
                                                                     A, B);
      assert(info == 0, "Solving a linear system of equations with function
    \"Matrices.C_solve2\" is not possible, because the system has either
    no or infinitely many solutions (A is singular).");
    end C_solve2;

## Modelica_LinearSystems2.WorkInProgress.Math.Matrices.C_RQ

RQ decomposition of a rectangular complex matrix without column pivoting (A = R*Q)

### Inputs

Type | Name | Default | Description
--- | --- | --- | ---
Complex | A[:, :] |  | Rectangular matrix with size(A,1) >= size(A,2)

### Outputs

Type | Name | Description
--- | --- | ---
Complex | R[size(A, 1), size(A, 1)] | Rectangular matrix with orthonormal columns such that Q*R=A[:,p]
Complex | Q[size(A, 1), size(A, 2)] | Square upper triangular matrix

### Modelica definition

    function C_RQ
      "RQ decomposition of a rectangular complex matrix without column pivoting (A = R*Q)"
      import Modelica_LinearSystems2.Math.Complex;

      input Complex A[:,:] "Rectangular matrix with size(A,1) >= size(A,2)";
      output Complex R[size(A, 1),size(A, 1)]
        "Rectangular matrix with orthonormal columns such that Q*R=A[:,p]";
      output Complex Q[size(A, 1),size(A, 2)] "Square upper triangular matrix";

    protected
      Integer nrow=size(A, 1);
      Integer ncol=size(A, 2);
      Complex RQ[size(A, 1),size(A, 2)];
      Complex tau[:];

    algorithm
      assert(nrow <= ncol, "\nInput matrix A[" + String(nrow) + "," + String(ncol)
         +"] has more rows as columns. This is not allowed when calling Modelica.Matrices.C_RQ(A).");

      if min(nrow,ncol) > 0 then
        (RQ, tau) := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zgerq2(
                                                                         A);

        for i in 1:nrow loop
          for ii in i:nrow loop
            R[i, ii] := RQ[i, ii+ncol-nrow];
          end for;
        end for;
        Q := Modelica_LinearSystems2.WorkInProgress.Math.LAPACK.zungrq(
                                                                 RQ, tau);
      else
        Q := fill(Complex(1), size(A, 1), 0);
        R := fill(Complex(0), 0,  0);
      end if;
    end C_RQ;