v0.16.0
Loading...
Searching...
No Matches
Classes | Typedefs | Enumerations | Functions | Variables
gradient.cpp File Reference
#include <boost/python.hpp>
#include <boost/python/def.hpp>
#include <boost/python/numpy.hpp>
#include <MoFEM.hpp>
#include <ElasticSpring.hpp>
#include <FluidLevel.hpp>
#include <CalculateTraction.hpp>
#include <NaturalDomainBC.hpp>
#include <NaturalBoundaryBC.hpp>
#include <HookeOps.hpp>
#include <ElasticPostProc.hpp>
#include <ObjectiveFunctionData.hpp>

Go to the source code of this file.

Classes

struct  DomainBCs
 [Define entities] More...
 
struct  BoundaryBCs
 Boundary conditions marker. More...
 
struct  PostProcEleByDim< 2 >
 
struct  PostProcEleByDim< 3 >
 
struct  Example
 [Example] More...
 
struct  OpStateSensitivity
 
struct  OpObjective
 
struct  OpAdJointObjective
 

Typedefs

using EntData = EntitiesFieldData::EntData
 Entity data for field operations.
 
using DomainEle = PipelineManager::ElementsAndOpsByDim< SPACE_DIM >::DomainEle
 Domain finite elements.
 
using BoundaryEle = PipelineManager::ElementsAndOpsByDim< SPACE_DIM >::BoundaryEle
 Boundary finite elements.
 
using DomainEleOp = DomainEle::UserDataOperator
 Domain element operators.
 
using BoundaryEleOp = BoundaryEle::UserDataOperator
 
using DomainRhsBCs = NaturalBC< DomainEleOp >::Assembly< A >::LinearForm< I >
 Domain RHS natural BCs.
 
using OpDomainRhsBCs = DomainRhsBCs::OpFlux< DomainBCs, 1, SPACE_DIM >
 Domain flux operator.
 
using BoundaryRhsBCs = NaturalBC< BoundaryEleOp >::Assembly< A >::LinearForm< I >
 Boundary RHS natural BCs.
 
using OpBoundaryRhsBCs = BoundaryRhsBCs::OpFlux< BoundaryBCs, 1, SPACE_DIM >
 Boundary flux operator.
 
using BoundaryLhsBCs = NaturalBC< BoundaryEleOp >::Assembly< A >::BiLinearForm< I >
 Boundary LHS natural BCs.
 
using OpBoundaryLhsBCs = BoundaryLhsBCs::OpFlux< BoundaryBCs, 1, SPACE_DIM >
 Boundary LHS flux operator.
 
using PostProcEleDomain = PostProcEleByDim< SPACE_DIM >::PostProcEleDomain
 
using SideEle = PostProcEleByDim< SPACE_DIM >::SideEle
 
using PostProcEleBdy = PostProcEleByDim< SPACE_DIM >::PostProcEleBdy
 
using DomainBaseOp = FormsIntegrators< DomainEleOp >::Assembly< A >::OpBase
 [Postprocess results]
 

Enumerations

enum  SensitivityMethod { DIRECT , ADJOINT }
 

Functions

Range get_range_from_block (MoFEM::Interface &m_field, const std::string block_name, int dim)
 
MoFEMErrorCode save_range (moab::Interface &moab, const std::string name, const Range r)
 
template<int DIM>
auto diff_symmetrize (FTensor::Number< DIM >)
 
int main (int argc, char *argv[])
 

Variables

constexpr int BASE_DIM = 1
 [Constants and material properties]
 
constexpr int SPACE_DIM
 [Define dimension]
 
constexpr AssemblyType A
 [Define dimension]
 
constexpr IntegrationType I
 Use Gauss quadrature for integration.
 
constexpr double young_modulus = 1
 [Material properties for linear elasticity]
 
constexpr double poisson_ratio = 0.3
 Poisson's ratio ν
 
constexpr double bulk_modulus_K
 Bulk modulus K = E/(3(1-2ν))
 
constexpr double shear_modulus_G
 Shear modulus G = E/(2(1+ν))
 
PetscBool is_plane_strain
 
SensitivityMethod derivative_type = ADJOINT
 
static char help [] = "...\n\n"
 [Finite difference check]
 

Typedef Documentation

◆ BoundaryEle

Boundary finite elements.

Definition at line 46 of file gradient.cpp.

◆ BoundaryEleOp

Boundary element operators

Definition at line 49 of file gradient.cpp.

◆ BoundaryLhsBCs

Boundary LHS natural BCs.

Definition at line 67 of file gradient.cpp.

◆ BoundaryRhsBCs

Boundary RHS natural BCs.

Definition at line 62 of file gradient.cpp.

◆ DomainBaseOp

[Postprocess results]

Examples
mofem/tutorials/vec-7_shape_optimisation/adjoint.cpp.

Definition at line 703 of file gradient.cpp.

◆ DomainEle

Domain finite elements.

Definition at line 44 of file gradient.cpp.

◆ DomainEleOp

