FractureMechanics::MWLSApprox Struct Reference

#include <users_modules/fracture_mechanics/src/MWLS.hpp>

Collaboration diagram for FractureMechanics::MWLSApprox:

Classes
struct	OpMWLSBase
struct	OpMWLSCalculateBaseCoeffcientsAtGaussPtsTmpl
struct	OpMWLSMaterialStressLhs_Dx
struct	OpMWLSMaterialStressLhs_DX
struct	OpMWLSMaterialStressRhs
	Integrate force vector. More...
struct	OpMWLSRhoAtGaussPts
	Evaluate density at integration points. More...
struct	OpMWLSRhoAtGaussUsingPrecalulatedCoeffs
struct	OpMWLSRhoPostProcess
struct	OpMWLSSpatialStressLhs_DX
struct	OpMWLSSpatialStressRhs
	Integrate force vector. More...
struct	OpMWLSStressAndErrorsAtGaussPts
	Evaluate stress at integration points. More...
struct	OpMWLSStressAtGaussPts
	Evaluate stress at integration points. More...
struct	OpMWLSStressAtGaussUsingPrecalulatedCoeffs
struct	OpMWLSStressPostProcess

Public Types
typedef ublas::vector< double, ublas::bounded_array< double, 9 > >	VecVals
typedef ublas::matrix< double, ublas::row_major, ublas::bounded_array< double, 27 > >	MatDiffVals
typedef ublas::matrix< double, ublas::row_major, ublas::bounded_array< double, 81 > >	MatDiffDiffVals
typedef OpMWLSCalculateBaseCoeffcientsAtGaussPtsTmpl< VolumeElementForcesAndSourcesCore >	OpMWLSCalculateBaseCoeffcientsAtGaussPts
typedef OpMWLSCalculateBaseCoeffcientsAtGaussPtsTmpl< FaceElementForcesAndSourcesCore >	OpMWLSCalculateBaseCoeffcientsForFaceAtGaussPts

Public Member Functions
	MWLSApprox (MoFEM::Interface &m_field, Vec F_lambda=PETSC_NULL, boost::shared_ptr< DofEntity > arc_length_dof=nullptr)
virtual	~MWLSApprox ()=default
bool	getUseGlobalBaseAtMaterialReferenceConfiguration () const
	Get the Use Local Base At Material Reference Configuration object. More...
MoFEMErrorCode	getMWLSOptions ()
	get options from command line for MWLS More...
MoFEMErrorCode	loadMWLSMesh (std::string file_name)
MoFEMErrorCode	getValuesToNodes (Tag th)
	get values form elements to nodes More...
MoFEMErrorCode	getInfluenceNodes (double material_coords[3])
	Get node in influence domain. More...
MoFEMErrorCode	calculateApproxCoeff (bool trim_influence_nodes, bool global_derivatives)
	Calculate approximation function coefficients. More...
MoFEMErrorCode	evaluateApproxFun (double eval_point[3])
	evaluate approximation function at material point More...
MoFEMErrorCode	evaluateFirstDiffApproxFun (double eval_point[3], bool global_derivatives)
	evaluate 1st derivative of approximation function at material point More...
MoFEMErrorCode	evaluateSecondDiffApproxFun (double eval_point[3], bool global_derivatives)
	evaluate 2nd derivative of approximation function at material point More...
MoFEMErrorCode	getTagData (Tag th)
MoFEMErrorCode	getErrorAtNode (Tag th, double &total_stress_at_node, double &total_stress_error_at_node, double &hydro_static_error_at_node, double &deviatoric_error_at_node, double &total_energy, double &total_energy_error, const double young_modulus, const double poisson_ratio)
	Get the norm of the error at nod. More...
const VecVals &	getDataApprox ()
const MatDiffVals &	getDiffDataApprox ()
const MatDiffDiffVals &	getDiffDiffDataApprox ()
MatrixDouble &	getInvAB ()
std::vector< EntityHandle > &	getInfluenceNodes ()
auto &	getShortenDm ()

Public Attributes
MoFEM::Interface &	mField
moab::Core	mwlsCore
moab::Interface &	mwlsMoab
Vec	F_lambda
boost::weak_ptr< DofEntity >	arcLengthDof
MatrixDouble	stressAtGaussPts
MatrixDouble	diffStressAtGaussPts
boost::shared_ptr< VectorDouble >	rhoAtGaussPts
boost::shared_ptr< MatrixDouble >	diffRhoAtGaussPts
boost::shared_ptr< MatrixDouble >	diffDiffRhoAtGaussPts
MatrixDouble	singularInitialDisplacement
MatrixDouble	singularCurrentDisplacement
MatrixDouble	mwlsMaterialPositions
MatrixDouble	approxBasePoint
VectorDouble3	approxPointCoords
Range	tetsInBlock
Range	trisInBlock
Tag	thMidNode
std::map< EntityHandle, std::vector< MatrixDouble > >	invABMap
	Store Coefficient of base functions at integration points. More...
std::map< EntityHandle, std::vector< std::vector< EntityHandle > > >	influenceNodesMap
	Influence nodes on elements on at integration points. More...
std::map< EntityHandle, std::vector< double > >	dmNodesMap
	Store weight function radius at the nodes. More...
boost::function< VectorDouble(const double, const double, const double)>	functionTestHack
	This is used to set up analytical function for testing approximation. More...

Private Member Functions
template<class T = FTensor::Tensor1<double, 3>>
MoFEMErrorCode	calculateBase (const T &t_coords)
template<class T = FTensor::Tensor1<double, 3>>
MoFEMErrorCode	calculateDiffBase (const T &t_coords, const double dm)
template<class T = FTensor::Tensor1<double, 3>>
MoFEMErrorCode	calculateDiffDiffBase (const T &t_coords, const double dm)
double	evalWeight (const double r)
double	evalDiffWeight (const double r)
double	evalDiffDiffWeight (const double r)

Private Attributes
boost::shared_ptr< BVHTree >	myTree
int	nbBasePolynomials
double	dmFactor
	Relative value of weight function radius. More...
double	shortenDm
int	maxElemsInDMFactor
int	maxThreeDepth
	Maximal three depths. More...
int	splitsPerDir
	Split per direction in the tree. More...
PetscBool	useGlobalBaseAtMaterialReferenceConfiguration
PetscBool	useNodalData
PetscBool	mwlsAnalyticalInternalStressTest
EntityHandle	treeRoot
double	maxEdgeL
Range	levelNodes
Range	levelEdges
Range	levelTets
std::vector< EntityHandle >	influenceNodes
std::map< EntityHandle, std::array< double, 3 > >	elemCentreMap
EntityHandle	nearestInfluenceNode
ublas::symmetric_matrix< double, ublas::lower >	A
ublas::symmetric_matrix< double, ublas::lower >	testA
MatrixDouble	invA
ublas::triangular_matrix< double, ublas::lower >	L
MatrixDouble	B
VectorDouble	baseFun
VectorDouble	approxFun
MatrixDouble	invAB
VectorDouble	baseFunInvA
ublas::symmetric_matrix< double, ublas::lower, ublas::row_major, DoubleAllocator >	diffA [3]
ublas::symmetric_matrix< double, ublas::lower, ublas::row_major, DoubleAllocator >	diffInvA [3]
MatrixDouble	diffB [3]
VectorDouble	diffBaseFun [3]
VectorDouble	diffDiffBaseFun [9]
MatrixDouble	diffApproxFun
MatrixDouble	diffDiffApproxFun
MatrixDouble	influenceNodesData
VecVals	outData
double	outHydroStaticError
VecVals	outDeviatoricDifference
double	outDeviatoricError
MatDiffVals	outDiffData
MatDiffDiffVals	outDiffDiffData
std::vector< EntityHandle >	treeTets
std::vector< EntityHandle >	leafsTetsVecLeaf
std::vector< double >	leafsTetsCentre
std::vector< EntityHandle >	leafsVec
std::vector< double >	leafsDist
std::vector< std::pair< double, EntityHandle > >	distLeafs

Detailed Description

Moving least weighted square approximation (MWLS)

\[ u^h(x) = p_k(x)a_k(x)\\ J(x)=\frac{1}{2} \sum_i w(x-x_i)\left(p_j(x_i)a_j(x)-u_i\right)^2 = \\ J(x)=\frac{1}{2} \sum_i w(x-x_i)\left( p_j(x_i)a_j(x) p_k(x_i)a_k(x) - 2p_j(x_i)a_j(x)u_i + u_i^2 \right) \\ \frac{\partial J}{\partial a_j} = \sum_i w(x_i-x)\left(p_j(x_i)p_k(x_i) a_k(x) -p_j(x_j)u_i\right)\\ A_{jk}(x) = \sum_i w(x_i-x)p_j(x_i)p_k(x_i)\\ B_{ji}(x) = w(x-x_i) p_j(x_i)\\ A_{jk}(x) a_k(x) = B_{ji}(x)u_i\\ a_k(x) = \left( A^{-1}_{kj}(x) B_{ji}(x) \right) u_i\\ u^h(x) = p_k(x) \left( A^{-1}_{kj}(x) B_{ji}(x) \right) u_i \\ u^h(x) = \left[ p_k(x) \left( A^{-1}_{kj}(x) B_{ji}(x) \right) \right] u_i \\ \phi_i = p_k(x) \left( A^{-1}_{kj}(x) B_{ji}(x) \right) \\ \phi_{i,l} = p_{k,l}(x) \left( A^{-1}_{kj}(x) B_{ji}(x) \right)+ p_k(x) \left( A^{-1}_{kj,l}(x) B_{ji}(x) + A^{-1}_{kj}(x) B_{ji,l}(x) \right)\\ \phi_{i,l} = p_{k,l}(x) \left( A^{-1}_{kj}(x) B_{ji}(x) \right)+ p_k(x) \left( A^{-1}_{kj,l}(x) B_{ji}(x) + A^{-1}_{kj}(x) B_{ji,l}(x) \right)\\ \left(\mathbf{A}\mathbf{A}^{-1}\right)_{,l} = \mathbf{0} \to \mathbf{A}^{-1}_{,l} = -\mathbf{A}^{-1} \mathbf{A}_{,l} \mathbf{A}^{-1}\\ A_{jk,l} = \sum_i w_{,l}(x_i-x)p_j(x_i)p_k(x_i) \\ B_{ij,l} = w_{,l}(x_i-x)p_j(x_i) \]

