Class documentation of Concepts

Loading...
Searching...
No Matches

#include <PETSc.hh>

Inheritance diagram for concepts::PETScMat:
concepts::Matrix< Real > concepts::Operator< F > concepts::OutputOperator

Public Types

typedef Real value_type
 
typedef Realtype< Real >::type r_type
 Real type of data type.
 
typedef Cmplxtype< Real >::type c_type
 Complex type of data type.
 
typedef std::conditional< std::is_same< typenameRealtype< Real >::type, Real >::value, typenameRealtype< Real >::type, typenameCmplxtype< Real >::type >::type d_type
 Data type, depending if F is real or complex.
 
typedef _Matrix_iterator< Real, Real &, Real * > iterator
 
typedef _Matrix_iterator< Real, const Real &, const Real * > const_iterator
 
typedeftype
 Type of data, e.g. matrix entries.
 

Public Member Functions

 PETScMat (const Space< Real > &spcX, const Space< Real > &spcY)
 
 PETScMat (const Space< Real > &spc, BilinearForm< Real, Real > &bf)
 
 PETScMat (const SparseMatrix< Real > &matrix)
 Constructor. Copies the matrix from matrix.
 
virtual void operator() (const Function< double > &fncY, Function< double > &fncX)
 
virtual void operator() (const Function< std::complex< double > > &fncY, Function< std::complex< double > > &fncX)
 
void operator() (const Vector< double > &fncY, Vector< double > &fncX)
 
virtual void transpMult (const Vector< double > &fncY, Vector< double > &fncX)
 
virtual void transpMult (const Vector< std::complex< double > > &fncY, Vector< std::complex< double > > &fncX)
 
virtual Real operator() (const uint i, const uint j) const
 Returns entry with indices i and j.
 
virtual Realoperator() (const uint i, const uint j)
 
virtual const Space< double > & spaceX () const
 
virtual const Space< double > & spaceY () const
 
 operator Mat ()
 Conversion operator. Returns the PETSc matrix.
 
void storeMatlab (const char *name) const
 
const uint nofRows () const
 Number of rows.
 
const uint nofCols () const
 Number of columns.
 
iterator begin (uint r=0)
 
const_iterator begin (uint r=0) const
 
iterator end ()
 Iterator, standing behind last element.
 
const_iterator end () const
 Constant iterator, standing behind last element.
 
virtual void set (const uint i, const uint j, const Real value, const bool use_threshold=false, const Real threshold_value=1e-8)
 
virtual void add (const uint i, const uint j, const Real value, const bool use_threshold=false, const Real threshold_value=1e-8)
 
virtual void operator() (const Function< r_type > &fncY, Function< Real > &fncX)=0
 Computes fncX = A(fncY) where A is this matrix.
 
virtual void operator() (const Function< c_type > &fncY, Function< c_type > &fncX)=0
 
virtual void operator() (const Function< r_type > &fncY, Function< F > &fncX)
 
virtual void operator() ()
 
virtual bool operator== (const Matrix< Real > &otherMat) const
 
virtual void transpMult (const Vector< r_type > &fncY, Vector< Real > &fncX)=0
 Computes fncX = AT fncY where A is this matrix.
 
virtual void transpMult (const Vector< c_type > &fncY, Vector< c_type > &fncX)=0
 
virtual const uint dimX () const
 
virtual const uint dimY () const
 
virtual void show_messages ()
 

Static Public Member Functions

static void assembly (Matrix< Real > &dest, const Space< G > &spc, const BilinearForm< Real, G > &bf, const Real threshold=0.0)
 
static void assembly (Matrix< Real > &dest, Scan< Element< G > > *sc, const BilinearForm< Real, G > &bf, const Real threshold=0.0)
 
static void assembly (Matrix< Real > &dest, const Sequence< ElementWithCell< G > * > seq, const BilinearForm< Real, G > &bf, const Real threshold=0.0)
 
static void assembly (Matrix< Real > &dest, const Space< G > &spc, const Sequence< bool > &seq, const BilinearForm< Real, G > &bf, const Real threshold=0.0)
 
static void assembly (Matrix< Real > &dest, const Space< G > &spcX, const Space< G > &spcY, const BilinearForm< Real, G > &bf, const Real threshold=0.0, const bool single=false)
 
static void assembly (Matrix< Real > &dest, const Sequence< ElementWithCell< G > * > seqX, const Space< G > &spcY, const BilinearForm< Real, G > &bf, const Real threshold=0.0)
 
static void assembly (Matrix< Real > &dest, const BilinearForm< Real, G > &bf, const ElementPairList< G > &pairs)
 
