1%%%%%%%%%%%%%%%%%%% 2% XLiFE++ is an extended library of finite elements written in C++ 3% Copyright (C) 2014 Lunéville, Eric; Kielbasiewicz, Nicolas; Lafranche, Yvon; Nguyen, Manh-Ha; Chambeyron, Colin 4% 5% This program is free software: you can redistribute it and/or modify 6% it under the terms of the GNU General Public License as published by 7% the Free Software Foundation, either version 3 of the License, or 8% (at your option) any later version. 9% This program is distributed in the hope that it will be useful, 10% but WITHOUT ANY WARRANTY; without even the implied warranty of 11% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 12% GNU General Public License for more details. 13% You should have received a copy of the GNU General Public License 14% along with this program. If not, see <http://www.gnu.org/licenses/>. 15%%%%%%%%%%%%%%%%%%% 16 17\section{Matrix} 18 19The purpose of the {\class Matrix} class is mainly to deal with complex or real \textbf{dense} 20matrices. In particular, this class is used in the definition of the user functions (see the 21section Function). This class is compliant with the {\class Vector} class. Although, it can 22deal with matrices of anything, it is only fully functional for real or complex matrices: 23\vspace{.2cm} 24\begin{lstlisting}[]{} 25Matrix<Real> rA; // an empty matrix 26Matrix<Real> rB(3,2); // a 3@$\times$@2 zeros matrix 27Matrix<Real> rC(3,2,1); // a 3@$\times$@2 ones matrix 28Vector<Real> w(3,2.5); // w=[2.5 2.5 2.5] 29Complex i(0,1); // the complex i 30Matrix<Real> cA(3,2,i); // a 3@$\times$@2 i matrix 31\end{lstlisting} 32\vspace{.2cm} 33 34It is possible to construct diagonal matrix from a {\class Vector} or a matrix from a {\class 35Vector} of {\class Vector}, to load (and save) a matrix from a file and to construct particular 36matrices ({\tt \_zeroMatrix, \_onesMatrix, \_idMatrix, \_hilbertMatrix}): 37\vspace{.2cm} 38\begin{lstlisting}[]{} 39Vector<Real> u(3,2.); // vector [2. 2. 2.] 40Matrix<Real> rA(u); // 3@$\times$@3 matrix with u as diagonal 41Matrix<Real> rB("mat.dat"); // matrix loaded from "mat.dat" file 42Matrix<Real> rO(3,_zeroMatrix); // a 3@$\times$@3 zeros matrix 43Matrix<Real> r1(3,_onesMatrix); // a 3@$\times$@3 ones matrix 44Matrix<Real> rI(3,_idMatrix); // a 3@$\times$@3 identity matrix 45Matrix<Real> rH(3,_hilbertMatrix); // the 3@$\times$@3 Hilbert matrix 46\end{lstlisting} 47\vspace{.2cm} 48Construction of complex matrix from real data are allowed (automatic cast). But the contrary 49is not. \\ 50 51There are some functions to access to the matrix properties: 52\vspace{.2cm} 53\begin{lstlisting}[]{} 54numberOfRows(), numberOfColumns() 55isSymmetric(), isSkewSymmetric(), isSelfAdjoint(), isSkewAdjoint() 56\end{lstlisting} 57\vspace{.2cm} 58ans some utilities to access to a coefficient, a row or a column or the diagonal of the matrix 59: 60\begin{lstlisting}[]{} 61Matrix<Real> A(2,2,1); //a 2@$\times$@2 ones matrix 62A(1,1)=2.; //change the coefficient A11 63Vector<Real> r=A.row(1); //first row of A 64r=A.column(2); //second column of A 65r=A.diag(); //diagonal of A 66A.column(1,r); //assign a vector to the first column 67A.row(2,r); //assign a vector to the second row 68A.diag(r); //assign a vector to the diagonal 69\end{lstlisting} 70\vspace{.2cm} 71All these functions support automatic cast from real to complex but not the contrary. \\ 72Advanced users can use member functions {\cmd begin()} and {\cmd end()} returning respectively 73iterators (or const iterators) to the beginning and the end of the Matrix. The data values of the matrix 74are stored according to the C convention, i.e. row-wise. \\ 75 76There are also generalized access tools either to extract submatrix ({\cmd get()} or {\cmd operator()}) or to set submatrix of matrix ({\cmd set()}): 77\vspace{.2cm} 78\begin{lstlisting}[]{} 79Matrix<Real> M(3,3); 80for(Number i=1;i<=3;i++) 81 for(Number j=1;j<=3;j++) 82 M(i,j)=i+j; //M=[2 3 4; 3 4 5; 4 5 6] 83Matrix<Real> N= M.get(2,3,2,3); //N=[4 5; 5 6] 84//N=M(2,3,2,3) gives the same 85Vector<Number> is(2); 86is(1)=1;is(2)=3; 87N= M.get(is,is); //N=[2 4; 4 6] 88//N=M(is,is) gives the same 89Matrix<Real> Z(2,2,0); //Z=[0 0; 0 0] 90M.set(1,2,1,2,Z); //M=[0 0 4; 0 0 5; 4 5 6] 91Matrix<Real> U(2,2,1); //U=[1 1; 1 1] 92M.set(is,is,Z); //M=[1 0 1; 0 0 5; 1 5 1] 93\end{lstlisting} 94\vspace{.2cm} 95Other syntaxes are proposed, see the developer's documentation.\\ 96\vspace{.2cm} 97 98Besides, the {\class Matrix} class proposes some transformations either as internal functions 99or external functions: 100\vspace{.2cm} 101\begin{lstlisting}[]{} 102Matrix<Real> A(2,2,1),B; 103Matrix<Complex> C(2,2,i),D; 104A.transpose(); // self transposition of A 105B=transpose(A); // transposition of A, A not changed 106C.adjoint(); // self transposition and conjugate C 107D=adjoint(A); // transpose and conjugate, C not changed 108B=diag(A); // from diagonal of A to a diagonal matrix 109A=real(C); // real part of C 110B=imag(C); // imaginary part of C 111D=conj(C); // conjugate of C 112D=cmplx(A); // forced casting from real to complex 113Real n2=norm2(A); //Frobenius norm 114Real ninf=norminfty(A);//infinite norm 115\end{lstlisting} 116\vspace{.2cm} 117 118Standard algebraic operations (+=,-=,*=,/=,+,-,*,/) are supported by the {\class Matrix} class. 119Some shortcuts are also possible, for instance a matrix plus a scalar, a scalar plus a matrix, 120... Automatic cast from real to complex is supported. There is no comparison operator.\\ 121\\ 122The {\class Matrix} proposes some solvers: 123\begin{lstlisting}[deletekeywords={[3] row}] 124Matrix<Real> A 125RealVector B; 126... 127// solve AX=B or AXs=Bxs using Gauss reduction 128gaussSolver(A, B, piv, row); 129gaussMultipleSolver(A, B, nbrhs, piv, row); 130//inverse of a square matrix 131RealMatrix invA=inverse(A); 132// QR factorization 133RealMatrix Q,R; 134qr(A,Q,R); 135// SVD factorization A= U S V*, if A a (m,n)-matrix, U is a (m,r)-matrix, V a (n,r)-matrix and S a r-vector where r=min(m,n) 136RealMatrix U, V; 137RealVector S; //singular values 138svdMat(A, U, S, V); 139\end{lstlisting} 140\vspace{.2cm} 141\begin{remark} 142 QR and SVD are available only if Eigen library is set on. 143\end{remark} 144\begin{remark} 145It is also possible to deal with matrix of matrices, for instance: 146\vspace{.2cm} 147\begin{lstlisting}[]{} 148Matrix<Real> ones(2,2,1); //a 2@$\times$@2 ones matrix 149Matrix<Matrix<Real> > A(2,2,ones); //a 2@$\times$@2 matrix of 2@$\times$@2 ones matrix 150\end{lstlisting} 151\vspace{.2cm} 152but all operations are not supported for such matrices!\\ 153\end{remark} 154 155To avoid explicit templates in user program, the following aliases are provided: 156\begin{itemize}[itemsep=-5pt] 157\item {\class RealMatrix} stands for {\class Matrix<Real>}, 158\item {\class ComplexMatrix} stands for {\class Matrix<Complex>}. 159\item {\class RealMatrices} stands for {\class Matrix<Matrix<Real> >}, 160\item {\class ComplexMatrices} stands for {\class Matrix<Matrix<Complex> >}. 161\end{itemize} 162 163