Definition at line 53 of file MWLS.hpp.

Member Typedef Documentation

MatDiffDiffVals

typedef ublas::matrix<double, ublas::row_major, ublas::bounded_array<double, 81> > FractureMechanics::MWLSApprox::MatDiffDiffVals

Definition at line 181 of file MWLS.hpp.

MatDiffVals

typedef ublas::matrix<double, ublas::row_major, ublas::bounded_array<double, 27> > FractureMechanics::MWLSApprox::MatDiffVals

Definition at line 178 of file MWLS.hpp.

OpMWLSCalculateBaseCoeffcientsAtGaussPts

typedef OpMWLSCalculateBaseCoeffcientsAtGaussPtsTmpl< VolumeElementForcesAndSourcesCore> FractureMechanics::MWLSApprox::OpMWLSCalculateBaseCoeffcientsAtGaussPts

Definition at line 296 of file MWLS.hpp.

OpMWLSCalculateBaseCoeffcientsForFaceAtGaussPts

typedef OpMWLSCalculateBaseCoeffcientsAtGaussPtsTmpl< FaceElementForcesAndSourcesCore> FractureMechanics::MWLSApprox::OpMWLSCalculateBaseCoeffcientsForFaceAtGaussPts

Definition at line 300 of file MWLS.hpp.

VecVals

typedef ublas::vector<double, ublas::bounded_array<double, 9> > FractureMechanics::MWLSApprox::VecVals

Definition at line 175 of file MWLS.hpp.

Constructor & Destructor Documentation

MWLSApprox()

FractureMechanics::MWLSApprox::MWLSApprox	(	MoFEM::Interface &	m_field,
		Vec	F_lambda = PETSC_NULL,
		boost::shared_ptr< DofEntity >	arc_length_dof = nullptr
	)

Definition at line 36 of file MWLS.cpp.

    : mField(m_field), mwlsMoab(mwlsCore), F_lambda(PETSC_NULL),
      arcLengthDof(arc_length_dof), nbBasePolynomials(1), dmFactor(1),
      maxElemsInDMFactor(2),
      useGlobalBaseAtMaterialReferenceConfiguration(PETSC_FALSE),
      useNodalData(PETSC_FALSE), mwlsAnalyticalInternalStressTest(PETSC_FALSE),
      approxPointCoords(3) {
 
  if (!LogManager::checkIfChannelExist("MWLSWorld")) {
    auto core_log = logging::core::get();
 
    core_log->add_sink(
        LogManager::createSink(LogManager::getStrmWorld(), "MWLSWorld"));
    core_log->add_sink(
        LogManager::createSink(LogManager::getStrmSync(), "MWLSSync"));
    core_log->add_sink(
        LogManager::createSink(LogManager::getStrmSelf(), "MWLSSelf"));
 
    LogManager::setLog("MWLSWorld");
    LogManager::setLog("MWLSSync");
    LogManager::setLog("MWLSSelf");
 
    MOFEM_LOG_TAG("MWLSWorld", "MWLS");
    MOFEM_LOG_TAG("MWLSSync", "MWLS");
    MOFEM_LOG_TAG("MWLSSelf", "MWLS");
  }
 
  MOFEM_LOG("MWLSWorld", Sev::noisy) << "MWLS created";
 
  rhoAtGaussPts = boost::make_shared<VectorDouble>();
  diffRhoAtGaussPts = boost::make_shared<MatrixDouble>();
  diffDiffRhoAtGaussPts = boost::make_shared<MatrixDouble>();
}

~MWLSApprox()

virtual FractureMechanics::MWLSApprox::~MWLSApprox ( )

virtualdefault

Member Function Documentation

calculateApproxCoeff()

MoFEMErrorCode FractureMechanics::MWLSApprox::calculateApproxCoeff	(	bool	trim_influence_nodes,
		bool	global_derivatives
	)

Calculate approximation function coefficients.

Parameters

material_coords

material position

Returns

error code

NOTE: This function remove form influenceNodes entities which are beyond influence radius.

