Example solution grad-div diffusion system with div-conforming (face) elements. More...

#include "TrilinosCouplings_config.h"
#include "TrilinosCouplings_Pamgen_Utils.hpp"
#include "Intrepid_FunctionSpaceTools.hpp"
#include "Intrepid_FieldContainer.hpp"
#include "Intrepid_CellTools.hpp"
#include "Intrepid_ArrayTools.hpp"
#include "Intrepid_HCURL_HEX_I1_FEM.hpp"
#include "Intrepid_HGRAD_HEX_C1_FEM.hpp"
#include "Intrepid_HDIV_HEX_I1_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_Vector.h"
#include "Epetra_LinearProblem.h"
#include "Epetra_Import.h"
#include "Epetra_Export.h"
#include "Teuchos_oblackholestream.hpp"
#include "Teuchos_RCP.hpp"
#include "Teuchos_BLAS.hpp"
#include "Teuchos_GlobalMPISession.hpp"
#include "Teuchos_ParameterList.hpp"
#include "Teuchos_XMLParameterListHelpers.hpp"
#include "Shards_CellTopology.hpp"
#include "EpetraExt_RowMatrixOut.h"
#include "EpetraExt_MultiVectorOut.h"
#include "EpetraExt_VectorOut.h"
#include "EpetraExt_MatrixMatrix.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_epetra_utils.h"
#include "ml_GradDiv.h"

Include dependency graph for example_GradDiv.cpp:

Classes
struct	fecomp

Macros
#define	ABS(x) ((x)>0?(x):-(x))

Typedefs
typedef Intrepid::FunctionSpaceTools	IntrepidFSTools

typedef Intrepid::RealSpaceTools< double >	IntrepidRSTools

typedef Intrepid::CellTools< double >	IntrepidCTools

Functions
int	Multiply_Abs (const Epetra_CrsMatrix &A, const Epetra_Vector &x, Epetra_Vector &y)
	Multiplies abs(A)x = y, where all non-zero entries of A are replaced with their absolute values value.

void	TestMultiLevelPreconditioner_GradDiv (char ProblemType[], Teuchos::ParameterList &MLList, Epetra_CrsMatrix &GradDiv, Epetra_CrsMatrix &D0clean, Epetra_CrsMatrix &D1clean, Epetra_CrsMatrix &FaceNode, Epetra_CrsMatrix &M1, Epetra_CrsMatrix &M1inv, Epetra_CrsMatrix &M2, Epetra_MultiVector &xh, Epetra_MultiVector &b, double &TotalErrorResidual, double &TotalErrorExactSol)
	ML Preconditioner.

int	evalu (double &uExact0, double &uExact1, double &uExact2, double &x, double &y, double &z)
	Exact solution evaluation.

double	evalDivu (double &x, double &y, double &z)
	Divergence of exact solution.

int	evalGradDivu (double &gradDivu0, double &gradDivu1, double &gradDivu2, double &x, double &y, double &z, double &kappa)
	Grad of Div of exact solution.

int	main (int argc, char *argv[])

int	Multiply_Ones (const Epetra_CrsMatrix &A, const Epetra_Vector &x, Epetra_Vector &y)
	Multiplies Ax = y, where all non-zero entries of A are replaced with the value 1.0.

void	solution_test (string msg, const Epetra_Operator &A, const Epetra_MultiVector &lhs, const Epetra_MultiVector &rhs, const Epetra_MultiVector &xexact, Epetra_Time &Time, double &TotalErrorExactSol, double &TotalErrorResidual)
	Compute ML solution residual.

int	evalCurlu (double &curlu0, double &curlu1, double &curlu2, double &x, double &y, double &z, double &mu)

Detailed Description

Example solution grad-div diffusion system with div-conforming (face) elements.

This example uses the following Trilinos packages:

Pamgen to generate a Hexahedral mesh.
Intrepid to build the discretization matrices and right-hand side.
Epetra to handle the global matrix and vector.
ML to solve the linear system.

For more details on the formulation see

P. Bochev, K. Peterson, C. Siefert, "Analysis and Computation of Compatible Least-Squares Methods for Div-Curl Equations", SIAM J. Numerical Analysis, vol 49, pp 159-191, 2011.

Grad-Div System:

grad kappa div F + mu F  = g  in Omega

Corresponding discrete linear system for edge element coeficients (x):
(Kc + Mc)x = b

Kd    - Hdiv stiffness matrix
Md    - Hcdiv mass 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_GradDiv.exe  inputfile.xml