static void assembly (Matrix< Real > &dest, const Space< G > &spcX, const Space< G > &spcY, const Sequence< bool > &seq, const BilinearForm< Real, G > &bf, const Real threshold=0.0, const bool single=false)
 
Timing Interface

These functions are used to get timings from class internal computations. The values are stored in a user defined concepts::InOutParameters structure in different arrays (see setTimings). These arrays can be grouped into a table for easier postprocessing with

table.addMap(concepts::ResultsTable::DOUBLE, "jacobian", output);
table.addMap(concepts::ResultsTable::DOUBLE, "whole_sumfact", output);
std::ofstream ofs("table.gnuplot");
ofs << std::setprecision(20);
table.print<concepts::ResultsTable::GNUPLOT>(ofs);
Set< F > makeSet(uint n, const F &first,...)
Definition set.hh:320
static void setTimings (InOutParameters *timings)
 
static bool timings ()
 

Protected Member Functions

std::ostream & info (std::ostream &os) const
 Returns information in an output stream.
 

Protected Attributes

uint dimX_
 Dimension of image space and the source space.
 
uint dimY_
 

Detailed Description

Interface to the sparse matrices from PETSc. Very usefull if you intend to use PETSc to solve your systems since this makes it possible to use the preconditioners of PETSc.

It is not possible tough, to combine different PETSc matrices. If you want to compute a linear combination of bilinear forms, then insert the linear combination of the bilinear forms into the constructor.

PETSc is a suite of data structures and routines for the scalable (parallel) solution of scientific applications modeled by partial differential equations. It employs the MPI standard for all message-passing communication.

Our interface to PETSc uses only its serial capablities.

See also
PETSc solvers in Concepts
documentation of PETSc on matrices
Satish Balay, Kris Buschelman, William D. Gropp, Dinesh Kaushik, Lois Curfman McInnes, and Barry F. Smith, PETSc home page, 2001.
Satish Balay, William D. Gropp, Lois Curfman McInnes, and Barry F. Smith. Efficient Management of Parallelism in Object Oriented Numerical Software Libraries. In E. Arge, A. M. Bruaset, and H. P. Langtangen, editors, Modern Software Tools in Scientific Computing, pages 163-202. Birkhauser Press, 1997.
Satish Balay, William D. Gropp, Lois Curfman McInnes, and Barry F. Smith. PETSc Users Manual. Technical Report ANL-95/11 - Revision 2.1.0, Argonne National Laboratory, 2001.
Author
Philipp Frauenfelder, 2001
Examples
hpFEM2d.cc.

Definition at line 148 of file PETSc.hh.

Member Typedef Documentation

◆ c_type

Complex type of data type.

Definition at line 45 of file matrix.hh.

◆ const_iterator

typedef _Matrix_iterator<Real , const Real &, const Real *> concepts::Matrix< Real >::const_iterator
inherited

Definition at line 51 of file matrix.hh.

◆ d_type

typedef std::conditional<std::is_same<typenameRealtype<Real >::type,Real >::value,typenameRealtype<Real >::type,typenameCmplxtype<Real >::type>::type concepts::Matrix< Real >::d_type
inherited

Data type, depending if F is real or complex.

Definition at line 48 of file matrix.hh.

◆ iterator

typedef _Matrix_iterator<Real , Real &, Real *> concepts::Matrix< Real >::iterator
inherited

Definition at line 50 of file matrix.hh.

◆ r_type

Real type of data type.

Definition at line 43 of file matrix.hh.

◆ type

template<class F >
typedef F concepts::Operator< F >::type
inherited

Type of data, e.g. matrix entries.

Definition at line 45 of file compositions.hh.

◆ value_type

typedef Real concepts::Matrix< Real >::value_type
inherited

Definition at line 41 of file matrix.hh.

Constructor & Destructor Documentation

◆ PETScMat() [1/2]

concepts::PETScMat::PETScMat ( const Space< Real > &  spcX,
const Space< Real > &  spcY 
)

Constructor. Creates an empty matrix.

◆ PETScMat() [2/2]

concepts::PETScMat::PETScMat ( const Space< Real > &  spc,
BilinearForm< Real, Real > &  bf 
)

Constructor. Computes the global matrix by assembling the element matrices.

Member Function Documentation

◆ add()

virtual void concepts::Matrix< Real >::add ( const uint  i,
const uint  j,
const Real  value,
const bool  use_threshold = false,
const Real  threshold_value = 1e-8 
)
virtualinherited

Addition operator Add the value to the entry (i,j) if not zero (for use_threshold = false), or if the absolute avlue is bigger than the threshold_value (for use_threshold = true)