Definition at line 576 of file MWLS.cpp.

                                                                         {
  MoFEMFunctionBegin;
 
  auto to_range = [&](auto ents) {
    Range r;
    r.insert_list(ents.begin(), ents.end());
    return r;
  };
 
  auto save_entities = [&](const std::string name, Range ents) {
    MoFEMFunctionBeginHot;
    if (!ents.empty()) {
      EntityHandle meshset;
      CHKERR mwlsMoab.create_meshset(MESHSET_SET, meshset);
      CHKERR mwlsMoab.add_entities(meshset, ents);
 
      CHKERR mwlsMoab.write_file(name.c_str(), "VTK", "", &meshset, 1);
      CHKERR mwlsMoab.delete_entities(&meshset, 1);
    }
    MoFEMFunctionReturnHot(0);
  };
 
  // Determine points in influence volume
  auto trim_nodes = [&](const double dm) {
    MoFEMFunctionBeginHot;
    FTensor::Index<'i', 3> i;
    FTensor::Tensor1<double, 3> t_approx_point(
        approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
    const double dm2 = shortenDm * shortenDm;
    VectorDouble node_coors_vec(3 * influenceNodes.size());
    CHKERR mwlsMoab.get_coords(&*influenceNodes.begin(), influenceNodes.size(),
                               &*node_coors_vec.begin());
    FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_node_coords(
        &node_coors_vec[0], &node_coors_vec[1], &node_coors_vec[2]);
    decltype(influenceNodes) influence_nodes_tmp;
    influence_nodes_tmp.reserve(influenceNodes.size());
    for (auto node : influenceNodes) {
      t_node_coords(i) -= t_approx_point(i);
      const double r2 = t_node_coords(i) * t_node_coords(i);
      if (r2 < dm2)
        influence_nodes_tmp.emplace_back(node);
      ++t_node_coords;
    }
    influenceNodes.swap(influence_nodes_tmp);
    MoFEMFunctionReturnHot(0);
  };
 
  auto eval_AB = [this, global_derivatives](const double dm) {
    MoFEMFunctionBegin;
    FTensor::Index<'i', 3> i;
 
    A.resize(nbBasePolynomials, nbBasePolynomials, false);
    B.resize(nbBasePolynomials, influenceNodes.size(), false);
    A.clear();
    B.clear();
    for (int d = 0; d != 3; d++) {
      diffA[d].resize(nbBasePolynomials, nbBasePolynomials, false);
      diffB[d].resize(nbBasePolynomials, influenceNodes.size(), false);
      diffA[d].clear();
      diffB[d].clear();
    }
 
    FTensor::Tensor1<double, 3> t_approx_point(
        approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
 
    VectorDouble nodes_coords(3 * influenceNodes.size());
    CHKERR mwlsMoab.get_coords(&*influenceNodes.begin(), influenceNodes.size(),
                               &*nodes_coords.begin());
    double *ptr = &*nodes_coords.begin();
    FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_node_coords(
        ptr, ptr + 1, ptr + 2);
 
    double min_distance = -1;
    nearestInfluenceNode = 0;
 
    for (int ii = 0; ii != influenceNodes.size(); ++ii) {
      t_node_coords(i) -= t_approx_point(i);
      t_node_coords(i) /= shortenDm;
      const double r = sqrt(t_node_coords(i) * t_node_coords(i));
 
      // Get the influence node that has minimum distance from the node in the
      // projected mesh
      if (!ii || r < min_distance) {
        min_distance = r;
        nearestInfluenceNode = influenceNodes[ii];
      }
 
      FTensor::Tensor1<double, 3> t_diff_r;
      if (r > 0) {
        t_diff_r(i) = (1. / r) * t_node_coords(i) * (-1 / shortenDm);
      } else {
        t_diff_r(i) = 0;
      }
      // Weights
      const double w = evalWeight(r);
      calculateBase<FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>>(
          t_node_coords);
 
      // iterate base functions
      // A_kj = sum_ii  w * b_k * b_j
      //
      // B_k = w * b_k
      //
      // diffA_kj = sum_ii diff_w * b_k * b_j + w * diff_b_k * b_j + w * b_k *
      // diff_b_j
      //
      // diffB =  diff_w * b_k + w * diff_bk
      for (int k = 0; k != nbBasePolynomials; k++) {
        const double wp = w * baseFun[k];
        for (int j = 0; j <= k; j++) {
          A(k, j) += wp * baseFun[j];
        }
        B(k, ii) = wp;
        if (global_derivatives) {
          const double diff_w_r = evalDiffWeight(r);
          const double diff_wp_r = diff_w_r * baseFun[k];
          for (int d = 0; d != 3; d++) {
            const double diff_wp_r_dx = diff_wp_r * t_diff_r(d);
            for (int j = 0; j <= k; j++) {
              diffA[d](k, j) += diff_wp_r_dx * baseFun[j];
            }
            diffB[d](k, ii) = diff_wp_r_dx;
          }
        }
      }
      ++t_node_coords;
    }
    if (nearestInfluenceNode == 0) {
      SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
              "Nearest node not found \n");
    }
 
    MoFEMFunctionReturn(0);
  };
 
  // Invert matrix A
  auto invert_A = [this]() {
    MoFEMFunctionBegin;
 
    auto solve = [&](const size_t nb_base_functions) {
      A.resize(nb_base_functions, nb_base_functions, true);
      L.resize(nb_base_functions, nb_base_functions, false);
      if (cholesky_decompose(A, L)) {
        return false;
      } else {
 
        MatrixDouble inv_a(nb_base_functions, nb_base_functions, false);
        for (int k = 0; k != nb_base_functions; k++) {
          ublas::matrix_column<MatrixDouble> mr(inv_a, k);
          mr(k) = 1;
          cholesky_solve(L, mr, ublas::lower());
        }
 
        if (nb_base_functions == nbBasePolynomials) {
          invA.swap(inv_a);
        } else {
          invA.resize(nbBasePolynomials, nbBasePolynomials, false);
          invA.clear();
          for (size_t i = 0; i != nb_base_functions; ++i)
            for (size_t j = 0; j != nb_base_functions; ++j)
              invA(i, j) = inv_a(i, j);
        }
 
        return true;
      }
    };
 
    auto throw_error = [&]() {
      MoFEMFunctionBegin;
      FTensor::Tensor1<double, 3> t_approx_point(
          approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
      cerr << "Point: " << t_approx_point << endl;
      cerr << "Matrix A: " << A << endl;
      SETERRQ1(mField.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
               "Failed to invert matrix - solution could be "
               "to increase dmFactor to bigger value, e.g. -mwls_dm 4, "
               "or reduce of number of base fuctions "
               "-mwls_number_of_base_functions 1. The number of nodes "
               "found in "
               "radius is %u",
               influenceNodes.size());
      MoFEMFunctionReturn(0);
    };
 
    auto solve_one = [&]() {
      MoFEMFunctionBegin;
      if (!solve(1))
        CHKERR throw_error();
      MoFEMFunctionReturn(0);
    };
 
    testA = A;
    if (!solve(nbBasePolynomials)) {
      if (nbBasePolynomials == 10) {
        if (!solve(4))
          CHKERR solve_one();
      } else if (nbBasePolynomials == 4) {
        CHKERR solve_one();
      } else {
        CHKERR throw_error();
      }
    }
 
    MoFEMFunctionReturn(0);
  };
 
  auto calculate_shape_functions_coefficients = [this]() {
    MoFEMFunctionBegin;
    // Calculate base functions (shape functions) coefficients
    invAB.resize(invA.size1(), B.size2(), false);
    noalias(invAB) = prod(invA, B);
    MoFEMFunctionReturn(0);
  };
 
  // CHKERR trim_nodes(shortenDm);
  // CHKERR save_entities("tree_tets_trim.vtk", to_range(influenceNodes));
 
  CHKERR eval_AB(shortenDm);
  CHKERR invert_A();
  CHKERR calculate_shape_functions_coefficients();
 
  MoFEMFunctionReturn(0);
}

calculateBase()

template<class T = FTensor::Tensor1<double, 3>>

MoFEMErrorCode FractureMechanics::MWLSApprox::calculateBase ( const T &  t_coords )

inlineprivate

Definition at line 675 of file MWLS.hpp.

                                                  {
    MoFEMFunctionBeginHot;
 
    baseFun.resize(nbBasePolynomials, false);
    baseFun.clear();
 
    baseFun[0] = 1;
    if (nbBasePolynomials > 1) {
      baseFun[1] = t_coords(0); // x
      baseFun[2] = t_coords(1); // y
      baseFun[3] = t_coords(2); // z
    }
 
    if (nbBasePolynomials > 4) {
      baseFun[4] = t_coords(0) * t_coords(0); // x^2
      baseFun[5] = t_coords(1) * t_coords(1); // y^2
      baseFun[6] = t_coords(2) * t_coords(2); // z^2
      baseFun[7] = t_coords(0) * t_coords(1); // xy
      baseFun[8] = t_coords(0) * t_coords(2); // xz
      baseFun[9] = t_coords(1) * t_coords(2); // yz
    }
 
    MoFEMFunctionReturnHot(0);
  }

calculateDiffBase()

template<class T = FTensor::Tensor1<double, 3>>

MoFEMErrorCode FractureMechanics::MWLSApprox::calculateDiffBase ( const T &  t_coords, const double  dm  )

inlineprivate

Definition at line 701 of file MWLS.hpp.

                                                                       {
    MoFEMFunctionBeginHot;
 
    // Local base functions
    for (int d = 0; d != 3; ++d) {
      diffBaseFun[d].resize(nbBasePolynomials, false);
      diffBaseFun[d].clear();
    }
 
    for (int d = 0; d != 3; d++) {
      diffBaseFun[d][0] = 0;
    }
 
    if (nbBasePolynomials > 1) {
      // derivative dx
      diffBaseFun[0][1] = 1. / dm;
      // derivative dy
      diffBaseFun[1][2] = 1. / dm;
      // derivative dz
      diffBaseFun[2][3] = 1. / dm;
    }
 
    if (nbBasePolynomials > 4) {
      // dx derivative
      diffBaseFun[0][4] = 2 * t_coords(0) / dm;
      diffBaseFun[0][5] = 0;
      diffBaseFun[0][6] = 0;
      diffBaseFun[0][7] = t_coords(1) / dm;
      diffBaseFun[0][8] = t_coords(2) / dm;
      diffBaseFun[0][9] = 0;
      // dy derivative
      diffBaseFun[1][4] = 0;
      diffBaseFun[1][5] = 2 * t_coords(1) / dm;
      diffBaseFun[1][6] = 0;
      diffBaseFun[1][7] = t_coords(0) / dm;
      diffBaseFun[1][8] = 0;
      diffBaseFun[1][9] = t_coords(2) / dm;
      // dz derivative
      diffBaseFun[2][4] = 0;
      diffBaseFun[2][5] = 0;
      diffBaseFun[2][6] = 2 * t_coords(2) / dm;
      diffBaseFun[2][7] = 0;
      diffBaseFun[2][8] = t_coords(0) / dm;
      diffBaseFun[2][9] = t_coords(1) / dm;
    }
    MoFEMFunctionReturnHot(0);
  }

calculateDiffDiffBase()

template<class T = FTensor::Tensor1<double, 3>>

MoFEMErrorCode FractureMechanics::MWLSApprox::calculateDiffDiffBase ( const T &  t_coords, const double  dm  )

inlineprivate

Definition at line 750 of file MWLS.hpp.

                                                                           {
    MoFEMFunctionBeginHot;
 
    // Local base functions
    for (int d = 0; d != 9; ++d) {
      diffDiffBaseFun[d].resize(nbBasePolynomials, false);
      diffDiffBaseFun[d].clear();
    }
 
    // if (nbBasePolynomials > 1) {
    // }
 
    if (nbBasePolynomials > 4) {
      const double dm2 = dm * dm;
      // set the non-zero values
      // dxd? derivative
      diffDiffBaseFun[0][4] = 2. / dm2; // dxdx
      diffDiffBaseFun[1][7] = 1. / dm2; // dxdy
      diffDiffBaseFun[2][8] = 1. / dm2; // dxdz
      // dyd? derivative
      diffDiffBaseFun[4][5] = 2. / dm2; // dydy
      diffDiffBaseFun[3][7] = 1. / dm2; // dydx
      diffDiffBaseFun[5][9] = 1. / dm2; // dydz
      // dzd? derivative
      diffDiffBaseFun[8][6] = 2. / dm2; // dzdz
      diffDiffBaseFun[6][8] = 1. / dm2; // dzdx
      diffDiffBaseFun[7][9] = 1. / dm2; // dzdy
    }
    MoFEMFunctionReturnHot(0);
  }

evalDiffDiffWeight()

double FractureMechanics::MWLSApprox::evalDiffDiffWeight ( const double  r )

inlineprivate

Definition at line 812 of file MWLS.hpp.

                                                   {
    const double rr = r * r;
    if (r < 1) {
      return -36 * rr + 48 * r - 12;
    } else {
      return 0;
    }
  }

evalDiffWeight()

double FractureMechanics::MWLSApprox::evalDiffWeight ( const double  r )

inlineprivate

Definition at line 797 of file MWLS.hpp.

                                               {
    const double rr = r * r;
    const double rrr = rr * r;
    if (r < 1) {
      return -12 * r + 24 * rr - 12 * rrr;
    } else {
      return 0;
    }
    // if(r<=0.5) {
    //   return -8*r+12*rr;
    // } else if(r<=1) {
    //   return -4+8*r-4*rr;
    // } else return 0;
  }

evaluateApproxFun()

MoFEMErrorCode FractureMechanics::MWLSApprox::evaluateApproxFun

(

double

eval_point[3]

)

evaluate approximation function at material point

Parameters

material_coords

material position

Returns

error code

Definition at line 801 of file MWLS.cpp.

                                                                 {
  MoFEMFunctionBegin;
  // Determine points in influence volume
 
  auto cal_base_functions_at_eval_point = [this, eval_point](const double dm) {
    MoFEMFunctionBegin;
 
    FTensor::Tensor1<double, 3> t_approx_point(
        approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
    FTensor::Tensor1<double, 3> t_eval_point(eval_point[0], eval_point[1],
                                             eval_point[2]);
 
    FTensor::Index<'i', 3> i;
    FTensor::Tensor1<double, 3> t_delta;
    t_delta(i) = t_eval_point(i) - t_approx_point(i);
    t_delta(i) /= dm;
 
    CHKERR calculateBase(t_delta);
    MoFEMFunctionReturn(0);
  };
 
  auto calculate_shape_functions = [this]() {
    MoFEMFunctionBegin;
    // Calculate approximation functions at eval point
    approxFun.resize(influenceNodes.size(), false);
    noalias(approxFun) = prod(baseFun, invAB);
    MoFEMFunctionReturn(0);
  };
 
  CHKERR cal_base_functions_at_eval_point(shortenDm);
  CHKERR calculate_shape_functions();
 
  MoFEMFunctionReturn(0);
}

evaluateFirstDiffApproxFun()

MoFEMErrorCode FractureMechanics::MWLSApprox::evaluateFirstDiffApproxFun	(	double	eval_point[3],
		bool	global_derivatives
	)

evaluate 1st derivative of approximation function at material point

Parameters

global_derivatives

decide whether global derivative will be calculated

Returns

error code

Definition at line 836 of file MWLS.cpp.

                                                                               {
  MoFEMFunctionBegin;
 
  auto cal_diff_base_functions_at_eval_point = [this,
                                                eval_point](const double dm) {
    MoFEMFunctionBegin;
 
    FTensor::Tensor1<double, 3> t_approx_point(
        approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
    FTensor::Tensor1<double, 3> t_eval_point(eval_point[0], eval_point[1],
                                             eval_point[2]);
 
    FTensor::Index<'i', 3> i;
    FTensor::Tensor1<double, 3> t_delta;
    t_delta(i) = t_eval_point(i) - t_approx_point(i);
    t_delta(i) /= dm;
 
    CHKERR calculateDiffBase(t_delta, dm);
    MoFEMFunctionReturn(0);
  };
 
  auto calculate_global_shape_functions_derivatives = [this]() {
    MoFEMFunctionBegin;
 
    // Calculate derivatives of base functions
    VectorDouble tmp;
    tmp.resize(nbBasePolynomials, false);
    baseFunInvA.resize(invA.size2(), false);
    noalias(baseFunInvA) = prod(baseFun, invA);
 
    for (int d = 0; d != 3; d++) {
      diffInvA[d].resize(nbBasePolynomials, nbBasePolynomials, false);
      diffApproxFun.resize(3, influenceNodes.size(), false);
    }
 
    for (int d = 0; d != 3; d++) {
      MatrixDouble a = prod(diffA[d], invA);
      noalias(diffInvA[d]) = -prod(invA, a);
    }
 
    // Calculate line derivatives for each coordinate direction
    for (int d = 0; d != 3; ++d) {
 
      ublas::matrix_row<MatrixDouble> mr(diffApproxFun, d);
 
      // b,x * A^-1 B
      noalias(mr) = prod(diffBaseFun[d], invAB);
 
      // b * (A^-1 A,x A^-1) * B
      noalias(tmp) = prod(baseFun, diffInvA[d]);
      mr += prod(tmp, B);
 
      // b * A^-1 * B,x
      mr += prod(baseFunInvA, diffB[d]);
    }
 
    MoFEMFunctionReturn(0);
  };
 
  auto calculate_local_shape_functions_derivatives = [this]() {
    MoFEMFunctionBegin;
 
    // Calculate derivatives of base functions
    diffApproxFun.resize(3, influenceNodes.size(), false);
    // Calculate line derivatives for each coordinate direction
    for (int d = 0; d != 3; ++d) {
      ublas::matrix_row<MatrixDouble> mr(diffApproxFun, d);
      // b,x * A^-1 B
      noalias(mr) = prod(diffBaseFun[d], invAB);
    }
 
    MoFEMFunctionReturn(0);
  };
 
  CHKERR cal_diff_base_functions_at_eval_point(shortenDm);
  if (global_derivatives)
    CHKERR calculate_global_shape_functions_derivatives();
  else
    CHKERR calculate_local_shape_functions_derivatives();
 
  MoFEMFunctionReturn(0);
}

evaluateSecondDiffApproxFun()

MoFEMErrorCode FractureMechanics::MWLSApprox::evaluateSecondDiffApproxFun	(	double	eval_point[3],
		bool	global_derivatives
	)

evaluate 2nd derivative of approximation function at material point

Parameters

global_derivatives

decide whether global derivative will be calculated

Returns

error code

Definition at line 921 of file MWLS.cpp.

                                                                 {
  MoFEMFunctionBegin;
 
  auto cal_diff_diff_base_functions_at_eval_point =
      [this, eval_point](const double dm) {
        MoFEMFunctionBegin;
 
        FTensor::Tensor1<double, 3> t_approx_point(
            approxPointCoords[0], approxPointCoords[1], approxPointCoords[2]);
        FTensor::Tensor1<double, 3> t_eval_point(eval_point[0], eval_point[1],
                                                 eval_point[2]);
 
        FTensor::Index<'i', 3> i;
        FTensor::Tensor1<double, 3> t_delta;
        t_delta(i) = t_eval_point(i) - t_approx_point(i);
        t_delta(i) /= dm;
 
        CHKERR calculateDiffDiffBase(t_delta, dm);
        MoFEMFunctionReturn(0);
      };
 
  auto calculate_global_shape_functions_second_derivatives = [this]() {
    MoFEMFunctionBegin;
    // to be implemented
    MoFEMFunctionReturn(0);
  };
 
  auto calculate_local_shape_functions_second_derivatives = [this]() {
    MoFEMFunctionBegin;
 
    // Calculate derivatives of base functions
    // for (int d = 0; d != 3; d++) {
    diffDiffApproxFun.resize(9, influenceNodes.size(), false);
    // }
    // Calculate line derivatives for each coordinate direction
    for (int n = 0; n != 3; ++n) {
      for (int d = 0; d != 3; ++d) {
        const int indx = d + 3 * n;
        ublas::matrix_row<MatrixDouble> mr(diffDiffApproxFun, indx);
        // b,xy * A^-1 B
        noalias(mr) = prod(diffDiffBaseFun[indx], invAB);
      }
    }
    MoFEMFunctionReturn(0);
  };
 
  CHKERR cal_diff_diff_base_functions_at_eval_point(shortenDm);
  if (global_derivatives) {
    SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented yet");
    // CHKERR calculate_global_shape_functions_second_derivatives();
  } else
    CHKERR calculate_local_shape_functions_second_derivatives();
 
  MoFEMFunctionReturn(0);
}

evalWeight()

double FractureMechanics::MWLSApprox::evalWeight ( const double  r )

inlineprivate

Definition at line 781 of file MWLS.hpp.

                                           {
    const double rr = r * r;
    const double rrr = rr * r;
    const double r4 = rrr * r;
    if (r < 1) {
      return 1 - 6 * rr + 8 * rrr - 3 * r4;
    } else {
      return 0;
    }
    // if(r<=0.5) {
    //   return 2./3-4*rr+4*rrr;
    // } else if(r<=1) {
    //   return 4./3-4*r+4*rr-(4./3)*rrr;
    // } else return 0;
  }

getDataApprox()

const MWLSApprox::VecVals & FractureMechanics::MWLSApprox::getDataApprox

(

)

Get values at point

Returns

error code

Definition at line 992 of file MWLS.cpp.

                                                   {
  noalias(outData) = prod(approxFun, influenceNodesData);
  return outData;
}

getDiffDataApprox()

const MWLSApprox::MatDiffVals & FractureMechanics::MWLSApprox::getDiffDataApprox

(

)

Get derivatives of values at points

Returns

error code

Definition at line 997 of file MWLS.cpp.

                                                           {
  noalias(outDiffData) = prod(diffApproxFun, influenceNodesData);
  return outDiffData;
}

getDiffDiffDataApprox()

const MWLSApprox::MatDiffDiffVals & FractureMechanics::MWLSApprox::getDiffDiffDataApprox

(

)

Get second derivatives of values at points

Returns

error code

Definition at line 1002 of file MWLS.cpp.

                                                                   {
  noalias(outDiffDiffData) = prod(diffDiffApproxFun, influenceNodesData);
  return outDiffDiffData;
}

getErrorAtNode()

MoFEMErrorCode FractureMechanics::MWLSApprox::getErrorAtNode	(	Tag	th,
		double &	total_stress_at_node,
		double &	total_stress_error_at_node,
		double &	hydro_static_error_at_node,
		double &	deviatoric_error_at_node,
		double &	total_energy,
		double &	total_energy_error,
		const double	young_modulus,
		const double	poisson_ratio
	)

Get the norm of the error at nod.

Parameters

th
total_stress_at_node
total_stress_error_at_node
hydro_static_error_at_node
deviatoric_error_at_node
total_energy
total_energy_error
young_modulus
poisson_ratio

Returns

MoFEMErrorCode

Definition at line 1007 of file MWLS.cpp.

                                                            {
  MoFEMFunctionBegin;
 
  int length;
  CHKERR mwlsMoab.tag_get_length(th, length);
 
  FTensor::Index<'i', 3> i;
  FTensor::Index<'j', 3> j;
  FTensor::Index<'k', 3> k;
  FTensor::Index<'l', 3> l;
 
  auto get_compliance_matrix = [&]() {
    FTensor::Ddg<double, 3, 3> t_C;
    t_C(i, j, k, l) = 0;
 
    const double c = poisson_ratio / young_modulus;
    t_C(0, 0, 0, 0) = 1. / young_modulus;
    t_C(0, 0, 1, 1) = -c;
    t_C(0, 0, 2, 2) = -c;
 
    t_C(1, 1, 0, 0) = -c;
    t_C(1, 1, 1, 1) = 1. / young_modulus;
    t_C(1, 1, 2, 2) = -c;
 
    t_C(2, 2, 0, 0) = -c;
    t_C(2, 2, 1, 1) = -c;
    t_C(2, 2, 2, 2) = 1. / young_modulus;
 
    const double d = (1 + poisson_ratio) / young_modulus;
    t_C(0, 1, 0, 1) = d;
    t_C(0, 2, 0, 2) = d;
    t_C(1, 2, 1, 2) = d;
 
    return t_C;
  };
 
  auto t_C = get_compliance_matrix();
 
  if (length != outData.size())
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "Tag length is not consistent with outData size. \nMost probably "
            "getTagData has not been invoked yet!\n");
 
  VectorDouble val_at_nearest_influence_node;
  val_at_nearest_influence_node.resize(length, false);
  val_at_nearest_influence_node.clear();
  CHKERR mwlsMoab.tag_get_data(th, &nearestInfluenceNode, 1,
                               &*val_at_nearest_influence_node.begin());
 
  // SIXX 0 (00) SIYY 1 (11) SIZZ 2 (22) SIXY 3 (01) SIXZ 4 (02) SIYZ 5 (12)
  // T* d00, T* d01, T* d02, T* d11, T* d12, T* d22
  auto get_symm_tensor = [](auto &m) {
    FTensor::Tensor2_symmetric<double *, 3> t(
 
        &m(0), &m(3), &m(4),
 
        &m(1), &m(5),
 
        &m(2)
 
    );
    return t;
  };
 
  auto t_total_stress_error = FTensor::Tensor2_symmetric<double, 3>();
  auto t_stress_at_nearset_node =
      get_symm_tensor(val_at_nearest_influence_node);
  auto t_mwls_stress = get_symm_tensor(outData);
 
  auto get_energy_norm = [&](auto &t_s) {
    return 0.5 * t_s(i, j) * (t_C(i, j, k, l) * t_s(k, l));
  };
  auto get_stress_nrm2 = [i, j](auto &t_s) { return t_s(i, j) * t_s(i, j); };
 
  total_stress_at_node = get_stress_nrm2(t_stress_at_nearset_node);
  total_energy = get_energy_norm(t_stress_at_nearset_node);
 
  t_total_stress_error(i, j) =
      t_mwls_stress(i, j) - t_stress_at_nearset_node(i, j);
  total_stress_error_at_node = get_stress_nrm2(t_total_stress_error);
  total_energy_error = get_energy_norm(t_total_stress_error);
 
  constexpr auto t_kd = FTensor::Kronecker_Delta_symmetric<int>();
  // calculate error of the hydro static part of the stress
  double hydro_static_difference =
      (t_kd(i, j) * t_total_stress_error(i, j)) / 3.;
  hydro_static_error_at_node =
      hydro_static_difference * hydro_static_difference;
 
  // calculate error of the deviator part of the stress
  auto t_dev_error = FTensor::Tensor2_symmetric<double, 3>();
  t_dev_error(i, j) =
      t_total_stress_error(i, j) - hydro_static_difference * t_kd(i, j);
  deviatoric_error_at_node = get_stress_nrm2(t_dev_error);
 
  MoFEMFunctionReturn(0);
}

getInfluenceNodes() [1/2]

std::vector<EntityHandle>& FractureMechanics::MWLSApprox::getInfluenceNodes ( )

inline

Definition at line 605 of file MWLS.hpp.

                                                      {
    return influenceNodes;
  }

getInfluenceNodes() [2/2]

MoFEMErrorCode FractureMechanics::MWLSApprox::getInfluenceNodes

(

double

material_coords[3]

)

Get node in influence domain.

Parameters

material_coords

material position

Returns

error code

Definition at line 373 of file MWLS.cpp.

                                                                      {
  constexpr double eps_overlap = 0.01;
  const size_t max_nb_elems = maxElemsInDMFactor * nbBasePolynomials + 1;
  MoFEMFunctionBegin;
 
    using DistEntPair = std::pair<double, EntityHandle>;
    using DistEntMap = multi_index_container<
        DistEntPair,
        indexed_by<
 
            ordered_non_unique<member<DistEntPair, double, &DistEntPair::first>>
 
            >>;
    DistEntMap dist_map;
 
    // FIXME: That is potentisl source of problems, should be sparate
    // function for setting approximation base point, or simply should be set
    // outside of that function. Other functions will behave depending how
    // approxPointCoords is set. And here is somehow hidden from the view,
    // i.e. easy to overlook this. TODO: create a new function like:
    // MoFEMErrorCode setApproxBasePt(double material_coords[3]) ?
 
    for (int dd : {0, 1, 2})
      approxPointCoords[dd] = material_coords[dd];
 
    if (myTree) {
      shortenDm = maxEdgeL * dmFactor; // Influence radius
 
      auto find_leafs = [&](const auto dm) {
        leafsVec.clear();
        leafsDist.clear();
        CHKERR myTree->distance_search(material_coords, shortenDm, leafsVec,
                                       1.0e-10, 1.0e-6, &leafsDist);
        distLeafs.clear();
        distLeafs.reserve(leafsDist.size());
        for (size_t i = 0; i != leafsVec.size(); ++i) {
          distLeafs.emplace_back(leafsDist[i], leafsVec[i]);
        }
        std::sort(
            distLeafs.begin(), distLeafs.end(),
            [](const auto &a, const auto &b) { return a.first < b.first; });
        return distLeafs;
      };
 
      auto get_influence_nodes = [&](const auto &tets) {
        std::vector<EntityHandle> influence_nodes;
        if (!useNodalData) {
          influence_nodes.resize(tets.size());
          CHKERR mwlsMoab.tag_get_data(thMidNode, &*tets.begin(), tets.size(),
                                       &*influence_nodes.begin());
        } else {
          Range nodes;
          CHKERR mwlsMoab.get_connectivity(&*tets.begin(), tets.size(), nodes,
                                           true);
          influence_nodes.clear();
          influence_nodes.reserve(nodes.size());
          for (auto n : nodes)
            influence_nodes.emplace_back(n);
        }
        return influence_nodes;
      };
 
      auto to_range = [&](auto ents) {
        Range r;
        r.insert_list(ents.begin(), ents.end());
        return r;
      };
 
      auto save_entities = [&](const std::string name, Range ents) {
        MoFEMFunctionBeginHot;
        if (!ents.empty()) {
          EntityHandle meshset;
          CHKERR mwlsMoab.create_meshset(MESHSET_SET, meshset);
          CHKERR mwlsMoab.add_entities(meshset, ents);
 
          CHKERR mwlsMoab.write_file(name.c_str(), "VTK", "", &meshset, 1);
          CHKERR mwlsMoab.delete_entities(&meshset, 1);
        }
        MoFEMFunctionReturnHot(0);
      };
 
      auto leafs = find_leafs(shortenDm);
 
      treeTets.clear();
      std::vector<EntityHandle> leafs_tets;
      moab::BoundBox box;
 
      auto get_dm2_from_influence_points = [&]() {
        FTensor::Index<'i', 3> i;
        VectorDouble node_coors_vec(3 * influenceNodes.size());
        CHKERR mwlsMoab.get_coords(&*influenceNodes.begin(),
                                   influenceNodes.size(),
                                   &*node_coors_vec.begin());
        FTensor::Tensor1<double, 3> t_approx_point(
            material_coords[0], material_coords[1], material_coords[2]);
        FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_node_coords(
            &node_coors_vec[0], &node_coors_vec[1], &node_coors_vec[2]);
        double dm2 = 0;
        for (auto node : influenceNodes) {
          t_node_coords(i) -= t_approx_point(i);
          const double r2 = t_node_coords(i) * t_node_coords(i);
          if (dm2 < r2)
            dm2 = r2;
          ++t_node_coords;
        }
        return dm2;
      };
 
      const double shorten_dm2 = shortenDm * shortenDm;
      double set_dm2 = shorten_dm2;
      auto min_max_dist_it = dist_map.get<0>().begin();
      FTensor::Index<'i', 3> i;
 
      for (auto l : leafs) {
        if (dist_map.size() < max_nb_elems ||
            ((l.first * l.first) <
             (set_dm2 + std::numeric_limits<double>::epsilon()))) {
 
          leafsTetsVecLeaf.clear();
          CHKERR mwlsMoab.get_entities_by_dimension(l.second, 3,
                                                    leafsTetsVecLeaf, false);
 
          leafsTetsCentre.resize(3 * leafsTetsVecLeaf.size());
          {
            FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 3> t_c{
                &leafsTetsCentre[0], &leafsTetsCentre[leafsTetsVecLeaf.size()],
                &leafsTetsCentre[2 * leafsTetsVecLeaf.size()]};
 
            for (auto tet : leafsTetsVecLeaf) {
              const auto &c = elemCentreMap[tet];
              for (auto ii : {0, 1, 2}) {
                t_c(ii) = c[ii] - material_coords[ii];
              }
              ++t_c;
            }
          }
 
          {
            FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 3> t_c{
                &leafsTetsCentre[0], &leafsTetsCentre[leafsTetsVecLeaf.size()],
                &leafsTetsCentre[2 * leafsTetsVecLeaf.size()]};
 
            for (auto tet : leafsTetsVecLeaf) {
 
              if (dist_map.size() > max_nb_elems)
                set_dm2 = std::min(min_max_dist_it->first, shorten_dm2);
 
              if (t_c(0) < set_dm2 && t_c(1) < set_dm2 && t_c(2) < set_dm2) {
                double dm2 = t_c(i) * t_c(i);
                if (dm2 < set_dm2) {
                  min_max_dist_it =
                      dist_map.get<0>().emplace(std::make_pair(dm2, tet)).first;
                  if (dist_map.size() > max_nb_elems) {
                    int nb = std::distance(dist_map.get<0>().begin(),
                                           min_max_dist_it);
                    for (; nb < max_nb_elems; ++nb)
                      ++min_max_dist_it;
                  }
                }
              }
 
              ++t_c;
            }
          }
        } else {
          break;
        }
      }
 
      treeTets.clear();
      treeTets.reserve(dist_map.size());
      const double dm2 = set_dm2 + std::numeric_limits<double>::epsilon();
      auto hi_it = dist_map.get<0>().upper_bound(dm2);
      for (auto it = dist_map.get<0>().begin(); it != hi_it; ++it)
        treeTets.emplace_back(it->second);
 
      if (treeTets.empty()) {
        for (auto &it : dist_map)
          MOFEM_LOG("MWLSWorld", Sev::error)
              << "dm map " << it.first << " " << it.second;
 
        MOFEM_LOG("MWLSWorld", Sev::error) << "leafs found " << leafs.size();
        MOFEM_LOG("MWLSWorld", Sev::error)
            << "dist map size " << dist_map.size();
        MOFEM_LOG("MWLSWorld", Sev::error) << "dm2 " << dm2;
        MOFEM_LOG("MWLSWorld", Sev::error) << "shorten dm is " << shortenDm;
        MOFEM_LOG_C("MWLSWorld", Sev::error, "point:  ( %g %g %g )",
                    approxPointCoords(0), approxPointCoords(1),
                    approxPointCoords(2));
        SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "No leafs found.");
      }
 
      influenceNodes = get_influence_nodes(treeTets);
      const auto dm2_from_influence_pts = get_dm2_from_influence_points();
      shortenDm = (1 + eps_overlap) * sqrt(dm2_from_influence_pts);
 
      // CHKERR save_entities("tree_tets.vtk", to_range(influenceNodes));
  } else {
    SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY, "kd-three not build");
  }
  MoFEMFunctionReturn(0);
}

getInvAB()

MatrixDouble& FractureMechanics::MWLSApprox::getInvAB ( )

inline

Definition at line 604 of file MWLS.hpp.

{ return invAB; }

getMWLSOptions()

MoFEMErrorCode FractureMechanics::MWLSApprox::getMWLSOptions

(

)

get options from command line for MWLS

Returns

error code

Definition at line 71 of file MWLS.cpp.

                                          {
  MoFEMFunctionBegin;
  ierr = PetscOptionsBegin(mField.get_comm(), "",
                           "Get options for moving least square approximation",
                           "none");
  CHKERRQ(ierr);
  CHKERR PetscOptionsScalar("-mwls_dm",
                            "edge length scaling for influence radius", "",
                            dmFactor, &dmFactor, PETSC_NULL);
  MOFEM_LOG_C("MWLSWorld", Sev::inform, "### Input parameter: -mwls_dm %6.4e",
              dmFactor);
 
  CHKERR PetscOptionsInt(
      "-mwls_number_of_base_functions",
      "1 = constant, 4 = linear, 10 = quadratic approximation", "",
      nbBasePolynomials, &nbBasePolynomials, PETSC_NULL);
  MOFEM_LOG_C("MWLSWorld", Sev::inform,
              "### Input parameter: -mwls_number_of_base_functions %d",
              nbBasePolynomials);
 
  CHKERR PetscOptionsInt("-mwls_max_elems_factor",
                         "Set max elements factor max_nb_elems = "
                         "maxElemsInDMFactor * nbBasePolynomials",
                         "", maxElemsInDMFactor, &maxElemsInDMFactor,
                         PETSC_NULL);
  MOFEM_LOG_C("MWLSWorld", Sev::inform,
              "### Input parameter: -mwls_max_elems_factor %d",
              maxElemsInDMFactor);
 
  CHKERR PetscOptionsBool(
      "-mwls_use_global_base_at_reference_configuration",
      "true local mwls base functions at reference configuration are used", "",
      useGlobalBaseAtMaterialReferenceConfiguration,
      &useGlobalBaseAtMaterialReferenceConfiguration, NULL);
  CHKERRQ(ierr);
 
  CHKERR PetscOptionsBool(
      "-mwls_use_nodal_data",
      "true nodal data values are used (without averaging to the midpoint)", "",
      useNodalData, &useNodalData, NULL);
  MOFEM_LOG_C("MWLSWorld", Sev::verbose,
              "### Input parameter: -mwls_use_nodal_data %d", useNodalData);
 
  CHKERR PetscOptionsBool(
      "-mwls_analytical_internal_stress",
      "use analytical strain field for internal stress mwls test", "",
      mwlsAnalyticalInternalStressTest, &mwlsAnalyticalInternalStressTest,
      NULL);
  MOFEM_LOG_C("MWLSWorld", Sev::verbose,
              "### Input parameter: -test_mwls_internal_stress_analytical %d",
              mwlsAnalyticalInternalStressTest);
 
  maxThreeDepth = 30;
  CHKERR PetscOptionsInt("-mwls_tree_depth", "maximal three depths", "",
                         maxThreeDepth, &maxThreeDepth, PETSC_NULL);
  MOFEM_LOG("MWLSWorld", Sev::verbose)
      << "Maximal MWLS three depths " << maxThreeDepth;
  splitsPerDir = 12;
  CHKERR PetscOptionsInt("-mwls_splits_per_dir", "splits per direction", "",
                         splitsPerDir, &splitsPerDir, PETSC_NULL);
  MOFEM_LOG("MWLSWorld", Sev::verbose)
      << "Splits per direction " << splitsPerDir;
 
  ierr = PetscOptionsEnd();
  CHKERRQ(ierr);
  MoFEMFunctionReturn(0);
}

getShortenDm()

auto& FractureMechanics::MWLSApprox::getShortenDm ( )

inline

Definition at line 608 of file MWLS.hpp.

{ return shortenDm; }

getTagData()

MoFEMErrorCode FractureMechanics::MWLSApprox::getTagData

(

Tag

th

)

Get tag data on influence nodes

Parameters

th

tag

Returns

error code

Definition at line 978 of file MWLS.cpp.

                                            {
  MoFEMFunctionBegin;
  int length;
  CHKERR mwlsMoab.tag_get_length(th, length);
  influenceNodesData.resize(influenceNodes.size(), length, false);
  CHKERR mwlsMoab.tag_get_data(th, &*influenceNodes.begin(),
                               influenceNodes.size(),
                               &influenceNodesData(0, 0));
  outData.resize(length, false);
  outDiffData.resize(3, length, false);
  outDiffDiffData.resize(9, length, false);
  MoFEMFunctionReturn(0);
}

getUseGlobalBaseAtMaterialReferenceConfiguration()

bool FractureMechanics::MWLSApprox::getUseGlobalBaseAtMaterialReferenceConfiguration ( ) const

inline

Get the Use Local Base At Material Reference Configuration object.

If true base is calculated at initail integration point coordinates and base function evaluated at another point. Local derivatives are uses, since weight function is at initial point.

Returns

true

false

Definition at line 76 of file MWLS.hpp.

                                                                {
    return (useGlobalBaseAtMaterialReferenceConfiguration == PETSC_TRUE);
  }

getValuesToNodes()

MoFEMErrorCode FractureMechanics::MWLSApprox::getValuesToNodes

(

Tag

th

)

get values form elements to nodes

Parameters

th	tag with approximated data
th	tag with approximated data

Returns

error code

Definition at line 275 of file MWLS.cpp.

                                                  {
  MoFEMFunctionBegin;
  if (useNodalData)
    MoFEMFunctionReturnHot(0);
 
  CHKERR mwlsMoab.tag_get_handle("MID_NODE", 1, MB_TYPE_HANDLE, thMidNode,
                                 MB_TAG_CREAT | MB_TAG_DENSE);
 
  int length;
  CHKERR mwlsMoab.tag_get_length(th, length);
 
  VectorDouble vals(length);
  vals.clear();
 
  if (mwlsAnalyticalInternalStressTest == PETSC_FALSE) {
 
    string name;
    CHKERR mwlsMoab.tag_get_name(th, name);
 
    std::function<MoFEMErrorCode(EntityHandle)> get_data_to_node;
    if (name.compare("RHO") == 0) {
      get_data_to_node = [this, th, &vals](EntityHandle tet) {
        MoFEMFunctionBeginHot;
        Range nodes;
        double avg_val = 0;
        CHKERR mwlsMoab.get_connectivity(&tet, 1, nodes, true);
        for (auto &n : nodes) {
          double tag_data;
          CHKERR mwlsMoab.tag_get_data(th, &n, 1, &tag_data);
          avg_val += tag_data;
        }
        avg_val /= nodes.size();
        vals[0] = avg_val;
        MoFEMFunctionReturnHot(0);
      };
    } else {
      get_data_to_node = [this, th, &vals](EntityHandle tet) {
        MoFEMFunctionBeginHot;
        CHKERR mwlsMoab.tag_get_data(th, &tet, 1, &*vals.begin());
        MoFEMFunctionReturnHot(0);
      };
    }
 
    // clean tag on nodes
    CHKERR mwlsMoab.tag_clear_data(th, levelNodes, &*vals.begin());
    for (auto tet : levelTets) {
 
      CHKERR get_data_to_node(tet);
      int num_nodes;
      const EntityHandle *conn;
      CHKERR mwlsMoab.get_connectivity(tet, conn, num_nodes, false);
      EntityHandle ho_node;
      EntityType tet_type = mwlsMoab.type_from_handle(tet);
      CHKERR mwlsMoab.high_order_node(tet, conn, tet_type, ho_node);
      CHKERR mwlsMoab.tag_set_data(th, &ho_node, 1, &*vals.begin());
      CHKERR mwlsMoab.tag_set_data(thMidNode, &tet, 1, &ho_node);
    }
 
  } else {
 
    PetscBool add_analytical_internal_stress_operators = PETSC_FALSE;
    CHKERR PetscOptionsGetBool(PETSC_NULL, "",
                               "-add_analytical_internal_stress_operators",
                               &add_analytical_internal_stress_operators, NULL);
 
    CHKERR mwlsMoab.tag_clear_data(th, levelNodes, &*vals.begin());
    for (auto tet : levelTets) {
      int num_nodes;
      const EntityHandle *conn;
      CHKERR mwlsMoab.get_connectivity(tet, conn, num_nodes, false);
      EntityHandle ho_node;
      EntityType tet_type = mwlsMoab.type_from_handle(tet);
      CHKERR mwlsMoab.high_order_node(tet, conn, tet_type, ho_node);
      CHKERR mwlsMoab.tag_set_data(thMidNode, &tet, 1, &ho_node);
 
      if (add_analytical_internal_stress_operators == PETSC_FALSE) {
        VectorDouble coords(3);
        CHKERR mwlsMoab.get_coords(&ho_node, 1, &coords[0]);
        FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 3> t_coords(
            &coords[0], &coords[1], &coords[2]);
        auto t_strain = CrackPropagation::analyticalStrainFunction(t_coords);
 
        constexpr double K = 166.666666666e9; // K = E / 3 / (1 - 2*nu),
                                              // E = 200 GPa, nu = 0.3
 
        // Only diagonal, vectorial notation
        vals[0] = 3 * K * t_strain(0, 0);
        vals[1] = 3 * K * t_strain(1, 1);
        vals[2] = 3 * K * t_strain(2, 2);
      }
 
      CHKERR mwlsMoab.tag_set_data(th, &ho_node, 1, &*vals.begin());
    }
  }
 
  MoFEMFunctionReturn(0);
}

loadMWLSMesh()

MoFEMErrorCode FractureMechanics::MWLSApprox::loadMWLSMesh

(

std::string

file_name

)

Load mesh with data and do basic calculation

Parameters

file_name

File name

Returns

Error code

Definition at line 139 of file MWLS.cpp.

                                                           {
  MoFEMFunctionBegin;
 
  CHKERR getMWLSOptions();
 
  const char *option;
  option = "";
  CHKERR mwlsMoab.load_file(file_name.c_str(), 0, option);
  CHKERR mwlsMoab.get_entities_by_dimension(0, 3, levelTets);
  if (levelTets.empty()) {
    SETERRQ(mField.get_comm(), MOFEM_DATA_INCONSISTENCY,
            "No 3d element in MWLS mesh");
  }
 
  CHKERR mwlsMoab.get_adjacencies(levelTets, 1, true, levelEdges,
                                  moab::Interface::UNION);
 
  // Create HO node on tetrahedral
  if (!useNodalData) {
    EntityHandle meshset;
    CHKERR mwlsMoab.create_meshset(MESHSET_SET, meshset);
    CHKERR mwlsMoab.add_entities(meshset, levelTets);
    CHKERR mwlsMoab.convert_entities(meshset, false, false, true);
    CHKERR mwlsMoab.delete_entities(&meshset, 1);
 
    // Get HO nodes in TET center
    std::vector<double> edge_length;
    for (auto tet : levelTets) {
      int num_nodes;
      const EntityHandle *conn;
      CHKERR mwlsMoab.get_connectivity(tet, conn, num_nodes, false);
      EntityHandle ho_node;
      EntityType tit_type = mwlsMoab.type_from_handle(tet);
      CHKERR mwlsMoab.high_order_node(tet, conn, tit_type, ho_node);
      levelNodes.insert(ho_node);
    }
 
  } else {
    Range tets;
    CHKERR mwlsMoab.get_entities_by_dimension(0, 3, tets);
    CHKERR mwlsMoab.get_connectivity(tets, levelNodes, true);
  }
 
  // Store edge nodes coordinates in FTensor
  double edge_node_coords[6];
  FTensor::Tensor1<double *, 3> t_node_edge[2] = {
      FTensor::Tensor1<double *, 3>(edge_node_coords, &edge_node_coords[1],
                                    &edge_node_coords[2]),
      FTensor::Tensor1<double *, 3>(&edge_node_coords[3], &edge_node_coords[4],
                                    &edge_node_coords[5])};
 
  FTensor::Index<'i', 3> i;
 
  // Get edge lengths
  maxEdgeL = 0;
  for (auto edge : levelEdges) {
    int num_nodes;
    const EntityHandle *conn;
    CHKERR mwlsMoab.get_connectivity(edge, conn, num_nodes, true);
    CHKERR mwlsMoab.get_coords(conn, num_nodes, edge_node_coords);
    t_node_edge[0](i) -= t_node_edge[1](i);
    double l = sqrt(t_node_edge[0](i) * t_node_edge[0](i));
    maxEdgeL = (maxEdgeL < l) ? l : maxEdgeL;
  }
 
  // Create tree and find maximal edge length. Maximal edge length is used
  // to find neighbours nodes to material point.
  myTree = boost::make_shared<BVHTree>(&mwlsMoab);
  {
    std::ostringstream opts;
    opts << "MAX_DEPTH=" << maxThreeDepth
        //  << ";MAX_PER_LEAF=" << maxElemsInDMFactor * nbBasePolynomials + 1
         << ";SPLITS_PER_DIR=" << splitsPerDir;
    FileOptions tree_opts(opts.str().c_str());
 
    MOFEM_LOG_CHANNEL("MWLSWorld");
    MOFEM_LOG_TAG("MWLSWorld", "MWLS");
    BOOST_LOG_SCOPED_THREAD_ATTR("Timeline", attrs::timer());
    MOFEM_LOG("MWLSWorld", Sev::verbose) << "Build tree ... ";
    CHKERR myTree->build_tree(levelTets, &treeRoot, &tree_opts);
    MOFEM_LOG("MWLSWorld", Sev::verbose) << "done";
  }
 
  // moab::BoundBox box;
  // CHKERR myTree->get_bounding_box(box);
  // cout << box << endl;
  // myTree->print();
 
  auto map_elems_positions = [&]() {
    std::map<EntityHandle, std::array<double, 3>> elem_coords_map;
    Range tets;
    CHKERR mwlsMoab.get_entities_by_dimension(0, 3, tets, false);
    MatrixDouble elem_coords;
    for (auto p = tets.pair_begin(); p != tets.pair_end(); ++p) {
      EntityHandle *conn;
      int count, verts_per_e;
      Range slot(p->first, p->second);
      CHKERR mwlsMoab.connect_iterate(slot.begin(), slot.end(), conn,
                                      verts_per_e, count);
      if (count != (int)slot.size())
        THROW_MESSAGE("Should be one continuos sequence");
      elem_coords.resize(verts_per_e, 3, false);
 
      for (auto t = 0; t != count; ++t) {
        CHKERR mwlsMoab.get_coords(conn, verts_per_e,
                                   &*elem_coords.data().begin());
        auto get_center = [&]() {
          std::array<double, 3> c{0, 0, 0};
          for (auto d : {0, 1, 2})
            for (auto n = 0; n < verts_per_e; ++n)
              c[d] += elem_coords(n, d);
          for (auto d : {0, 1, 2})
            c[d] /= verts_per_e;
          return c;
        };
        elem_coords_map.emplace(std::make_pair(p->first + t, get_center()));
        conn += verts_per_e;
      }
    }
    return elem_coords_map;
  };
 
  elemCentreMap = map_elems_positions();
 
  // Set tag for MWLS stress when SSLV116 test is run
  if (mwlsAnalyticalInternalStressTest) {
    Tag th;
    std::array<double, 9> def;
    def.fill(0);
    CHKERR mwlsMoab.tag_get_handle("MED_SSLV116", 9, MB_TYPE_DOUBLE, th,
                                   MB_TAG_CREAT | MB_TAG_SPARSE, def.data());
  }
 
  MoFEMFunctionReturn(0);
}

Member Data Documentation

A

ublas::symmetric_matrix<double, ublas::lower> FractureMechanics::MWLSApprox::A

private

Definition at line 641 of file MWLS.hpp.

approxBasePoint

MatrixDouble FractureMechanics::MWLSApprox::approxBasePoint

Definition at line 213 of file MWLS.hpp.

approxFun

VectorDouble FractureMechanics::MWLSApprox::approxFun

private

Definition at line 647 of file MWLS.hpp.

approxPointCoords

VectorDouble3 FractureMechanics::MWLSApprox::approxPointCoords

Definition at line 214 of file MWLS.hpp.

arcLengthDof

boost::weak_ptr<DofEntity> FractureMechanics::MWLSApprox::arcLengthDof

Definition at line 60 of file MWLS.hpp.

B

MatrixDouble FractureMechanics::MWLSApprox::B

private

Definition at line 645 of file MWLS.hpp.

baseFun

VectorDouble FractureMechanics::MWLSApprox::baseFun

private

Definition at line 646 of file MWLS.hpp.

baseFunInvA

VectorDouble FractureMechanics::MWLSApprox::baseFunInvA

private

Definition at line 649 of file MWLS.hpp.

diffA

ublas::symmetric_matrix<double, ublas::lower, ublas::row_major, DoubleAllocator> FractureMechanics::MWLSApprox::diffA[3]

private

Definition at line 653 of file MWLS.hpp.

diffApproxFun

MatrixDouble FractureMechanics::MWLSApprox::diffApproxFun

private

Definition at line 660 of file MWLS.hpp.

diffB

MatrixDouble FractureMechanics::MWLSApprox::diffB[3]

private

Definition at line 657 of file MWLS.hpp.

diffBaseFun

VectorDouble FractureMechanics::MWLSApprox::diffBaseFun[3]

private

Definition at line 658 of file MWLS.hpp.

diffDiffApproxFun

MatrixDouble FractureMechanics::MWLSApprox::diffDiffApproxFun

private

Definition at line 661 of file MWLS.hpp.

diffDiffBaseFun

VectorDouble FractureMechanics::MWLSApprox::diffDiffBaseFun[9]

private

Definition at line 659 of file MWLS.hpp.

diffDiffRhoAtGaussPts

boost::shared_ptr<MatrixDouble> FractureMechanics::MWLSApprox::diffDiffRhoAtGaussPts

Definition at line 208 of file MWLS.hpp.

diffInvA

ublas::symmetric_matrix<double, ublas::lower, ublas::row_major, DoubleAllocator> FractureMechanics::MWLSApprox::diffInvA[3]

private

Definition at line 656 of file MWLS.hpp.

diffRhoAtGaussPts

boost::shared_ptr<MatrixDouble> FractureMechanics::MWLSApprox::diffRhoAtGaussPts

Definition at line 207 of file MWLS.hpp.

diffStressAtGaussPts

MatrixDouble FractureMechanics::MWLSApprox::diffStressAtGaussPts

Definition at line 203 of file MWLS.hpp.

distLeafs

std::vector<std::pair<double, EntityHandle> > FractureMechanics::MWLSApprox::distLeafs

private

Definition at line 826 of file MWLS.hpp.

dmFactor

double FractureMechanics::MWLSApprox::dmFactor

private

Relative value of weight function radius.

Definition at line 617 of file MWLS.hpp.

dmNodesMap

std::map<EntityHandle, std::vector<double> > FractureMechanics::MWLSApprox::dmNodesMap

Store weight function radius at the nodes.

Definition at line 236 of file MWLS.hpp.

elemCentreMap

std::map<EntityHandle, std::array<double, 3> > FractureMechanics::MWLSApprox::elemCentreMap

private

Definition at line 637 of file MWLS.hpp.

F_lambda

Vec FractureMechanics::MWLSApprox::F_lambda

Definition at line 58 of file MWLS.hpp.

functionTestHack

boost::function<VectorDouble(const double, const double, const double)> FractureMechanics::MWLSApprox::functionTestHack

This is used to set up analytical function for testing approximation.

Definition at line 602 of file MWLS.hpp.

influenceNodes

std::vector<EntityHandle> FractureMechanics::MWLSApprox::influenceNodes

private

Definition at line 636 of file MWLS.hpp.

influenceNodesData

MatrixDouble FractureMechanics::MWLSApprox::influenceNodesData

private

Definition at line 663 of file MWLS.hpp.

influenceNodesMap

std::map<EntityHandle, std::vector<std::vector<EntityHandle> > > FractureMechanics::MWLSApprox::influenceNodesMap

Influence nodes on elements on at integration points.

Definition at line 230 of file MWLS.hpp.

invA

MatrixDouble FractureMechanics::MWLSApprox::invA

private

Definition at line 643 of file MWLS.hpp.

invAB

MatrixDouble FractureMechanics::MWLSApprox::invAB

private

Definition at line 648 of file MWLS.hpp.

invABMap

std::map<EntityHandle, std::vector<MatrixDouble> > FractureMechanics::MWLSApprox::invABMap

Store Coefficient of base functions at integration points.

Definition at line 223 of file MWLS.hpp.

L

ublas::triangular_matrix<double, ublas::lower> FractureMechanics::MWLSApprox::L

private

Definition at line 644 of file MWLS.hpp.

leafsDist

std::vector<double> FractureMechanics::MWLSApprox::leafsDist

private

Definition at line 825 of file MWLS.hpp.

leafsTetsCentre

std::vector<double> FractureMechanics::MWLSApprox::leafsTetsCentre

private

Definition at line 823 of file MWLS.hpp.

leafsTetsVecLeaf

std::vector<EntityHandle> FractureMechanics::MWLSApprox::leafsTetsVecLeaf

private

Definition at line 822 of file MWLS.hpp.

leafsVec

std::vector<EntityHandle> FractureMechanics::MWLSApprox::leafsVec

private

Definition at line 824 of file MWLS.hpp.

levelEdges

Range FractureMechanics::MWLSApprox::levelEdges

private

Definition at line 635 of file MWLS.hpp.

levelNodes

Range FractureMechanics::MWLSApprox::levelNodes

private

Definition at line 635 of file MWLS.hpp.

levelTets

Range FractureMechanics::MWLSApprox::levelTets

private

Definition at line 635 of file MWLS.hpp.

maxEdgeL

double FractureMechanics::MWLSApprox::maxEdgeL

private

Definition at line 633 of file MWLS.hpp.

maxElemsInDMFactor

int FractureMechanics::MWLSApprox::maxElemsInDMFactor

private

Maximal number of elements in the horizon is calculated from max_nb_elems = maxElemsInDMFactor *nbBasePolynomials. Value is set by command line -mwls_max_elems_factor

Definition at line 620 of file MWLS.hpp.

maxThreeDepth

int FractureMechanics::MWLSApprox::maxThreeDepth

private

Maximal three depths.

Definition at line 625 of file MWLS.hpp.

mField

MoFEM::Interface& FractureMechanics::MWLSApprox::mField

Definition at line 55 of file MWLS.hpp.

mwlsAnalyticalInternalStressTest

PetscBool FractureMechanics::MWLSApprox::mwlsAnalyticalInternalStressTest

private

Definition at line 630 of file MWLS.hpp.

mwlsCore

moab::Core FractureMechanics::MWLSApprox::mwlsCore

Definition at line 56 of file MWLS.hpp.

mwlsMaterialPositions

MatrixDouble FractureMechanics::MWLSApprox::mwlsMaterialPositions

Definition at line 212 of file MWLS.hpp.

mwlsMoab

moab::Interface& FractureMechanics::MWLSApprox::mwlsMoab

Definition at line 57 of file MWLS.hpp.

myTree

boost::shared_ptr<BVHTree> FractureMechanics::MWLSApprox::myTree

private

Definition at line 611 of file MWLS.hpp.

nbBasePolynomials

int FractureMechanics::MWLSApprox::nbBasePolynomials

private

Number of base functions (1 - linear, 4 - linear, 10 - quadratic)

Definition at line 614 of file MWLS.hpp.

nearestInfluenceNode

EntityHandle FractureMechanics::MWLSApprox::nearestInfluenceNode

private

Definition at line 639 of file MWLS.hpp.

outData

VecVals FractureMechanics::MWLSApprox::outData

private

Definition at line 664 of file MWLS.hpp.

outDeviatoricDifference

VecVals FractureMechanics::MWLSApprox::outDeviatoricDifference

private

Definition at line 668 of file MWLS.hpp.

outDeviatoricError

double FractureMechanics::MWLSApprox::outDeviatoricError

private

Definition at line 669 of file MWLS.hpp.

outDiffData

MatDiffVals FractureMechanics::MWLSApprox::outDiffData

private

Definition at line 671 of file MWLS.hpp.

outDiffDiffData

MatDiffDiffVals FractureMechanics::MWLSApprox::outDiffDiffData

private

Definition at line 672 of file MWLS.hpp.

outHydroStaticError

double FractureMechanics::MWLSApprox::outHydroStaticError

private

Definition at line 666 of file MWLS.hpp.

rhoAtGaussPts

boost::shared_ptr<VectorDouble> FractureMechanics::MWLSApprox::rhoAtGaussPts

Definition at line 206 of file MWLS.hpp.

shortenDm

double FractureMechanics::MWLSApprox::shortenDm

private

Radius of weight function used to calculate base function.

Definition at line 618 of file MWLS.hpp.

singularCurrentDisplacement

MatrixDouble FractureMechanics::MWLSApprox::singularCurrentDisplacement

Definition at line 211 of file MWLS.hpp.

singularInitialDisplacement

MatrixDouble FractureMechanics::MWLSApprox::singularInitialDisplacement

Definition at line 210 of file MWLS.hpp.

splitsPerDir

int FractureMechanics::MWLSApprox::splitsPerDir

private

Split per direction in the tree.

Definition at line 626 of file MWLS.hpp.

stressAtGaussPts

MatrixDouble FractureMechanics::MWLSApprox::stressAtGaussPts

Definition at line 202 of file MWLS.hpp.

testA

ublas::symmetric_matrix<double, ublas::lower> FractureMechanics::MWLSApprox::testA

private

Definition at line 642 of file MWLS.hpp.

tetsInBlock

Range FractureMechanics::MWLSApprox::tetsInBlock

Definition at line 215 of file MWLS.hpp.

thMidNode

Tag FractureMechanics::MWLSApprox::thMidNode

Definition at line 217 of file MWLS.hpp.

treeRoot

EntityHandle FractureMechanics::MWLSApprox::treeRoot

private

Definition at line 632 of file MWLS.hpp.

treeTets

std::vector<EntityHandle> FractureMechanics::MWLSApprox::treeTets

private

Definition at line 821 of file MWLS.hpp.

trisInBlock

Range FractureMechanics::MWLSApprox::trisInBlock

Definition at line 216 of file MWLS.hpp.

useGlobalBaseAtMaterialReferenceConfiguration

PetscBool FractureMechanics::MWLSApprox::useGlobalBaseAtMaterialReferenceConfiguration

private

Definition at line 628 of file MWLS.hpp.

useNodalData

PetscBool FractureMechanics::MWLSApprox::useNodalData

private

Definition at line 629 of file MWLS.hpp.

---

The documentation for this struct was generated from the following files:

• users_modules/fracture_mechanics/src/MWLS.hpp
• users_modules/fracture_mechanics/src/impl/MWLS.cpp