![]() |
v0.15.0 |
Topology optimization using adjoint method with Python objective functions. More...
#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>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 | ObjectiveFunctionData |
| Abstract interface for Python-defined objective functions. More... | |
| struct | Example |
| [Example] More... | |
| struct | OpAdJointTestOp< SPACE_DIM, IntegrationType::GAUSS, DomainBaseOp > |
| Adjoint test operator for elastic problems. More... | |
| struct | OpAdJointGradTimesSymTensor< SPACE_DIM, IntegrationType::GAUSS, DomainBaseOp > |
| struct | OpAdJointObjective |
| Operator for computing objective function and gradients in topology optimization. More... | |
| struct | ObjectiveFunctionDataImpl |
| Implementation of ObjectiveFunctionData interface using Python integration. More... | |
Functions | |
| boost::shared_ptr< ObjectiveFunctionData > | create_python_objective_function (std::string py_file) |
| Factory function to create Python-integrated objective function interface. | |
| 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[]) |
| Main function for topology optimization tutorial using adjoint method. | |
Variables | |
| constexpr int | BASE_DIM = 1 |
| [Constants and material properties] | |
| constexpr int | SPACE_DIM |
| [Define dimension] | |
| constexpr AssemblyType | A = AssemblyType::PETSC |
| [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 = young_modulus / (3 * (1 - 2 * poisson_ratio)) |
| Bulk modulus K = E/(3(1-2ν)) | |
| constexpr double | shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio)) |
| Shear modulus G = E/(2(1+ν)) | |
| PetscBool | is_plane_strain = PETSC_FALSE |
| static char | help [] = "...\n\n" |
| [calculateGradient] | |
Topology optimization using adjoint method with Python objective functions.
This tutorial demonstrates:
Definition in file adjoint.cpp.
Boundary finite elements.
Definition at line 48 of file adjoint.cpp.
Boundary element operators
Definition at line 50 of file adjoint.cpp.
| using BoundaryLhsBCs = NaturalBC<BoundaryEleOp>::Assembly<A>::BiLinearForm<I> |
Boundary LHS natural BCs.
Definition at line 62 of file adjoint.cpp.
| using BoundaryRhsBCs = NaturalBC<BoundaryEleOp>::Assembly<A>::LinearForm<I> |
Boundary RHS natural BCs
Definition at line 60 of file adjoint.cpp.
| using DomainBaseOp = FormsIntegrators<DomainEleOp>::Assembly<A>::OpBase |
[Postprocess results]
[calculateGradient]
Definition at line 1544 of file adjoint.cpp.
Domain finite elements.
Definition at line 47 of file adjoint.cpp.
Domain element operators.
Definition at line 49 of file adjoint.cpp.
| using DomainRhsBCs = NaturalBC<DomainEleOp>::Assembly<A>::LinearForm<I> |
Domain RHS natural BCs.
Definition at line 58 of file adjoint.cpp.
| using EntData = EntitiesFieldData::EntData |
Entity data for field operations.
Definition at line 46 of file adjoint.cpp.
| using OpBoundaryLhsBCs = BoundaryLhsBCs::OpFlux<BoundaryBCs, 1, SPACE_DIM> |
Boundary LHS flux operator.
Definition at line 63 of file adjoint.cpp.
| using OpBoundaryRhsBCs = BoundaryRhsBCs::OpFlux<BoundaryBCs, 1, SPACE_DIM> |
Boundary flux operator.
Definition at line 61 of file adjoint.cpp.
| using OpDomainRhsBCs = DomainRhsBCs::OpFlux<DomainBCs, 1, SPACE_DIM> |
Domain flux operator.
Definition at line 59 of file adjoint.cpp.
Definition at line 81 of file adjoint.cpp.
Definition at line 79 of file adjoint.cpp.
| using SideEle = PostProcEleByDim<SPACE_DIM>::SideEle |
Definition at line 80 of file adjoint.cpp.
| boost::shared_ptr< ObjectiveFunctionData > create_python_objective_function | ( | std::string | py_file | ) |
Factory function to create Python-integrated objective function interface.
Creates and initializes an ObjectiveFunctionDataImpl instance that bridges MoFEM finite element computations with Python-defined objective functions. This enables flexible objective function definition for topology optimization while maintaining computational efficiency.
The Python file must define specific functions with correct signatures:
| py_file | Path to Python file containing objective function definitions |
| MoFEM | exception if Python initialization fails |
Definition at line 2667 of file adjoint.cpp.
|
inline |
Definition at line 1656 of file adjoint.cpp.
| Range get_range_from_block | ( | MoFEM::Interface & | m_field, |
| const std::string | block_name, | ||
| int | dim | ||
| ) |
Definition at line 3290 of file adjoint.cpp.
| int main | ( | int | argc, |
| char * | argv[] | ||
| ) |
Main function for topology optimization tutorial using adjoint method.
This tutorial demonstrates structural topology optimization using:
Workflow:
The adjoint method enables efficient gradient computation with cost independent of the number of design variables, making it suitable for large-scale topology optimization problems.
Required input files:
| argc | Command line argument count |
| argv | Command line argument values |
[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 2304 of file adjoint.cpp.
| MoFEMErrorCode save_range | ( | moab::Interface & | moab, |
| const std::string | name, | ||
| const Range | r | ||
| ) |
Definition at line 3329 of file adjoint.cpp.
|
constexpr |
|
constexpr |
[Constants and material properties]
Dimension of the base functions
Definition at line 25 of file adjoint.cpp.
|
constexpr |
Bulk modulus K = E/(3(1-2ν))
Definition at line 39 of file adjoint.cpp.
|
static |
[calculateGradient]
Definition at line 2272 of file adjoint.cpp.
|
constexpr |
Use Gauss quadrature for integration.
Definition at line 33 of file adjoint.cpp.
| PetscBool is_plane_strain = PETSC_FALSE |
Flag for plane strain vs plane stress in 2D
Definition at line 42 of file adjoint.cpp.
|
constexpr |
Poisson's ratio ν
Definition at line 38 of file adjoint.cpp.
|
constexpr |
Shear modulus G = E/(2(1+ν))
Definition at line 40 of file adjoint.cpp.
|
constexpr |
[Define dimension]
Space dimension of problem (2D or 3D), set at compile time
Definition at line 28 of file adjoint.cpp.
|
constexpr |
[Material properties for linear elasticity]
[Adjoint operator declarations]
Young's modulus E
Definition at line 37 of file adjoint.cpp.