◆ assembly() [1/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const BilinearForm< Real , G > &  bf,
const ElementPairList< G > &  pairs 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator uses the element pairs taken from pairs. For every two elements found in a ElementPair in pairs, the bilinear form is evaluated and the result assembled into dest.

◆ assembly() [2/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Sequence< ElementWithCell< G > * >  seq,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator does not compute element matrices for two different elements. The elements are taken from the element sequence seq.

◆ assembly() [3/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Sequence< ElementWithCell< G > * >  seqX,
const Space< G > &  spcY,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator computes also the element matrices for two different elements. The elements are taken from the element sequence seqX for the test space, and the trial space spcY is fully taken into account.

◆ assembly() [4/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Space< G > &  spc,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator does not compute element matrices for two different elements. The elements are taken from the space spc.

◆ assembly() [5/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Space< G > &  spc,
const Sequence< bool > &  seq,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0 
)
staticinherited

Assembly operator for dest using the bilinear form bf on space spc. This assembly operator does not compute element matrices for two different elements. The elements are computing on the cells that are flagged by seq.

◆ assembly() [6/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Space< G > &  spcX,
const Space< G > &  spcY,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0,
const bool  single = false 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator computes also the element matrices for two different elements (coming from test space spcX and trial space spcY).

◆ assembly() [7/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
const Space< G > &  spcX,
const Space< G > &  spcY,
const Sequence< bool > &  seq,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0,
const bool  single = false 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator computes also the element matrices for two different elements (coming from test space spcX and trial space spcY). The elements are computing on the cells that are flagged by seq.

◆ assembly() [8/8]

static void concepts::Matrix< Real >::assembly ( Matrix< Real > &  dest,
Scan< Element< G > > *  sc,
const BilinearForm< Real , G > &  bf,
const Real  threshold = 0.0 
)
staticinherited

Assembly operator for dest using the bilinear form bf. This assembly operator does not compute element matrices for two different elements. The elements are taken from the space scanner sc.

◆ begin() [1/2]

iterator concepts::Matrix< Real >::begin ( uint  r = 0)
inlineinherited

Iterator over the elements, standing at position (r,0).

Might be implemented differently for derived classes.

Definition at line 64 of file matrix.hh.

◆ begin() [2/2]

const_iterator concepts::Matrix< Real >::begin ( uint  r = 0) const
inlineinherited

Constant iterator over the elements, standing at position (r,0)

Might be implemented differently for derived classes.

Definition at line 71 of file matrix.hh.

◆ dimX()

template<class F >
virtual const uint concepts::Operator< F >::dimX ( ) const
inlinevirtualinherited

Returns the size of the image space of the operator (number of rows of the corresponding matrix)

Examples
hpFEM2d-simple.cc, hpFEM2d.cc, and matfileTutorial.cc.

Definition at line 93 of file compositions.hh.

◆ dimY()

template<class F >
virtual const uint concepts::Operator< F >::dimY ( ) const
inlinevirtualinherited

Returns the size of the source space of the operator (number of columns of the corresponding matrix)

Examples
matfileTutorial.cc.

Definition at line 98 of file compositions.hh.

◆ end() [1/2]

iterator concepts::Matrix< Real >::end ( )
inlineinherited

Iterator, standing behind last element.

Definition at line 66 of file matrix.hh.

◆ end() [2/2]

const_iterator concepts::Matrix< Real >::end ( ) const
inlineinherited

Constant iterator, standing behind last element.

Definition at line 73 of file matrix.hh.

◆ info()

std::ostream & concepts::PETScMat::info ( std::ostream &  os) const
protectedvirtual

Returns information in an output stream.

Reimplemented from concepts::Operator< F >.

◆ nofCols()

const uint concepts::Matrix< Real >::nofCols ( ) const
inlineinherited

Number of columns.

Definition at line 58 of file matrix.hh.

◆ nofRows()

const uint concepts::Matrix< Real >::nofRows ( ) const
inlineinherited

Number of rows.

Definition at line 56 of file matrix.hh.

◆ operator Mat()

concepts::PETScMat::operator Mat ( )
inline

Conversion operator. Returns the PETSc matrix.

Definition at line 187 of file PETSc.hh.

◆ operator()() [1/5]

◆ operator()() [2/5]

virtual void concepts::Matrix< Real >::operator() ( const Function< c_type > &  fncY,
Function< c_type > &  fncX 
)
pure virtualinherited

Application operator for complex function fncY.

Computes fncX = A(fncY) where A is this operator. fncX becomes complex.

In derived classes its enough to implement the operator() for complex Operator's. If a real counterpart is not implemented, the function fncY is splitted into real and imaginary part and the application operator for real functions is called for each. Then the result is combined.

If in a derived class the operator() for complex Operator's is not implemented, a exception is thrown from here.

Reimplemented from concepts::Operator< F >.

◆ operator()() [3/5]

template<class F >
virtual void concepts::Operator< F >::operator() ( const Function< r_type > &  fncY,
Function< F > &  fncX 
)
virtualinherited

Application operator for real function fncY.

Computes fncX = A(fncY) where A is this operator.

fncX becomes the type of the operator, for real data it becomes real, for complex data it becomes complex.

In derived classes its enough to implement the operator() for real Operator's. If a complex counterpart is not implemented, the function fncY is transformed to a complex function and then the application operator for complex functions is called.

If in a derived class the operator() for real Operator's is not implemented, a exception is thrown from here.

Reimplemented in aglowav::C2W< F >, aglowav2::C2W< F >, aglowav::W2C< F >, aglowav2::W2C< F >, aglowav::C2_tl2< F >, aglowav::C2tl2< F >, aglowav::CGt2< F >, aglowav::ComposeN< F >, aglowav2::Operator00< F >, bem::D< F >, bem::D_1< F >, concepts::TrivExtendRestrict< F >, sparseqr::GivensRotations< F >, vectorial::BlockOperator< F >, concepts::TrivExtendRestrict< Real >, sparseqr::GivensRotations< Real >, concepts::AfterIteration< F >, concepts::Compose< F, H >, concepts::DDSolver< F, G >, concepts::Multiple< F >, concepts::LiCoI< F >, concepts::LiCo< F >, concepts::DenseMatrix< F >, concepts::DiagonalMatrix< F >, concepts::Permutation< F >, concepts::Matrix< F >, concepts::Matrix< F::type >, and concepts::DiagonalSolver< F >.

◆ operator()() [4/5]

virtual Real & concepts::PETScMat::operator() ( const uint  i,
const uint  j 
)
virtual

Should return and allow access to entry with indices i and j. This is not possible with the interface of PETSc. Therefore, this method throws an exception.

Exceptions
MissingFeature

Implements concepts::Matrix< Real >.

◆ operator()() [5/5]

virtual Real concepts::PETScMat::operator() ( const uint  i,
const uint  j 
) const
virtual

Returns entry with indices i and j.

Implements concepts::Matrix< Real >.

◆ operator==()

virtual bool concepts::Matrix< Real >::operator== ( const Matrix< Real > &  otherMat) const
inlinevirtualinherited

Definition at line 179 of file matrix.hh.

◆ set()

virtual void concepts::Matrix< Real >::set ( const uint  i,
const uint  j,
const Real  value,
const bool  use_threshold = false,
const Real  threshold_value = 1e-8 
)
virtualinherited

Affectation operator Affet the value to the entry (i,j) if not zero (for use_threshold = false), or if the absolute value is bigger than the threshold_value (for use_threshold = true)

◆ setTimings()

static void concepts::Matrix< Real >::setTimings ( InOutParameters timings)
staticinherited

Sets the class to store the timing values in. Additionally, the timeCntr_ is reset to 0. This counter is used to fill in the values into the arrays listed below in subsequent calls. The following timings are taken and stored in timings:

  • evaluation of bilinear form in bilinear_form
  • application of T matrix in tmatrix_apply
  • assembling into global matrix in global_assembly

◆ show_messages()

template<class F >
virtual void concepts::Operator< F >::show_messages ( )
inlinevirtualinherited

Reimplemented in concepts::MumpsOverlap< F >.

Definition at line 100 of file compositions.hh.

◆ spaceX()

virtual const Space< double > & concepts::PETScMat::spaceX ( ) const
inlinevirtual

Definition at line 183 of file PETSc.hh.

◆ spaceY()

virtual const Space< double > & concepts::PETScMat::spaceY ( ) const
inlinevirtual

Definition at line 184 of file PETSc.hh.

◆ timings()

static bool concepts::Matrix< Real >::timings ( )
staticinherited

Returns true if the class is able to do timings. The ability to do timings depends on a compiler switch in matrix.cc file.

Member Data Documentation

◆ dimX_

template<class F >
uint concepts::Operator< F >::dimX_
protectedinherited

Dimension of image space and the source space.

Definition at line 104 of file compositions.hh.

◆ dimY_

template<class F >
uint concepts::Operator< F >::dimY_
protectedinherited

Definition at line 104 of file compositions.hh.


The documentation for this class was generated from the following file: