v0.13.1
Classes | Typedefs | Functions | Variables
free_surface.cpp File Reference
#include <MoFEM.hpp>
#include <BasicFiniteElements.hpp>
#include <FreeSurfaceOps.hpp>

Go to the source code of this file.

Classes

struct  ElementsAndOps< DIM >
 
struct  ElementsAndOps< 2 >
 
struct  ExtractParentType< PARENT >
 
struct  FreeSurface
 
struct  Monitor
 Monitor solution. More...
 

Typedefs

using DomainEle = ElementsAndOps< SPACE_DIM >::DomainEle
 
using DomianParentEle = ElementsAndOps< SPACE_DIM >::DomianParentEle
 
using DomainEleOp = DomainEle::UserDataOperator
 
using BoundaryEle = ElementsAndOps< SPACE_DIM >::BoundaryEle
 
using BoundaryEleOp = BoundaryEle::UserDataOperator
 
using BoundaryParentEle = ElementsAndOps< SPACE_DIM >::BoundaryParentEle
 
using EntData = EntitiesFieldData::EntData
 
using PostProcEle = PostProcBrokenMeshInMoab< DomainEle >
 
using PostProcEdgeEle = PostProcBrokenMeshInMoab< BoundaryEle >
 
using AssemblyDomainEleOp = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::OpBase
 
using AssemblyBoundaryEleOp = FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::OpBase
 
using OpDomainMassU = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< BASE_DIM, U_FIELD_DIM >
 
using OpDomainMassH = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< BASE_DIM, 1 >
 
using OpDomainMassP = OpDomainMassH
 
using OpDomainSourceU = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< BASE_DIM, U_FIELD_DIM >
 
using OpDomainSourceH = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< BASE_DIM, 1 >
 
using OpBaseTimesScalar = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1, 1 >
 
using OpMixScalarTimesDiv = FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMixScalarTimesDiv< SPACE_DIM, coord_type >
 

Functions

template<int T>
constexpr int powof2 ()
 
int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 
constexpr int BASE_DIM = 1
 
constexpr int SPACE_DIM = 2
 
constexpr int U_FIELD_DIM = SPACE_DIM
 
constexpr CoordinateTypes coord_type
 select coordinate system <CARTESIAN, CYLINDRICAL>; More...
 
FTensor::Index< 'i', SPACE_DIMi
 
FTensor::Index< 'j', SPACE_DIMj
 
FTensor::Index< 'k', SPACE_DIMk
 
FTensor::Index< 'l', SPACE_DIMl
 
constexpr auto t_kd = FTensor::Kronecker_Delta_symmetric<int>()
 
constexpr int order = 3
 approximation order More...
 
constexpr double a0 = 0.98
 
constexpr double rho_m = 0.998
 
constexpr double mu_m = 0.0101
 
constexpr double rho_p = 0.0012
 
constexpr double mu_p = 0.000182
 
constexpr double lambda = 7.4
 
constexpr double W = 0.25
 
constexpr double h = 0.025
 
constexpr double eta = h
 
constexpr double eta2 = eta * eta
 
constexpr double md = 1e-2
 
constexpr double eps = 1e-12
 
constexpr double tol = std::numeric_limits<float>::epsilon()
 
constexpr double rho_ave = (rho_p + rho_m) / 2
 
constexpr double rho_diff = (rho_p - rho_m) / 2
 
constexpr double mu_ave = (mu_p + mu_m) / 2
 
constexpr double mu_diff = (mu_p - mu_m) / 2
 
const double kappa = (3. / (4. * std::sqrt(2. * W))) * (lambda / eta)
 
auto integration_rule = [](int, int, int) { return 2 * order; }
 
auto cylindrical
 
auto my_max = [](const double x) { return (x - 1 + std::abs(x + 1)) / 2; }
 
auto my_min = [](const double x) { return (x + 1 - std::abs(x - 1)) / 2; }
 
auto cut_off = [](const double h) { return my_max(my_min(h)); }
 
auto d_cut_off
 
auto phase_function
 
auto d_phase_function_h
 
auto get_f = [](const double h) { return 4 * W * h * (h * h - 1); }
 
auto get_f_dh = [](const double h) { return 4 * W * (3 * h * h - 1); }
 
