TrilinosCouplings  Development
 All Classes Namespaces Files Functions Pages
Typedefs | Functions
example_Poisson.cpp File Reference

Example solution of a Poisson equation on a hexahedral mesh using nodal (Hgrad) elements. More...

#include "TrilinosCouplings_config.h"
#include "TrilinosCouplings_Pamgen_Utils.hpp"
#include "Intrepid_FunctionSpaceTools.hpp"
#include "Intrepid_CellTools.hpp"
#include "Intrepid_ArrayTools.hpp"
#include "Intrepid_HGRAD_HEX_C1_FEM.hpp"
#include "Intrepid_RealSpaceTools.hpp"
#include "Intrepid_DefaultCubatureFactory.hpp"
#include "Intrepid_Utils.hpp"
#include "Epetra_Time.h"
#include "Epetra_Map.h"
#include "Epetra_SerialComm.h"
#include "Epetra_FECrsMatrix.h"
#include "Epetra_FEVector.h"
#include "Epetra_Import.h"
#include "Teuchos_oblackholestream.hpp"
#include "Teuchos_RCP.hpp"
#include "Teuchos_BLAS.hpp"
#include "Teuchos_GlobalMPISession.hpp"
#include "Teuchos_XMLParameterListHelpers.hpp"
#include "Shards_CellTopology.hpp"
#include "EpetraExt_RowMatrixOut.h"
#include "EpetraExt_MultiVectorOut.h"
#include "create_inline_mesh.h"
#include "pamgen_im_exodusII_l.h"
#include "pamgen_im_ne_nemesisI_l.h"
#include "pamgen_extras.h"
#include "AztecOO.h"
#include "ml_MultiLevelPreconditioner.h"
#include "ml_epetra_utils.h"
#include "Sacado_No_Kokkos.hpp"
Include dependency graph for example_Poisson.cpp:

Typedefs

typedef Sacado::Fad::SFad
< double, 3 > 
Fad3
 
typedef
Intrepid::FunctionSpaceTools 
IntrepidFSTools
 
typedef
Intrepid::RealSpaceTools
< double > 
IntrepidRSTools
 
typedef Intrepid::CellTools
< double > 
IntrepidCTools
 

Functions

int TestMultiLevelPreconditionerLaplace (char ProblemType[], Teuchos::ParameterList &MLList, Epetra_CrsMatrix &A, const Epetra_MultiVector &xexact, Epetra_MultiVector &b, Epetra_MultiVector &uh, double &TotalErrorResidual, double &TotalErrorExactSol)
 ML Preconditioner. More...
 
template<typename Scalar >
const Scalar exactSolution (const Scalar &x, const Scalar &y, const Scalar &z)
 User-defined exact solution. More...
 
template<typename Scalar >
void materialTensor (Scalar material[][3], const Scalar &x, const Scalar &y, const Scalar &z)
 User-defined material tensor. More...
 
template<typename Scalar >
void exactSolutionGrad (Scalar gradExact[3], const Scalar &x, const Scalar &y, const Scalar &z)
 Computes gradient of the exact solution. Requires user-defined exact solution. More...
 
template<typename Scalar >
const Scalar sourceTerm (Scalar &x, Scalar &y, Scalar &z)
 Computes source term: f = -div(A.grad u). Requires user-defined exact solution and material tensor. More...
 
template<class ArrayOut , class ArrayIn >
void evaluateMaterialTensor (ArrayOut &worksetMaterialValues, const ArrayIn &evaluationPoints)
 Computation of the material tensor at array of points in physical space. More...
 
template<class ArrayOut , class ArrayIn >
void evaluateSourceTerm (ArrayOut &sourceTermValues, const ArrayIn &evaluationPoints)
 Computation of the source term at array of points in physical space. More...
 
template<class ArrayOut , class ArrayIn >
void evaluateExactSolution (ArrayOut &exactSolutionValues, const ArrayIn &evaluationPoints)
 Computation of the exact solution at array of points in physical space. More...
 
template<class ArrayOut , class ArrayIn >
void evaluateExactSolutionGrad (ArrayOut &exactSolutionGradValues, const ArrayIn &evaluationPoints)
 Computation of the gradient of the exact solution at array of points in physical space. More...
 
int main (int argc, char *argv[])
 

Detailed Description

Example solution of a Poisson equation on a hexahedral mesh using nodal (Hgrad) elements.

Example solution of a Poisson equation on a quad mesh using nodal (Hgrad) elements.

   This example uses the following Trilinos packages:
 Poisson system:

        div A grad u = f in Omega
                   u = g on Gamma

  where
         A is a symmetric, positive definite material tensor
         f is a given source term

 Corresponding discrete linear system for nodal coefficients(x):

             Kx = b

        K - HGrad stiffness matrix
        b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe
Example driver requires input file named Poisson.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.
   This example uses the following Trilinos packages:
 Poisson system:

        div A grad u = f in Omega
                   u = g on Gamma

  where
         A is a symmetric, positive definite material tensor
         f is a given source term

 Corresponding discrete linear system for nodal coefficients(x):

             Kx = b

        K - HGrad stiffness matrix
        b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe
Example driver requires input file named Poisson2D.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.

This example uses the following Trilinos packages:

Poisson system:

div A grad u = f in Omega
u = g on Gamma

where
A is a symmetric, positive definite material tensor
f is a given source term

Corresponding discrete linear system for nodal coefficients(x):

Kx = b

K - HGrad stiffness matrix
b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe
Example driver requires input file named Poisson2D.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.

This example uses the following Trilinos packages:

Poisson system:

div A grad u = f in Omega
u = g on Gamma