inputfile.xml (optional)  -  xml input file containing Pamgen mesh description
and material parameters for each Pamgen block,
if not present code attempts to read DivLSFEMin.xml.

Function Documentation

◆ evalDivu()

double evalDivu	(	double &	x,
		double &	y,
		double &	z
	)

Divergence of exact solution.

Parameters

x	[in] x coordinate
y	[in] y coordinate
z	[in] z coordinate

Returns: Value of the divergence of exact solution at (x,y,z)

◆ evalGradDivu()

int evalGradDivu	(	double &	gradDivu0,
		double &	gradDivu1,
		double &	gradDivu2,
		double &	x,
		double &	y,
		double &	z,
		double &	kappa
	)

Grad of Div of exact solution.

Parameters

gradDivu0	[out] first component of grad div of exact solution
gradDivu1	[out] second component of grad div of exact solution
gradDivu2	[out] third component of grad div of exact solution
x	[in] x coordinate
y	[in] y coordinate
z	[in] z coordinate
mu	[in] material parameter

◆ evalu()

int evalu	(	double &	uExact0,
		double &	uExact1,
		double &	uExact2,
		double &	x,
		double &	y,
		double &	z
	)

Exact solution evaluation.

Parameters

uExact0	[out] first component of exact solution at (x,y,z)
uExact1	[out] second component of exact solution at (x,y,z)
uExact2	[out] third component of exact solution at (x,y,z)
x	[in] x coordinate
y	[in] y coordinate
z	[in] z coordinate

◆ Multiply_Abs()

int Multiply_Abs	(	const Epetra_CrsMatrix &	A,
		const Epetra_Vector &	x,
		Epetra_Vector &	y
	)

Multiplies abs(A)x = y, where all non-zero entries of A are replaced with their absolute values value.

Parameters

A	[in] matrix
x	[in] vector
y	[out] vector

◆ Multiply_Ones()

int Multiply_Ones	(	const Epetra_CrsMatrix &	A,
		const Epetra_Vector &	x,
		Epetra_Vector &	y
	)

Multiplies Ax = y, where all non-zero entries of A are replaced with the value 1.0.

Parameters

A	[in] matrix
x	[in] vector
y	[in] vector

◆ solution_test()

void solution_test	(	string	msg,
		const Epetra_Operator &	A,
		const Epetra_MultiVector &	lhs,
		const Epetra_MultiVector &	rhs,
		const Epetra_MultiVector &	xexact,
		Epetra_Time &	Time,
		double &	TotalErrorExactSol,
		double &	TotalErrorResidual
	)

Compute ML solution residual.

Parameters

A	[in] discrete operator
lhs	[in] solution vector
rhs	[in] right hand side vector
Time	[in] elapsed time for output
TotalErrorResidual	[out] error residual
TotalErrorExactSol	[out] error in xh (not an appropriate measure for H(div) basis functions)

Referenced by TestMultiLevelPreconditioner_GradDiv().

◆ TestMultiLevelPreconditioner_GradDiv()

void TestMultiLevelPreconditioner_GradDiv	(	char	ProblemType[],
		Teuchos::ParameterList &	MLList,
		Epetra_CrsMatrix &	GradDiv,
		Epetra_CrsMatrix &	D0clean,
		Epetra_CrsMatrix &	D1clean,
		Epetra_CrsMatrix &	FaceNode,
		Epetra_CrsMatrix &	M1,
		Epetra_CrsMatrix &	M1inv,
		Epetra_CrsMatrix &	M2,
		Epetra_MultiVector &	xh,
		Epetra_MultiVector &	b,
		double &	TotalErrorResidual,
		double &	TotalErrorExactSol
	)

ML Preconditioner.

Parameters

ProblemType	[in] problem type
MLList	[in] ML parameter list
GradDiv	[in] H(curl) stiffness matrix
D0clean	[in] Edge to node incidence matrix
D1clean	[in] Face to edge incidence matrix
FaceNode	[in] Face to node incidence matrix
M1inv	[in] H(curl) mass matrix inverse
M2	[in] H(div) mass matrix
xh	[out] solution vector
b	[in] right-hand-side vector
TotalErrorResidual	[out] error residual
TotalErrorExactSol	[out] error in xh

References solution_test().

Classes

Macros

Typedefs

Functions

Detailed Description

Function Documentation

◆ evalDivu()

◆ evalGradDivu()

◆ evalu()

◆ Multiply_Abs()

◆ Multiply_Ones()

◆ solution_test()

◆ TestMultiLevelPreconditioner_GradDiv()