auto get_M0 = [](auto h) { return md; }
 
auto get_M0_dh = [](auto h) { return 0; }
 
auto get_M2
 
auto get_M2_dh
 
auto get_M3
 
auto get_M3_dh
 
auto get_M = [](auto h) { return get_M0(h); }
 
auto get_M_dh = [](auto h) { return get_M0_dh(h); }
 
auto get_D
 
auto kernel_oscillation
 
auto kernel_eye
 
auto init_h
 
auto bit = [](auto b) { return BitRefLevel().set(b); }
 
auto marker = [](auto b) { return BitRefLevel().set(BITREFLEVEL_SIZE - b); }
 
auto save_range
 
auto get_dofs_ents
 

Typedef Documentation

◆ AssemblyBoundaryEleOp

Definition at line 48 of file free_surface.cpp.

◆ AssemblyDomainEleOp

Definition at line 46 of file free_surface.cpp.

◆ BoundaryEle

Definition at line 37 of file free_surface.cpp.

◆ BoundaryEleOp

Definition at line 38 of file free_surface.cpp.

◆ BoundaryParentEle

Definition at line 39 of file free_surface.cpp.

◆ DomainEle

Definition at line 34 of file free_surface.cpp.

◆ DomainEleOp

Definition at line 36 of file free_surface.cpp.

◆ DomianParentEle

Definition at line 35 of file free_surface.cpp.

◆ EntData

Examples
free_surface.cpp.

Definition at line 41 of file free_surface.cpp.

◆ OpBaseTimesScalar

Examples
free_surface.cpp.

Definition at line 62 of file free_surface.cpp.

◆ OpDomainMassH

Examples
free_surface.cpp.

Definition at line 53 of file free_surface.cpp.

◆ OpDomainMassP

Examples
free_surface.cpp.

Definition at line 55 of file free_surface.cpp.

◆ OpDomainMassU

Definition at line 51 of file free_surface.cpp.

◆ OpDomainSourceH

using OpDomainSourceH = FormsIntegrators<DomainEleOp>::Assembly< PETSC>::LinearForm<GAUSS>::OpSource<BASE_DIM, 1>

Definition at line 59 of file free_surface.cpp.

◆ OpDomainSourceU

using OpDomainSourceU = FormsIntegrators<DomainEleOp>::Assembly< PETSC>::LinearForm<GAUSS>::OpSource<BASE_DIM, U_FIELD_DIM>

Definition at line 57 of file free_surface.cpp.

◆ OpMixScalarTimesDiv

Examples
free_surface.cpp.

Definition at line 65 of file free_surface.cpp.

◆ PostProcEdgeEle

Definition at line 44 of file free_surface.cpp.

◆ PostProcEle

Definition at line 43 of file free_surface.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)
Examples
free_surface.cpp.

Definition at line 859 of file free_surface.cpp.

859 {
860
861 // Initialisation of MoFEM/PETSc and MOAB data structures
862 const char param_file[] = "param_file.petsc";
864
865 // Add logging channel for example
866 auto core_log = logging::core::get();
867 core_log->add_sink(LogManager::createSink(LogManager::getStrmWorld(), "FS"));
868 LogManager::setLog("FS");
869 MOFEM_LOG_TAG("FS", "free surface");
870
871 try {
872
873 //! [Register MoFEM discrete manager in PETSc]
874 DMType dm_name = "DMMOFEM";
875 CHKERR DMRegister_MoFEM(dm_name);
876 //! [Register MoFEM discrete manager in PETSc
877
878 //! [Create MoAB]
879 moab::Core mb_instance; ///< mesh database
880 moab::Interface &moab = mb_instance; ///< mesh database interface
881 //! [Create MoAB]
882
883 //! [Create MoFEM]
884 MoFEM::Core core(moab); ///< finite element database
885 MoFEM::Interface &m_field = core; ///< finite element database insterface
886 //! [Create MoFEM]
887
888 //! [FreeSurface]
889 FreeSurface ex(m_field);
890 CHKERR ex.runProblem();
891 //! [FreeSurface]
892 }
894
896}
std::string param_file
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
#define CHKERR
Inline error check.
Definition: definitions.h:535
static char help[]
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:47
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
CoreTmp< 0 > Core
Definition: Core.hpp:1086
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
Core (interface) class.
Definition: Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
Deprecated interface functions.