Domain element operators.

Definition at line 48 of file gradient.cpp.

◆ DomainRhsBCs

Domain RHS natural BCs.

Definition at line 58 of file gradient.cpp.

◆ EntData

Entity data for field operations.

Definition at line 42 of file gradient.cpp.

◆ OpBoundaryLhsBCs

Boundary LHS flux operator.

Definition at line 69 of file gradient.cpp.

◆ OpBoundaryRhsBCs

Boundary flux operator.

Definition at line 64 of file gradient.cpp.

◆ OpDomainRhsBCs

Domain flux operator.

Definition at line 60 of file gradient.cpp.

◆ PostProcEleBdy

Definition at line 89 of file gradient.cpp.

◆ PostProcEleDomain

Definition at line 87 of file gradient.cpp.

◆ SideEle

Definition at line 88 of file gradient.cpp.

Enumeration Type Documentation

◆ SensitivityMethod

Enumerator
DIRECT 
ADJOINT 

Definition at line 94 of file gradient.cpp.

94{ DIRECT, ADJOINT };
@ DIRECT
Definition gradient.cpp:94
@ ADJOINT
Definition gradient.cpp:94

Function Documentation

◆ diff_symmetrize()

template<int DIM>
auto diff_symmetrize ( FTensor::Number< DIM >  )
inline
Examples
mofem/tutorials/vec-7_shape_optimisation/adjoint.cpp.

Definition at line 706 of file gradient.cpp.

706 {
707
708 FTensor::Index<'i', DIM> i;
709 FTensor::Index<'j', DIM> j;
710 FTensor::Index<'k', DIM> k;
711 FTensor::Index<'l', DIM> l;
712
714
715 t_diff(i, j, k, l) = 0;
716 t_diff(0, 0, 0, 0) = 1;
717 t_diff(1, 1, 1, 1) = 1;
718
719 t_diff(1, 0, 1, 0) = 0.5;
720 t_diff(1, 0, 0, 1) = 0.5;
721
722 t_diff(0, 1, 0, 1) = 0.5;
723 t_diff(0, 1, 1, 0) = 0.5;
724
725 if constexpr (DIM == 3) {
726 t_diff(2, 2, 2, 2) = 1;
727
728 t_diff(2, 0, 2, 0) = 0.5;
729 t_diff(2, 0, 0, 2) = 0.5;
730 t_diff(0, 2, 0, 2) = 0.5;
731 t_diff(0, 2, 2, 0) = 0.5;
732
733 t_diff(2, 1, 2, 1) = 0.5;
734 t_diff(2, 1, 1, 2) = 0.5;
735 t_diff(1, 2, 1, 2) = 0.5;
736 t_diff(1, 2, 2, 1) = 0.5;
737 }
738
739 return t_diff;
740};
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k

◆ get_range_from_block()

Range get_range_from_block ( MoFEM::Interface m_field,
const std::string  block_name,
int  dim 
)

Definition at line 1385 of file gradient.cpp.

1386 {
1387 Range r;
1388
1389 auto mesh_mng = m_field.getInterface<MeshsetsManager>();
1390 auto bcs = mesh_mng->getCubitMeshsetPtr(
1391
1392 std::regex((boost::format("%s(.*)") % block_name).str())
1393
1394 );
1395
1396 for (auto bc : bcs) {
1397 Range faces;
1398 CHK_MOAB_THROW(bc->getMeshsetIdEntitiesByDimension(m_field.get_moab(), dim,
1399 faces, true),
1400 "get meshset ents");
1401 r.merge(faces);
1402 }
1403
1404 for (auto dd = dim - 1; dd >= 0; --dd) {
1405 if (dd >= 0) {
1406 Range ents;
1407 CHK_MOAB_THROW(m_field.get_moab().get_adjacencies(r, dd, false, ents,
1408 moab::Interface::UNION),
1409 "get adjs");
1410 r.merge(ents);
1411 } else {
1412 Range verts;
1413 CHK_MOAB_THROW(m_field.get_moab().get_connectivity(r, verts),
1414 "get verts");
1415 r.merge(verts);
1416 }
1418 m_field.getInterface<CommInterface>()->synchroniseEntities(r), "comm");
1419 }
1420
1421 return r;
1422};
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition ddTensor0.hpp:33
int r
Definition sdf.py:205
Managing BitRefLevels.
MoFEMErrorCode synchroniseEntities(Range &ent, std::map< int, Range > *received_ents, int verb=DEFAULT_VERBOSITY)
synchronize entity range on processors (collective)
virtual moab::Interface & get_moab()=0
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ main()

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

[Register MoFEM discrete manager in PETSc]

