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

#include <unistd.h>
#include <vector>
#include <map>
#include "TrilinosCouplings_config.h"
#include "TrilinosCouplings_TpetraIntrepidPoissonExample.hpp"
#include "TrilinosCouplings_IntrepidPoissonExampleHelpers.hpp"
#include <Teuchos_CommandLineProcessor.hpp>
#include <Teuchos_DefaultMpiComm.hpp>
#include <Teuchos_TimeMonitor.hpp>
#include "Intrepid_FunctionSpaceTools.hpp"
#include "Intrepid_CellTools.hpp"
#include "Intrepid_ArrayTools.hpp"
#include "Intrepid_Basis.hpp"
#include "Intrepid_HGRAD_HEX_C1_FEM.hpp"
#include "Intrepid_HGRAD_TET_C1_FEM.hpp"
#include "Intrepid_HGRAD_QUAD_C1_FEM.hpp"
#include "Intrepid_HGRAD_TRI_C1_FEM.hpp"
#include "Intrepid_RealSpaceTools.hpp"
#include "Intrepid_DefaultCubatureFactory.hpp"
#include "Intrepid_Utils.hpp"
#include "Tpetra_Map.hpp"
#include "Tpetra_FECrsMatrix.hpp"
#include "Tpetra_FECrsGraph.hpp"
#include "Tpetra_FEMultiVector.hpp"
#include "Tpetra_Import.hpp"
#include "Tpetra_Assembly_Helpers.hpp"
#include "MatrixMarket_Tpetra.hpp"
#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 "MueLu.hpp"
#include "MueLu_ParameterListInterpreter.hpp"
#include "MueLu_CreateTpetraPreconditioner.hpp"
#include "Sacado_No_Kokkos.hpp"
#include "Ionit_Initializer.h"
#include "Ioss_SubSystem.h"
#include "stk_io/IossBridge.hpp"
#include "stk_io/StkMeshIoBroker.hpp"
#include "stk_util/parallel/Parallel.hpp"
#include "stk_mesh/base/Types.hpp"
#include "stk_mesh/base/MetaData.hpp"
#include "stk_mesh/base/BulkData.hpp"
#include "stk_mesh/base/ForEachEntity.hpp"
#include "stk_mesh/base/GetEntities.hpp"
#include "stk_mesh/base/Selector.hpp"
#include "stk_mesh/base/Bucket.hpp"
#include "stk_mesh/base/Entity.hpp"
#include "stk_mesh/base/Field.hpp"
#include "TrilinosCouplings_Statistics.hpp"

Typedefs
typedef shards::CellTopology	ShardsCellTopology

typedef Intrepid::FunctionSpaceTools	IntrepidFSTools

typedef Intrepid::RealSpaceTools< double >	IntrepidRSTools

typedef Intrepid::CellTools< double >	IntrepidCTools

typedef Intrepid::FieldContainer< double >	IntrepidFieldContainer

using	Tpetra_CrsMatrix = Tpetra::CrsMatrix<>

using	Tpetra_MultiVector = Tpetra::MultiVector<>

using	Tpetra_Vector = Tpetra::Vector<>

Functions
template<typename Scalar >
const Scalar	exactSolution (const Scalar &x, const Scalar &y, const Scalar &z)
	User-defined exact solution.

template<typename Scalar >
void	materialTensor (Scalar material[][3], const Scalar &x, const Scalar &y, const Scalar &z)
	User-defined material tensor.

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.

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.

template<class ArrayOut , class ArrayIn >
void	evaluateMaterialTensor (ArrayOut &worksetMaterialValues, const ArrayIn &evaluationPoints)
	Computation of the material tensor at array of points in physical space.

template<class ArrayOut , class ArrayIn >
void	evaluateSourceTerm (ArrayOut &sourceTermValues, const ArrayIn &evaluationPoints)
	Computation of the source term at array of points in physical space.

template<class ArrayOut , class ArrayIn >
void	evaluateExactSolution (ArrayOut &exactSolutionValues, const ArrayIn &evaluationPoints)
	Computation of the exact solution at array of points in physical space.

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.

int	TestMultiLevelPreconditioner (char ProblemType[], Teuchos::ParameterList &MLList, Teuchos::RCP< Tpetra_CrsMatrix > &A, Teuchos::RCP< Tpetra_Vector > &xexact, Teuchos::RCP< Tpetra_MultiVector > &b, Teuchos::RCP< Tpetra_MultiVector > &uh, Teuchos::RCP< Tpetra_MultiVector > &coords, double &TotalErrorResidual, double &TotalErrorExactSol)

void	getBasis (Teuchos::RCP< Intrepid::Basis< double, IntrepidFieldContainer > > &basis, const ShardsCellTopology &cellTopology, int order)
	Simple factory that chooses basis function based on cell topology.

int	getDimension (const ShardsCellTopology &cellTopology)

void	mesh_read_write (const std::string &type, const std::string &working_directory, const std::string &filename, stk::io::StkMeshIoBroker &broker, int db_integer_size, stk::io::HeartbeatType hb_type)

template<class crs_matrix_type , class multivector_type1 , class multivector_type2 , class solution_type >
void	Apply_Dirichlet_BCs (std::vector< int > &BCNodes, crs_matrix_type &A, multivector_type1 &x, multivector_type2 &b, solution_type &solution_values)

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

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

Variables
int	spaceDim

Detailed Description

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

This example requires a hexahedral or tetrahedral mesh in Exodus format with a nodeset containing boundary nodes. STK is used to read the mesh and populate a mesh database, Intrepid is used to build the stiffness matrix and right-hand side, and ML is used to solve the resulting linear system.

 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 C. Siefert.

Remarks: Usage:
./example_Poisson_stk --help; Example requires a hexahedral or tetrahedral mesh in Exodus format

Function Documentation

◆ evaluateExactSolution()

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

◆ evaluateExactSolutionGrad()

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

◆ evaluateMaterialTensor()

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

◆ evaluateSourceTerm()

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

◆ exactSolution()

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)

Referenced by exactSolutionGrad(), and sourceTerm().

◆ exactSolutionGrad()

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().

Referenced by sourceTerm().

◆ getBasis()

void getBasis	(	Teuchos::RCP< Intrepid::Basis< double, IntrepidFieldContainer > > &	basis,
		const ShardsCellTopology &	cellTopology,
		int	order
	)

Simple factory that chooses basis function based on cell topology.

Parameters

cellTopology	[in] Shards cell topology
order	[in] basis function order, currently unused
basis	[out] pointer to Intrepid basis

Returns: Intrepid basis

◆ materialTensor()

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).

◆ sourceTerm()

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().

Typedefs

Functions

Variables

Detailed Description

Function Documentation

◆ evaluateExactSolution()

◆ evaluateExactSolutionGrad()

◆ evaluateMaterialTensor()

◆ evaluateSourceTerm()

◆ exactSolution()

◆ exactSolutionGrad()

◆ getBasis()

◆ materialTensor()

◆ sourceTerm()