◆ powof2()

template<int T>
constexpr int powof2 ( )
constexpr
Examples
free_surface.cpp.

Definition at line 87 of file free_surface.cpp.

87 {
88 if constexpr (T == 0)
89 return 1;
90 else
91 return powof2<T - 1>() * 2;
92};
constexpr int powof2()
const double T

Variable Documentation

◆ a0

constexpr double a0 = 0.98
constexpr
Examples
free_surface.cpp.

Definition at line 79 of file free_surface.cpp.

◆ BASE_DIM

constexpr int BASE_DIM = 1
constexpr
Examples
free_surface.cpp.

Definition at line 18 of file free_surface.cpp.

◆ bit

auto bit = [](auto b) { return BitRefLevel().set(b); }
Examples
free_surface.cpp.

Definition at line 204 of file free_surface.cpp.

◆ coord_type

constexpr CoordinateTypes coord_type
constexpr
Initial value:
=
@ CARTESIAN
Definition: definitions.h:115

select coordinate system <CARTESIAN, CYLINDRICAL>;

Examples
free_surface.cpp.

Definition at line 21 of file free_surface.cpp.

◆ cut_off

auto cut_off = [](const double h) { return my_max(my_min(h)); }
Examples
free_surface.cpp.

Definition at line 123 of file free_surface.cpp.

◆ cylindrical

auto cylindrical
Initial value:
= [](const double r) {
if constexpr (coord_type == CYLINDRICAL)
return 2 * M_PI * r;
else
return 1.;
}
@ CYLINDRICAL
Definition: definitions.h:117
constexpr CoordinateTypes coord_type
select coordinate system <CARTESIAN, CYLINDRICAL>;
const double r
rate factor
Examples
free_surface.cpp.

Definition at line 113 of file free_surface.cpp.

◆ d_cut_off

auto d_cut_off
Initial value:
= [](const double h) {
if (h >= -1 && h < 1)
return 1.;
else
return 0.;
}
constexpr double h
Examples
free_surface.cpp.

Definition at line 124 of file free_surface.cpp.

◆ d_phase_function_h

auto d_phase_function_h
Initial value:
= [](const double h, const double diff) {
return diff * d_cut_off(h);
}
auto d_cut_off
Examples
free_surface.cpp.

Definition at line 135 of file free_surface.cpp.

◆ eps

constexpr double eps = 1e-12
constexpr
Examples
free_surface.cpp.

Definition at line 101 of file free_surface.cpp.

◆ eta

constexpr double eta = h
constexpr
Examples
free_surface.cpp, and plot_base.cpp.

Definition at line 96 of file free_surface.cpp.

◆ eta2

constexpr double eta2 = eta * eta
constexpr
Examples
free_surface.cpp.

Definition at line 97 of file free_surface.cpp.

◆ get_D

auto get_D
Initial value:
= [](const double A) {
t_D(i, j, k, l) = A * ((t_kd(i, k) ^ t_kd(j, l)) / 4.);
return t_D;
}
FTensor::Index< 'j', SPACE_DIM > j
FTensor::Index< 'k', SPACE_DIM > k
FTensor::Index< 'i', SPACE_DIM > i
constexpr auto t_kd
FTensor::Index< 'l', SPACE_DIM > l
double A
Examples
free_surface.cpp.

Definition at line 176 of file free_surface.cpp.

◆ get_dofs_ents

auto get_dofs_ents
Initial value:
= [](auto dm) {
auto prb_ptr = getProblemPtr(dm);
std::vector<EntityHandle> ents_vec;
ents_vec.reserve(prb_ptr->numeredRowDofsPtr->size());
for (auto dof : *prb_ptr->numeredRowDofsPtr) {
ents_vec.push_back(dof->getEnt());
}
std::sort(ents_vec.begin(), ents_vec.end());
auto it = std::unique(ents_vec.begin(), ents_vec.end());
r.insert_list(ents_vec.begin(), it);
return r;
}
auto getProblemPtr(DM dm)
get problem pointer from DM
Definition: DMMoFEM.hpp:941
Examples
free_surface.cpp.

Definition at line 218 of file free_surface.cpp.

◆ get_f