where
A is a symmetric, positive definite material tensor
f is a given source term

Corresponding discrete linear system for nodal coefficients(x):

Kx = b

K - HGrad stiffness matrix
b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe
Example driver requires input file named Poisson2D.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.

This example uses the following Trilinos packages:

Poisson system:

div A grad u = f in Omega
u = g on Gamma

where
A is a symmetric, positive definite material tensor
f is a given source term

Corresponding discrete linear system for nodal coefficients(x):

Kx = b

K - HGrad stiffness matrix
b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe
Example driver requires input file named Poisson2D.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.
   This example uses the following Trilinos packages:
 Poisson system:

        div A grad u = f in Omega
                   u = g on Gamma

  where
         A is a symmetric, positive definite material tensor
         f is a given source term

 Corresponding discrete linear system for nodal coefficients(x):

             Kx = b

        K - HGrad stiffness matrix
        b - right hand side vector
Author
Created by P. Bochev, D. Ridzal, K. Peterson, D. Hensinger, C. Siefert.
Remarks
Usage:
./TrilinosCouplings_examples_scaling_example_Poisson.exe [--meshfile filename]
Example driver requires input file named Poisson.xml with Pamgen formatted mesh description and settings for Isorropia (a version is included in the Trilinos repository with this driver).
The exact solution (u) and material tensor (A) are set in the functions "exactSolution" and "materialTensor" and may be modified by the user.

Function Documentation

template<class ArrayOut , class ArrayIn >
void evaluateExactSolution ( ArrayOut &  exactSolutionValues,
const ArrayIn &  evaluationPoints 
)

Computation of the exact solution at array of points in physical space.

Parameters
exactSolutionValues[out] Rank-2 (C,P) array with the values of the exact solution
evaluationPoints[in] Rank-3 (C,P,D) array with the evaluation points in physical frame
template<class ArrayOut , class ArrayIn >
void evaluateExactSolutionGrad ( ArrayOut &  exactSolutionGradValues,
const ArrayIn &  evaluationPoints 
)

Computation of the gradient of the exact solution at array of points in physical space.

Parameters
exactSolutionGradValues[out] Rank-3 (C,P,D) array with the values of the gradient of the exact solution
evaluationPoints[in] Rank-3 (C,P,D) array with the evaluation points in physical frame
template<class ArrayOut , class ArrayIn >
void evaluateMaterialTensor ( ArrayOut &  worksetMaterialValues,
const ArrayIn &  evaluationPoints 
)

Computation of the material tensor at array of points in physical space.

Parameters
worksetMaterialValues[out] Rank-2, 3 or 4 array with dimensions (C,P), (C,P,D) or (C,P,D,D) with the values of the material tensor
evaluationPoints[in] Rank-3 (C,P,D) array with the evaluation points in physical frame
template<class ArrayOut , class ArrayIn >
void evaluateSourceTerm ( ArrayOut &  sourceTermValues,
const ArrayIn &  evaluationPoints 
)

Computation of the source term at array of points in physical space.

Parameters
sourceTermValues[out] Rank-2 (C,P) array with the values of the source term
evaluationPoints[in] Rank-3 (C,P,D) array with the evaluation points in physical frame
template<typename Scalar >
const Scalar exactSolution ( const Scalar &  x,
const Scalar &  y,
const Scalar &  z 
)

User-defined exact solution.

Parameters
x[in] x-coordinate of the evaluation point
y[in] y-coordinate of the evaluation point
z[in] z-coordinate of the evaluation point
Returns
Value of the exact solution at (x,y,z)
template<typename Scalar >
void exactSolutionGrad ( Scalar  gradExact[3],
const Scalar &  x,
const Scalar &  y,
const Scalar &  z 
)

Computes gradient of the exact solution. Requires user-defined exact solution.

Parameters
gradExact[out] gradient of the exact solution evaluated at (x,y,z)
x[in] x-coordinate of the evaluation point
y[in] y-coordinate of the evaluation point
z[in] z-coordinate of the evaluation point

References exactSolution().

template<typename Scalar >
void materialTensor ( Scalar  material[][3],
const Scalar &  x,
const Scalar &  y,
const Scalar &  z 
)

User-defined material tensor.

Parameters
material[out] 3 x 3 material tensor evaluated at (x,y,z)
x[in] x-coordinate of the evaluation point
y[in] y-coordinate of the evaluation point
z[in] z-coordinate of the evaluation point
Warning
Symmetric and positive definite tensor is required for every (x,y,z).
template<typename Scalar >
const Scalar sourceTerm ( Scalar &  x,
Scalar &  y,
Scalar &  z 
)

Computes source term: f = -div(A.grad u). Requires user-defined exact solution and material tensor.

Parameters
x[in] x-coordinate of the evaluation point
y[in] y-coordinate of the evaluation point
z[in] z-coordinate of the evaluation point
Returns
Source term corresponding to the user-defined exact solution evaluated at (x,y,z)

References exactSolution(), and exactSolutionGrad().

int TestMultiLevelPreconditionerLaplace ( char  ProblemType[],
Teuchos::ParameterList &  MLList,
Epetra_CrsMatrix &  A,
const Epetra_MultiVector &  xexact,
Epetra_MultiVector &  b,
Epetra_MultiVector &  uh,
double &  TotalErrorResidual,
double &  TotalErrorExactSol 
)

ML Preconditioner.

Parameters
ProblemType[in] problem type
MLList[in] ML parameter list
A[in] discrete operator matrix
xexact[in] exact solution
b[in] right-hand-side vector
uh[out] solution vector
TotalErrorResidual[out] error residual
TotalErrorExactSol[out] error in uh