[Register MoFEM discrete manager in PETSc

[Create MoAB]

< mesh database

< mesh database interface

[Create MoAB]

[Create MoFEM]

< finite element database

< finite element database interface

[Create MoFEM]

[Example]

[Example]

Definition at line 1332 of file gradient.cpp.

1332 {
1333
1334 // Initialize Python environment for objective function interface
1335 Py_Initialize();
1336 np::initialize();
1337
1338 // Initialize MoFEM/PETSc and MOAB data structures
1339 const char param_file[] = "param_file.petsc";
1340 MoFEM::Core::Initialize(&argc, &argv, param_file, help);
1341
1342 auto core_log = logging::core::get();
1343 core_log->add_sink(
1345
1346 core_log->add_sink(
1347 LogManager::createSink(LogManager::getStrmSync(), "FieldEvaluator"));
1348 LogManager::setLog("FieldEvaluator");
1349 MOFEM_LOG_TAG("FieldEvaluator", "field_eval");
1350
1351 try {
1352
1353 //! [Register MoFEM discrete manager in PETSc]
1354 DMType dm_name = "DMMOFEM";
1355 CHKERR DMRegister_MoFEM(dm_name);
1356 DMType dm_name_mg = "DMMOFEM_MG";
1358 //! [Register MoFEM discrete manager in PETSc
1359
1360 //! [Create MoAB]
1361 moab::Core mb_instance; ///< mesh database
1362 moab::Interface &moab = mb_instance; ///< mesh database interface
1363 //! [Create MoAB]
1364
1365 //! [Create MoFEM]
1366 MoFEM::Core core(moab); ///< finite element database
1367 MoFEM::Interface &m_field = core; ///< finite element database interface
1368 //! [Create MoFEM]
1369
1370 //! [Example]
1371
1372 Example ex(m_field);
1373 CHKERR ex.runProblem();
1374 //! [Example]
1375 }
1377
1379
1380 if (Py_FinalizeEx() < 0) {
1381 exit(120);
1382 }
1383}
#define CATCH_ERRORS
Catch errors.
#define CHKERR
Inline error check.
static char help[]
[Finite difference check]
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition DMMoFEM.cpp:43
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
[Example]
Definition plastic.cpp:217
Core (interface) class.
Definition Core.hpp:83
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition Core.cpp:68
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition Core.cpp:123
Deprecated interface functions.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.

◆ save_range()

MoFEMErrorCode save_range ( moab::Interface &  moab,
const std::string  name,
const Range  r 
)

Definition at line 1424 of file gradient.cpp.

1425 {
1427 auto out_meshset = get_temp_meshset_ptr(moab);
1428 CHKERR moab.add_entities(*out_meshset, r);
1429 CHKERR moab.write_file(name.c_str(), "VTK", "", out_meshset->get_ptr(), 1);
1431};
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temporary meshset.

Variable Documentation

◆ A

constexpr AssemblyType A
constexpr
Initial value:
=
AssemblyType::PETSC

[Define dimension]

Use PETSc for matrix/vector assembly

Definition at line 23 of file gradient.cpp.

◆ BASE_DIM

constexpr int BASE_DIM = 1
constexpr

[Constants and material properties]

Dimension of the base functions

Definition at line 16 of file gradient.cpp.

◆ bulk_modulus_K

constexpr double bulk_modulus_K
constexpr
Initial value:
=
(3 * (1 - 2 * poisson_ratio))
constexpr double poisson_ratio
Poisson's ratio ν
Definition gradient.cpp:30
constexpr double young_modulus
[Material properties for linear elasticity]
Definition gradient.cpp:29

Bulk modulus K = E/(3(1-2ν))

Definition at line 31 of file gradient.cpp.

◆ derivative_type

SensitivityMethod derivative_type = ADJOINT

Definition at line 96 of file gradient.cpp.

◆ help

char help[] = "...\n\n"
static

[Finite difference check]

Definition at line 1330 of file gradient.cpp.

◆ I

constexpr IntegrationType I
constexpr
Initial value:
=
IntegrationType::GAUSS

Use Gauss quadrature for integration.

Examples
mofem/tutorials/vec-7_shape_optimisation/adjoint.cpp.

Definition at line 25 of file gradient.cpp.

◆ is_plane_strain

PetscBool is_plane_strain
Initial value:
=
PETSC_FALSE

Flag for plane strain vs plane stress in 2D

Definition at line 37 of file gradient.cpp.

◆ poisson_ratio

constexpr double poisson_ratio = 0.3
constexpr

Poisson's ratio ν

Definition at line 30 of file gradient.cpp.

◆ shear_modulus_G

constexpr double shear_modulus_G
constexpr
Initial value:
=

Shear modulus G = E/(2(1+ν))

Definition at line 34 of file gradient.cpp.

◆ SPACE_DIM

constexpr int SPACE_DIM
constexpr
Initial value:
=
#define EXECUTABLE_DIMENSION
Definition plastic.cpp:13

[Define dimension]

Space dimension of problem (2D or 3D), set at compile time

Examples
mofem/tutorials/vec-7_shape_optimisation/adjoint.cpp.

Definition at line 19 of file gradient.cpp.

◆ young_modulus

constexpr double young_modulus = 1
constexpr

[Material properties for linear elasticity]

Young's modulus E

Definition at line 29 of file gradient.cpp.