auto get_f = [](const double h) { return 4 * W * h * (h * h - 1); }
Examples
free_surface.cpp.

Definition at line 139 of file free_surface.cpp.

◆ get_f_dh

auto get_f_dh = [](const double h) { return 4 * W * (3 * h * h - 1); }
Examples
free_surface.cpp.

Definition at line 140 of file free_surface.cpp.

◆ get_M

auto get_M = [](auto h) { return get_M0(h); }
Examples
free_surface.cpp.

Definition at line 173 of file free_surface.cpp.

◆ get_M0

auto get_M0 = [](auto h) { return md; }
Examples
free_surface.cpp.

Definition at line 142 of file free_surface.cpp.

◆ get_M0_dh

auto get_M0_dh = [](auto h) { return 0; }
Examples
free_surface.cpp.

Definition at line 143 of file free_surface.cpp.

◆ get_M2

auto get_M2
Initial value:
= [](auto h_tmp) {
const double h = cut_off(h_tmp);
return md * (1 - h * h);
}
auto cut_off
constexpr double md
Examples
free_surface.cpp.

Definition at line 145 of file free_surface.cpp.

◆ get_M2_dh

auto get_M2_dh
Initial value:
= [](auto h_tmp) {
const double h = cut_off(h_tmp);
return -md * 2 * h * d_cut_off(h_tmp);
}
Examples
free_surface.cpp.

Definition at line 150 of file free_surface.cpp.

◆ get_M3

auto get_M3
Initial value:
= [](auto h_tmp) {
const double h = cut_off(h_tmp);
const double h2 = h * h;
const double h3 = h2 * h;
if (h >= 0)
return md * (2 * h3 - 3 * h2 + 1);
else
return md * (-2 * h3 - 3 * h2 + 1);
}
Examples
free_surface.cpp.

Definition at line 155 of file free_surface.cpp.

◆ get_M3_dh

auto get_M3_dh
Initial value:
= [](auto h_tmp) {
const double h = cut_off(h_tmp);
if (h >= 0)
return md * (6 * h * (h - 1)) * d_cut_off(h_tmp);
else
return md * (-6 * h * (h + 1)) * d_cut_off(h_tmp);
}
Examples
free_surface.cpp.

Definition at line 165 of file free_surface.cpp.

◆ get_M_dh

auto get_M_dh = [](auto h) { return get_M0_dh(h); }
Examples
free_surface.cpp.

Definition at line 174 of file free_surface.cpp.

◆ h

constexpr double h = 0.025
constexpr
Examples
free_surface.cpp.

Definition at line 95 of file free_surface.cpp.

◆ help

char help[] = "...\n\n"
static
Examples
free_surface.cpp.

Definition at line 14 of file free_surface.cpp.

◆ i

Examples
free_surface.cpp.

Definition at line 68 of file free_surface.cpp.

◆ init_h

auto init_h
Initial value:
= [](double r, double y, double theta) {
return kernel_eye(r, y, theta);
}
auto kernel_eye
Examples
free_surface.cpp, and shallow_wave.cpp.

Definition at line 200 of file free_surface.cpp.

◆ integration_rule

auto integration_rule = [](int, int, int) { return 2 * order; }

◆ j

Examples
free_surface.cpp.

Definition at line 69 of file free_surface.cpp.

◆ k

Examples
free_surface.cpp.

Definition at line 70 of file free_surface.cpp.

◆ kappa

const double kappa = (3. / (4. * std::sqrt(2. * W))) * (lambda / eta)
Examples
free_surface.cpp.

Definition at line 109 of file free_surface.cpp.

◆ kernel_eye

auto kernel_eye
Initial value:
= [](double r, double y, double) {
constexpr double y0 = 0.5;
constexpr double R = 0.4;
y -= y0;
const double d = std::sqrt(r * r + y * y);
return tanh((R - d) / (eta * std::sqrt(2)));
}
constexpr double eta
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
Examples
free_surface.cpp.

Definition at line 192 of file free_surface.cpp.

◆ kernel_oscillation

auto kernel_oscillation
Initial value:
= [](double r, double y, double) {
constexpr int n = 3;
constexpr double R = 0.0125;
constexpr double A = R * 0.2;
const double theta = atan2(r, y);
const double w = R + A * cos(n * theta);
const double d = std::sqrt(r * r + y * y);
return tanh((w - d) / (eta * std::sqrt(2)));
}
FTensor::Index< 'n', SPACE_DIM > n
double w(const double g, const double t)
Examples
free_surface.cpp.

