complex_matrices.h
Go to the documentation of this file.
1 //LIC// ====================================================================
2 //LIC// This file forms part of oomph-lib, the object-oriented,
3 //LIC// multi-physics finite-element library, available
4 //LIC// at http://www.oomph-lib.org.
5 //LIC//
6 //LIC// Version 1.0; svn revision $LastChangedRevision: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 //This header file contains classes and inline function definitions for
31 //matrices of complex numbers and their derived types
32 
33 //Include guards to prevent multiple inclusion of the header
34 #ifndef OOMPH_COMPLEX_MATRICES_HEADER
35 #define OOMPH_COMPLEX_MATRICES_HEADER
36 
37 // Config header generated by autoconfig
38 #ifdef HAVE_CONFIG_H
39 #include <oomph-lib-config.h>
40 #endif
41 
42 #ifdef OOMPH_HAS_MPI
43 #include "mpi.h"
44 #endif
45 
46 //Of course we need the complex type
47 #include <complex>
48 
49 //Also need the standard matrices header
50 #include "matrices.h"
51 
52 namespace oomph
53 {
54 
55 //===================================================================
56 /// \short Abstract base class for matrices of complex doubles -- adds
57 /// abstract interfaces for solving, LU decomposition and
58 /// multiplication by vectors.
59 //===================================================================
61 {
62  public:
63 
64  /// (Empty) constructor.
66 
67  /// Broken copy constructor
69  {
70  BrokenCopy::broken_copy("ComplexMatrixBase");
71  }
72 
73  /// Broken assignment operator
75  {
76  BrokenCopy::broken_assign("ComplexMatrixBase");
77  }
78 
79  /// Return the number of rows of the matrix
80  virtual unsigned long nrow() const=0;
81 
82  /// Return the number of columns of the matrix
83  virtual unsigned long ncol() const=0;
84 
85  /// virtual (empty) destructor
86  virtual ~ComplexMatrixBase() { }
87 
88  /// \short Round brackets to give access as a(i,j) for read only
89  /// (we're not providing a general interface for component-wise write
90  /// access since not all matrix formats allow efficient direct access!)
91  virtual std::complex<double> operator()(const unsigned long &i,
92  const unsigned long &j)
93  const=0;
94 
95  /// \short LU decomposition of the matrix using the appropriate
96  /// linear solver. Return the sign of the determinant
97  virtual int ludecompose()
98  {
99  throw OomphLibError(
100  "ludecompose() has not been written for this matrix class\n",
101  OOMPH_CURRENT_FUNCTION,
102  OOMPH_EXCEPTION_LOCATION);
103 
104  /// Dummy return
105  return 1;
106  }
107 
108  /// \short LU backsubstitue a previously LU-decomposed matrix;
109  /// The passed rhs will be over-written with the solution vector
110  virtual void lubksub(Vector<std::complex<double> > &rhs)
111  {
112  throw OomphLibError(
113  "lubksub() has not been written for this matrix class\n",
114  OOMPH_CURRENT_FUNCTION,
115  OOMPH_EXCEPTION_LOCATION);
116  }
117 
118 
119  /// \short Complete LU solve (replaces matrix by its LU decomposition
120  /// and overwrites RHS with solution). The default should not need
121  /// to be over-written
122  virtual void solve(Vector<std::complex<double> > &rhs);
123 
124  /// \short Complete LU solve (Nothing gets overwritten!). The default should
125  /// not need to be overwritten
126  virtual void solve(const Vector<std::complex<double> > &rhs,
127  Vector<std::complex<double> > &soln);
128 
129  /// \short Find the residual, i.e. r=b-Ax the residual
130  virtual void residual(const Vector<std::complex<double> > &x,
131  const Vector<std::complex<double> > &b,
132  Vector<std::complex<double> > &residual)=0;
133 
134  /// \short Find the maximum residual r=b-Ax -- generic version, can be
135  /// overloaded for specific derived classes where the
136  /// max. can be determined "on the fly"
137  virtual double max_residual(const Vector<std::complex<double> > &x,
138  const Vector<std::complex<double> > &rhs)
139  {
140  unsigned long n=rhs.size();
142  residual(x,rhs,res);
143  double ans=0.0;
144  for (unsigned long i=0;i<n;i++)
145  {
146  ans = std::max(ans,std::abs(res[i]));
147  }
148  return ans;
149  }
150 
151  /// \short Multiply the matrix by the vector x: soln=Ax.
152  virtual void multiply(const Vector<std::complex<double> > &x,
153  Vector<std::complex<double> > &soln)=0;
154 
155  /// \short Multiply the transposed matrix by the vector x: soln=A^T x
156  virtual void multiply_transpose(const Vector<std::complex<double> > &x,
157  Vector<std::complex<double> > &soln)=0;
158 
159 };
160 
161 
162 
163 //=================================================================
164 /// \short Class of matrices containing double complex, and stored as a
165 /// DenseMatrix<complex<double> >, but with solving functionality inherited
166 /// from the abstract ComplexMatrix class.
167 //=================================================================
169  public DenseMatrix<std::complex<double> >
170  {
171  private:
172 
173  /// Pointer to storage for the index of permutations in the LU solve
175 
176  /// Pointer to storage for the LU decomposition
177  std::complex<double>* LU_factors;
178 
179  /// Boolean flag used to decided whether the LU decomposition is stored
180  /// separately, or not
182 
183  /// Function that deletes the storage for the LU_factors, if it is
184  /// not the same as the matrix storage
185  void delete_lu_factors();
186 
187  public:
188 
189  /// Empty constructor, simply assign the lengths N and M to 0
190  DenseComplexMatrix(): DenseMatrix<std::complex<double> >(),
191  Index(0), LU_factors(0),
192  Overwrite_matrix_storage(false) {}
193 
194  /// Constructor to build a square n by n matrix.
195  DenseComplexMatrix(const unsigned long &n) :
196  DenseMatrix<std::complex<double> >(n), Index(0), LU_factors(0),
197  Overwrite_matrix_storage(false) {}
198 
199  /// Constructor to build a matrix with n rows and m columns.
200  DenseComplexMatrix(const unsigned long &n, const unsigned long &m) :
201  DenseMatrix<std::complex<double> >(n,m), Index(0), LU_factors(0),
202  Overwrite_matrix_storage(false) {}
203 
204  /// \short Constructor to build a matrix with n rows and m columns,
205  /// with initial value initial_val
206  DenseComplexMatrix(const unsigned long &n, const unsigned long &m,
207  const std::complex<double> &initial_val) :
208  DenseMatrix<std::complex<double> >(n,m,initial_val),
209  Index(0), LU_factors(0),
210  Overwrite_matrix_storage(false) {}
211 
212 
213  /// Broken copy constructor
215  {
216  BrokenCopy::broken_copy("DenseComplexMatrix");
217  }
218 
219  /// Broken assignment operator
221  {
222  BrokenCopy::broken_assign("DenseComplexMatrix");
223  }
224 
225 
226  /// Return the number of rows of the matrix
227  inline unsigned long nrow() const
229 
230  /// Return the number of columns of the matrix
231  inline unsigned long ncol() const
233 
234  /// \short Overload the const version of the round-bracket access operator
235  /// for read-only access.
236  inline std::complex<double> operator()(const unsigned long &i,
237  const unsigned long &j)
238  const {return DenseMatrix<std::complex<double> >::get_entry(i,j);}
239 
240  /// \short Overload the non-const version of the round-bracket access
241  /// operator for read-write access
242  inline std::complex<double>& operator()(const unsigned long &i,
243  const unsigned long &j)
245 
246  /// Destructor, delete the storage for Index vector and LU storage (if any)
247  virtual ~DenseComplexMatrix();
248 
249  /// \short Overload the LU decomposition.
250  /// For this storage scheme the matrix storage will be over-writeen
251  /// by its LU decomposition
252  int ludecompose();
253 
254  /// \short Overload the LU back substitution. Note that the rhs will be
255  /// overwritten with the solution vector
256  void lubksub(Vector<std::complex<double> >&rhs);
257 
258  /// \short Find the residual of Ax=b, ie r=b-Ax for the
259  /// "solution" x.
260  void residual(const Vector<std::complex<double> >& x,
261  const Vector<std::complex<double> >& rhs,
262  Vector<std::complex<double> >& residual);
263 
264  /// \short Multiply the matrix by the vector x: soln=Ax
265  void multiply(const Vector<std::complex<double> > &x,
266  Vector<std::complex<double> > &soln);
267 
268  /// \short Multiply the transposed matrix by the vector x: soln=A^T x
269  void multiply_transpose(const Vector<std::complex<double> > &x,
270  Vector<std::complex<double> > &soln);
271 
272 };
273 
274 
275 //=================================================================
276 /// \short A class for compressed row matrices
277 //=================================================================
278 class CRComplexMatrix : public CRMatrix<std::complex<double> >,
279  public ComplexMatrixBase
280 {
281  private:
282 
283  /// \short Storage for the LU factors as required by SuperLU
284  void* F_factors;
285 
286  /// \short Info flag for the SuperLU solver
287  int Info;
288 
289  public:
290 
291  /// Default constructor
292  CRComplexMatrix() : CRMatrix<std::complex<double> >(), F_factors(0), Info(0)
293  {
294  // By default SuperLU solves linear systems quietly
296  }
297 
298  /// \short Constructor: Pass vector of values, vector of column indices,
299  /// vector of row starts and number of columns (can be suppressed
300  /// for square matrices)
301  CRComplexMatrix(const Vector<std::complex<double> >& value,
302  const Vector<int>& column_index,
303  const Vector<int>& row_start,
304  const unsigned long &n,
305  const unsigned long &m) :
306  CRMatrix<std::complex<double> >(value,column_index,row_start,n,m),
307  F_factors(0), Info(0)
308  {
309  // By default SuperLU solves linear systems quietly
311  }
312 
313 
314  /// Broken copy constructor
316  {
317  BrokenCopy::broken_copy("CRComplexMatrix");
318  }
319 
320 
321  /// Broken assignment operator
323  {
324  BrokenCopy::broken_assign("CRComplexMatrix");
325  }
326 
327  /// Destructor: Kill the LU decomposition if it has been computed
329 
330  /// \short Set flag to indicate that stats are to be displayed during
331  /// solution of linear system with SuperLU
333 
334  // \short Set flag to indicate that stats are not to be displayed during
335  /// the solve
337 
338  /// Return the number of rows of the matrix
339  inline unsigned long nrow() const
341 
342  /// Return the number of columns of the matrix
343  inline unsigned long ncol() const
345 
346  /// Overload the round-bracket access operator for read-only access
347  inline std::complex<double> operator()(const unsigned long &i,
348  const unsigned long &j)
349  const {return CRMatrix<std::complex<double> >::get_entry(i,j);}
350 
351  /// \short LU decomposition using SuperLU
352  int ludecompose();
353 
354  /// \short LU back solve for given RHS
355  void lubksub(Vector<std::complex<double> > &rhs);
356 
357  /// \short LU clean up (perhaps this should happen in the destructor)
358  void clean_up_memory();
359 
360  /// \short Find the residual to x of Ax=b, i.e. r=b-Ax
361  void residual(const Vector<std::complex<double> > &x,
362  const Vector<std::complex<double> > &b,
363  Vector<std::complex<double> > &residual);
364 
365  /// \short Multiply the matrix by the vector x: soln=Ax
366  void multiply(const Vector<std::complex<double> > &x,
367  Vector<std::complex<double> > &soln);
368 
369 
370  /// \short Multiply the transposed matrix by the vector x: soln=A^T x
371  void multiply_transpose(const Vector<std::complex<double> > &x,
372  Vector<std::complex<double> > &soln);
373 
374  protected:
375 
376  /// \short Flag to indicate if stats are to be displayed during
377  /// solution of linear system with SuperLU
379 
380 };
381 
382 
383 
384 //=================================================================
385 /// \short A class for compressed column matrices that store doubles
386 //=================================================================
388  public CCMatrix<std::complex<double> >
389 {
390  private:
391 
392  /// \short Storage for the LU factors as required by SuperLU
393  void* F_factors;
394 
395  /// \short Info flag for the SuperLU solver
396  int Info;
397 
398  protected:
399 
400  /// \short Flag to indicate if stats are to be displayed during
401  /// solution of linear system with SuperLU
403 
404  public:
405 
406  /// Default constructor
407  CCComplexMatrix() : CCMatrix<std::complex<double> >(), F_factors(0), Info(0)
408  {
409  // By default SuperLU solves linear systems quietly
411  }
412 
413  /// \short Constructor: Pass vector of values, vector of row indices,
414  /// vector of column starts and number of rows (can be suppressed
415  /// for square matrices). Number of nonzero entries is read
416  /// off from value, so make sure the vector has been shrunk
417  /// to its correct length.
418  CCComplexMatrix(const Vector<std::complex<double> >& value,
419  const Vector<int>& row_index,
420  const Vector<int>& column_start,
421  const unsigned long &n,
422  const unsigned long &m) :
423  CCMatrix<std::complex<double> >(value,row_index,column_start,n,m),
424  F_factors(0), Info(0)
425  {
426  // By default SuperLU solves linear systems quietly
428  }
429 
430  /// Broken copy constructor
432  {
433  BrokenCopy::broken_copy("CCComplexMatrix");
434  }
435 
436  /// Broken assignment operator
438  {
439  BrokenCopy::broken_assign("CCComplexMatrix");
440  }
441 
442  /// Destructor: Kill the LU factors if they have been setup.
444 
445  /// \short Set flag to indicate that stats are to be displayed during
446  /// solution of linear system with SuperLU
448 
449  // \short Set flag to indicate that stats are not to be displayed during
450  /// the solve
452 
453  /// Return the number of rows of the matrix
454  inline unsigned long nrow() const
456 
457  /// Return the number of columns of the matrix
458  inline unsigned long ncol() const
460 
461  /// \short Overload the round-bracket access operator to provide
462  /// read-only (const) access to the data
463  inline std::complex<double> operator()(const unsigned long &i,
464  const unsigned long &j)
465  const {return CCMatrix<std::complex<double> >::get_entry(i,j);}
466 
467  /// \short LU decomposition using SuperLU
468  int ludecompose();
469 
470  /// \short LU back solve for given RHS
471  void lubksub(Vector<std::complex<double> > &rhs);
472 
473  /// \short LU clean up (perhaps this should happen in the destructor)
474  void clean_up_memory();
475 
476  /// \short Find the residulal to x of Ax=b, ie r=b-Ax
477  void residual(const Vector<std::complex<double> > &x,
478  const Vector<std::complex<double> > &b,
479  Vector<std::complex<double> > &residual);
480 
481 
482  /// \short Multiply the matrix by the vector x: soln=Ax
483  void multiply(const Vector<std::complex<double> > &x,
484  Vector<std::complex<double> > &soln);
485 
486 
487  /// \short Multiply the transposed matrix by the vector x: soln=A^T x
488  void multiply_transpose(const Vector<std::complex<double> > &x,
489  Vector<std::complex<double> > &soln);
490 
491 };
492 }
493 
494 #endif
virtual double max_residual(const Vector< std::complex< double > > &x, const Vector< std::complex< double > > &rhs)
Find the maximum residual r=b-Ax – generic version, can be overloaded for specific derived classes wh...
void residual(const Vector< std::complex< double > > &x, const Vector< std::complex< double > > &b, Vector< std::complex< double > > &residual)
Find the residual to x of Ax=b, i.e. r=b-Ax.
void lubksub(Vector< std::complex< double > > &rhs)
LU back solve for given RHS.
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
int * row_start()
Access to C-style row_start array.
Definition: matrices.h:783
virtual void multiply_transpose(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)=0
Multiply the transposed matrix by the vector x: soln=A^T x.
Abstract base class for matrices of complex doubles – adds abstract interfaces for solving...
CRComplexMatrix(const Vector< std::complex< double > > &value, const Vector< int > &column_index, const Vector< int > &row_start, const unsigned long &n, const unsigned long &m)
Constructor: Pass vector of values, vector of column indices, vector of row starts and number of colu...
unsigned long nrow() const
Return the number of rows of the matrix.
DenseComplexMatrix(const unsigned long &n, const unsigned long &m)
Constructor to build a matrix with n rows and m columns.
bool Doc_stats_during_solve
Flag to indicate if stats are to be displayed during solution of linear system with SuperLU...
void multiply_transpose(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the transposed matrix by the vector x: soln=A^T x.
void residual(const Vector< std::complex< double > > &x, const Vector< std::complex< double > > &b, Vector< std::complex< double > > &residual)
Find the residulal to x of Ax=b, ie r=b-Ax.
cstr elem_len * i
Definition: cfortran.h:607
int ludecompose()
LU decomposition using SuperLU.
virtual ~ComplexMatrixBase()
virtual (empty) destructor
Vector< long > * Index
Pointer to storage for the index of permutations in the LU solve.
virtual ~CCComplexMatrix()
Destructor: Kill the LU factors if they have been setup.
int * row_index()
Access to C-style row index array.
Definition: matrices.h:2492
void multiply_transpose(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the transposed matrix by the vector x: soln=A^T x.
ComplexMatrixBase()
(Empty) constructor.
CCComplexMatrix(const Vector< std::complex< double > > &value, const Vector< int > &row_index, const Vector< int > &column_start, const unsigned long &n, const unsigned long &m)
Constructor: Pass vector of values, vector of row indices, vector of column starts and number of rows...
A class for compressed row matrices.
int ludecompose()
LU decomposition using SuperLU.
bool Doc_stats_during_solve
Flag to indicate if stats are to be displayed during solution of linear system with SuperLU...
void multiply_transpose(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the transposed matrix by the vector x: soln=A^T x.
virtual unsigned long nrow() const =0
Return the number of rows of the matrix.
DenseComplexMatrix(const DenseComplexMatrix &matrix)
Broken copy constructor.
int ludecompose()
Overload the LU decomposition. For this storage scheme the matrix storage will be over-writeen by its...
void disable_doc_stats()
the solve
CRComplexMatrix()
Default constructor.
std::complex< double > get_entry(const unsigned long &i, const unsigned long &j) const
Access function that will be called by the read-only round-bracket operator (const) ...
Definition: matrices.h:745
void clean_up_memory()
LU clean up (perhaps this should happen in the destructor)
void enable_doc_stats()
Set flag to indicate that stats are to be displayed during solution of linear system with SuperLU...
virtual ~CRComplexMatrix()
Destructor: Kill the LU decomposition if it has been computed.
std::complex< double > operator()(const unsigned long &i, const unsigned long &j) const
Overload the round-bracket access operator to provide read-only (const) access to the data...
void * F_factors
Storage for the LU factors as required by SuperLU.
int Info
Info flag for the SuperLU solver.
CCComplexMatrix(const CCComplexMatrix &matrix)
Broken copy constructor.
void clean_up_memory()
LU clean up (perhaps this should happen in the destructor)
unsigned long ncol() const
Return the number of columns of the matrix.
void multiply(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the matrix by the vector x: soln=Ax.
std::complex< double > operator()(const unsigned long &i, const unsigned long &j) const
Overload the const version of the round-bracket access operator for read-only access.
void disable_doc_stats()
the solve
void residual(const Vector< std::complex< double > > &x, const Vector< std::complex< double > > &rhs, Vector< std::complex< double > > &residual)
Find the residual of Ax=b, ie r=b-Ax for the "solution" x.
void operator=(const ComplexMatrixBase &)
Broken assignment operator.
void operator=(const DenseComplexMatrix &)
Broken assignment operator.
virtual std::complex< double > operator()(const unsigned long &i, const unsigned long &j) const =0
Round brackets to give access as a(i,j) for read only (we're not providing a general interface for co...
std::complex< double > get_entry(const unsigned long &i, const unsigned long &j) const
Access function that will be called by the read-only round-bracket operator (const) ...
Definition: matrices.h:2448
A class for compressed column matrices: a sparse matrix format The class is passed as the MATRIX_TYPE...
Definition: matrices.h:2377
virtual void residual(const Vector< std::complex< double > > &x, const Vector< std::complex< double > > &b, Vector< std::complex< double > > &residual)=0
Find the residual, i.e. r=b-Ax the residual.
int Info
Info flag for the SuperLU solver.
std::complex< double > & entry(const unsigned long &i, const unsigned long &j)
The access function that will be called by the read-write round-bracket operator. ...
Definition: matrices.h:470
ComplexMatrixBase(const ComplexMatrixBase &matrix)
Broken copy constructor.
DenseComplexMatrix(const unsigned long &n, const unsigned long &m, const std::complex< double > &initial_val)
Constructor to build a matrix with n rows and m columns, with initial value initial_val.
std::complex< double > * LU_factors
Pointer to storage for the LU decomposition.
void operator=(const CCComplexMatrix &)
Broken assignment operator.
virtual void multiply(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)=0
Multiply the matrix by the vector x: soln=Ax.
void * F_factors
Storage for the LU factors as required by SuperLU.
unsigned long nrow() const
Return the number of rows of the matrix.
CRComplexMatrix(const CRComplexMatrix &matrix)
Broken copy constructor.
unsigned long ncol() const
Return the number of columns of the matrix.
A class for compressed row matrices, a sparse storage format Once again the recursive template trick ...
Definition: matrices.h:675
int * column_index()
Access to C-style column index array.
Definition: matrices.h:789
void lubksub(Vector< std::complex< double > > &rhs)
Overload the LU back substitution. Note that the rhs will be overwritten with the solution vector...
virtual ~DenseComplexMatrix()
Destructor, delete the storage for Index vector and LU storage (if any)
void lubksub(Vector< std::complex< double > > &rhs)
LU back solve for given RHS.
void broken_assign(const std::string &class_name)
Issue error message and terminate execution.
virtual void lubksub(Vector< std::complex< double > > &rhs)
LU backsubstitue a previously LU-decomposed matrix; The passed rhs will be over-written with the solu...
unsigned long nrow() const
Return the number of rows of the matrix.
std::complex< double > operator()(const unsigned long &i, const unsigned long &j) const
Overload the round-bracket access operator for read-only access.
DenseComplexMatrix(const unsigned long &n)
Constructor to build a square n by n matrix.
A class for compressed column matrices that store doubles.
Class for dense matrices, storing all the values of the matrix as a pointer to a pointer with assorte...
Definition: communicator.h:50
void multiply(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the matrix by the vector x: soln=Ax.
int * column_start()
Access to C-style column_start array.
Definition: matrices.h:2486
DenseComplexMatrix()
Empty constructor, simply assign the lengths N and M to 0.
CCComplexMatrix()
Default constructor.
void operator=(const CRComplexMatrix &)
Broken assignment operator.
std::complex< double > & operator()(const unsigned long &i, const unsigned long &j)
Overload the non-const version of the round-bracket access operator for read-write access...
virtual void solve(Vector< std::complex< double > > &rhs)
Complete LU solve (replaces matrix by its LU decomposition and overwrites RHS with solution)...
std::complex< double > get_entry(const unsigned long &i, const unsigned long &j) const
The access function the will be called by the read-only (const version) round-bracket operator...
Definition: matrices.h:480
std::complex< double > * value()
Access to C-style value array.
Definition: matrices.h:618
virtual int ludecompose()
LU decomposition of the matrix using the appropriate linear solver. Return the sign of the determinan...
unsigned long ncol() const
Return the number of columns of the matrix.
Class of matrices containing double complex, and stored as a DenseMatrix<complex<double> >...
void multiply(const Vector< std::complex< double > > &x, Vector< std::complex< double > > &soln)
Multiply the matrix by the vector x: soln=Ax.
void enable_doc_stats()
Set flag to indicate that stats are to be displayed during solution of linear system with SuperLU...
virtual unsigned long ncol() const =0
Return the number of columns of the matrix.