Definition at line 182 of file free_surface.cpp.

◆ l

Examples
free_surface.cpp.

Definition at line 71 of file free_surface.cpp.

◆ lambda

constexpr double lambda = 7.4
constexpr

◆ marker

auto marker = [](auto b) { return BitRefLevel().set(BITREFLEVEL_SIZE - b); }
Examples
free_surface.cpp.

Definition at line 205 of file free_surface.cpp.

◆ md

constexpr double md = 1e-2
constexpr
Examples
free_surface.cpp.

Definition at line 100 of file free_surface.cpp.

◆ mu_ave

constexpr double mu_ave = (mu_p + mu_m) / 2
constexpr
Examples
free_surface.cpp.

Definition at line 106 of file free_surface.cpp.

◆ mu_diff

constexpr double mu_diff = (mu_p - mu_m) / 2
constexpr
Examples
free_surface.cpp.

Definition at line 107 of file free_surface.cpp.

◆ mu_m

constexpr double mu_m = 0.0101
constexpr
Examples
free_surface.cpp.

Definition at line 81 of file free_surface.cpp.

◆ mu_p

constexpr double mu_p = 0.000182
constexpr
Examples
free_surface.cpp.

Definition at line 83 of file free_surface.cpp.

◆ my_max

auto my_max = [](const double x) { return (x - 1 + std::abs(x + 1)) / 2; }
Examples
free_surface.cpp.

Definition at line 121 of file free_surface.cpp.

◆ my_min

auto my_min = [](const double x) { return (x + 1 - std::abs(x - 1)) / 2; }
Examples
free_surface.cpp.

Definition at line 122 of file free_surface.cpp.

◆ order

constexpr int order = 3
constexpr

approximation order

Examples
free_surface.cpp.

Definition at line 76 of file free_surface.cpp.

◆ phase_function

auto phase_function
Initial value:
= [](const double h, const double diff, const double ave) {
return diff * cut_off(h) + ave;
}
Examples
free_surface.cpp.

Definition at line 131 of file free_surface.cpp.

◆ rho_ave

constexpr double rho_ave = (rho_p + rho_m) / 2
constexpr
Examples
free_surface.cpp.

Definition at line 104 of file free_surface.cpp.

◆ rho_diff

constexpr double rho_diff = (rho_p - rho_m) / 2
constexpr
Examples
free_surface.cpp.

Definition at line 105 of file free_surface.cpp.

◆ rho_m

constexpr double rho_m = 0.998
constexpr
Examples
free_surface.cpp.

Definition at line 80 of file free_surface.cpp.

◆ rho_p

constexpr double rho_p = 0.0012
constexpr
Examples
free_surface.cpp.

Definition at line 82 of file free_surface.cpp.

◆ save_range

auto save_range
Initial value:
= [](moab::Interface &moab, const std::string name,
const Range r) {
EntityHandle out_meshset;
CHKERR moab.create_meshset(MESHSET_SET, out_meshset);
CHKERR moab.add_entities(out_meshset, r);
CHKERR moab.write_file(name.c_str(), "VTK", "", &out_meshset, 1);
CHKERR moab.delete_entities(&out_meshset, 1);
}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
Examples
free_surface.cpp.

Definition at line 207 of file free_surface.cpp.

◆ SPACE_DIM

constexpr int SPACE_DIM = 2
constexpr
Examples
free_surface.cpp.

Definition at line 19 of file free_surface.cpp.

◆ t_kd

constexpr auto t_kd = FTensor::Kronecker_Delta_symmetric<int>()
constexpr

◆ tol

constexpr double tol = std::numeric_limits<float>::epsilon()
constexpr
Examples
free_surface.cpp.

Definition at line 102 of file free_surface.cpp.

◆ U_FIELD_DIM

constexpr int U_FIELD_DIM = SPACE_DIM
constexpr
Examples
free_surface.cpp.

Definition at line 20 of file free_surface.cpp.

◆ W

constexpr double W = 0.25
constexpr
Examples
free_surface.cpp.

Definition at line 85 of file free_surface.cpp.