v0.15.0
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Attributes | List of all members
OpAdJointObjective Struct Reference

Operator for computing objective function and gradients in topology optimization. More...

Inheritance diagram for OpAdJointObjective:
[legend]
Collaboration diagram for OpAdJointObjective:
[legend]

Public Types

using OP = ForcesAndSourcesCore::UserDataOperator
 

Public Member Functions

 OpAdJointObjective (boost::shared_ptr< ObjectiveFunctionData > python_ptr, boost::shared_ptr< HookeOps::CommonData > comm_ptr, boost::shared_ptr< MatrixDouble > jac_ptr, boost::shared_ptr< MatrixDouble > diff_jac, boost::shared_ptr< VectorDouble > cof_vals, boost::shared_ptr< MatrixDouble > d_grad_ptr, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< double > glob_objective_ptr, boost::shared_ptr< double > glob_objective_grad_ptr)
 
MoFEMErrorCode doWork (int side, EntityType type, EntitiesFieldData::EntData &data)
 Compute objective function contributions at element level.
 

Private Attributes

boost::shared_ptr< ObjectiveFunctionDatapythonPtr
 
boost::shared_ptr< HookeOps::CommonDatacommPtr
 
boost::shared_ptr< MatrixDoublejacPtr
 
boost::shared_ptr< MatrixDoublediffJacPtr
 
boost::shared_ptr< VectorDoublecofVals
 
boost::shared_ptr< MatrixDoubledGradPtr
 
boost::shared_ptr< MatrixDoubleuPtr
 
boost::shared_ptr< doubleglobObjectivePtr
 
boost::shared_ptr< doubleglobObjectiveGradPtr
 

Detailed Description

Operator for computing objective function and gradients in topology optimization.

This operator interfaces with Python-defined objective functions to:

  1. Evaluate objective function f(u, x) at current design point
  2. Compute gradients ∂f/∂u and ∂f/∂x for adjoint sensitivity analysis

The objective function typically depends on:

This enables flexible objective function definition (compliance, stress constraints, volume constraints, etc.) through Python scripting.

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

Definition at line 1801 of file adjoint.cpp.

Member Typedef Documentation

◆ OP

using OpAdJointObjective::OP = ForcesAndSourcesCore::UserDataOperator

Definition at line 1802 of file adjoint.cpp.

Constructor & Destructor Documentation

◆ OpAdJointObjective()

OpAdJointObjective::OpAdJointObjective ( boost::shared_ptr< ObjectiveFunctionData python_ptr,
boost::shared_ptr< HookeOps::CommonData comm_ptr,
boost::shared_ptr< MatrixDouble jac_ptr,
boost::shared_ptr< MatrixDouble diff_jac,
boost::shared_ptr< VectorDouble cof_vals,
boost::shared_ptr< MatrixDouble d_grad_ptr,
boost::shared_ptr< MatrixDouble u_ptr,
boost::shared_ptr< double glob_objective_ptr,
boost::shared_ptr< double glob_objective_grad_ptr 
)
inline

Definition at line 1804 of file adjoint.cpp.

1813 : OP(NOSPACE, OP::OPSPACE), pythonPtr(python_ptr), commPtr(comm_ptr),
1814 jacPtr(jac_ptr), diffJacPtr(diff_jac), cofVals(cof_vals),
1815 dGradPtr(d_grad_ptr), uPtr(u_ptr), globObjectivePtr(glob_objective_ptr),
1816 globObjectiveGradPtr(glob_objective_grad_ptr) {}
@ NOSPACE
Definition definitions.h:83
boost::shared_ptr< MatrixDouble > dGradPtr
Definition adjoint.cpp:1947
boost::shared_ptr< double > globObjectiveGradPtr
Definition adjoint.cpp:1951
boost::shared_ptr< ObjectiveFunctionData > pythonPtr
Definition adjoint.cpp:1942
boost::shared_ptr< double > globObjectivePtr
Definition adjoint.cpp:1950
ForcesAndSourcesCore::UserDataOperator OP
Definition adjoint.cpp:1802
boost::shared_ptr< HookeOps::CommonData > commPtr
Definition adjoint.cpp:1943
boost::shared_ptr< VectorDouble > cofVals
Definition adjoint.cpp:1946
boost::shared_ptr< MatrixDouble > uPtr
Definition adjoint.cpp:1948
boost::shared_ptr< MatrixDouble > diffJacPtr
Definition adjoint.cpp:1945
boost::shared_ptr< MatrixDouble > jacPtr
Definition adjoint.cpp:1944

Member Function Documentation

◆ doWork()

MoFEMErrorCode OpAdJointObjective::doWork ( int  side,
EntityType  type,
EntitiesFieldData::EntData data 
)
inline

Compute objective function contributions at element level.

Evaluates Python objective function with current displacement and stress state, and accumulates global objective value and gradients.

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

Definition at line 1824 of file adjoint.cpp.

1825 {
1827
1828 // Define tensor indices for calculations
1833
1836
1837 constexpr auto symm_size = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
1838
1839 auto t_diff_symm = diff_symmetrize(FTensor::Number<SPACE_DIM>());
1840
1841 auto nb_gauss_pts = getGaussPts().size2();
1842 auto objective_ptr = boost::make_shared<VectorDouble>(nb_gauss_pts);
1843 auto objective_dstress =
1844 boost::make_shared<MatrixDouble>(symm_size, nb_gauss_pts);
1845 auto objective_dstrain =
1846 boost::make_shared<MatrixDouble>(symm_size, nb_gauss_pts);
1847 auto obj_grad = boost::make_shared<MatrixDouble>(SPACE_DIM, nb_gauss_pts);
1848
1849 auto evaluate_python = [&]() {
1851 auto &coords = OP::getCoordsAtGaussPts();
1852 CHKERR pythonPtr->evalObjectiveFunction(
1853 coords, uPtr, commPtr->getMatCauchyStress(), commPtr->getMatStrain(),
1854 objective_ptr);
1855 CHKERR pythonPtr->evalObjectiveGradientStress(
1856 coords, uPtr, commPtr->getMatCauchyStress(), commPtr->getMatStrain(),
1857 objective_dstress);
1858 CHKERR pythonPtr->evalObjectiveGradientStrain(
1859 coords, uPtr, commPtr->getMatCauchyStress(), commPtr->getMatStrain(),
1860 objective_dstrain);
1861
1862 auto t_grad_u =
1863 getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(commPtr->matGradPtr));
1864 auto t_D =
1865 getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*(commPtr->matDPtr));
1866 auto t_jac = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(jacPtr));
1867 auto t_diff_jac = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(diffJacPtr));
1868 auto t_cof = getFTensor0FromVec(*(cofVals));
1869 auto t_d_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(dGradPtr));
1870
1871 auto t_obj = getFTensor0FromVec(*objective_ptr);
1872 auto t_obj_dstress =
1873 getFTensor2SymmetricFromMat<SPACE_DIM>(*objective_dstress);
1874 auto t_obj_dstrain =
1875 getFTensor2SymmetricFromMat<SPACE_DIM>(*objective_dstrain);
1876
1877 auto vol = OP::getMeasure();
1878 auto t_w = getFTensor0IntegrationWeight();
1879 for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
1880
1881 auto t_det = determinantTensor(t_jac);
1883 CHKERR invertTensor(t_jac, t_det, t_inv_jac);
1884
1886 t_diff_inv_jac(i, j) =
1887 -(t_inv_jac(i, I) * t_diff_jac(I, J)) * t_inv_jac(J, j);
1889 t_diff_grad(i, j) = t_grad_u(i, k) * t_diff_inv_jac(k, j);
1890
1892 t_d_strain(i, j) = t_diff_symm(i, j, k, l) * (
1893
1894 t_d_grad(k, l)
1895
1896 +
1897
1898 t_diff_grad(k, l)
1899
1900 );
1901
1902 auto alpha = t_w * vol;
1903
1904 (*globObjectivePtr) += alpha * t_obj;
1905 (*globObjectiveGradPtr) +=
1906 alpha *
1907 (
1908
1909 t_obj_dstress(i, j) * (t_D(i, j, k, l) * t_d_strain(k, l))
1910
1911 +
1912
1913 t_obj_dstrain(i, j) * t_d_strain(i, j)
1914
1915 +
1916
1917 t_obj * t_cof
1918
1919 );
1920
1921 ++t_w;
1922 ++t_jac;
1923 ++t_diff_jac;
1924 ++t_cof;
1925
1926 ++t_obj;
1927 ++t_obj_dstress;
1928 ++t_obj_dstrain;
1929
1930 ++t_grad_u;
1931 ++t_d_grad;
1932 }
1934 };
1935
1936 CHKERR evaluate_python();
1937
1939 }
#define FTENSOR_INDEX(DIM, I)
auto diff_symmetrize(FTensor::Number< DIM >)
Definition adjoint.cpp:1656
constexpr int SPACE_DIM
[Define dimension]
Definition adjoint.cpp:28
constexpr IntegrationType I
Use Gauss quadrature for integration.
Definition adjoint.cpp:33
#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()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
static MoFEMErrorCode invertTensor(FTensor::Tensor2< T1, DIM, DIM > &t, T2 &det, FTensor::Tensor2< T3, DIM, DIM > &inv_t)
static auto determinantTensor(FTensor::Tensor2< T, DIM, DIM > &t)
Calculate the determinant of a tensor of rank DIM.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.

Member Data Documentation

◆ cofVals

boost::shared_ptr<VectorDouble> OpAdJointObjective::cofVals
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1946 of file adjoint.cpp.

◆ commPtr

boost::shared_ptr<HookeOps::CommonData> OpAdJointObjective::commPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1943 of file adjoint.cpp.

◆ dGradPtr

boost::shared_ptr<MatrixDouble> OpAdJointObjective::dGradPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1947 of file adjoint.cpp.

◆ diffJacPtr

boost::shared_ptr<MatrixDouble> OpAdJointObjective::diffJacPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1945 of file adjoint.cpp.

◆ globObjectiveGradPtr

boost::shared_ptr<double> OpAdJointObjective::globObjectiveGradPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1951 of file adjoint.cpp.

◆ globObjectivePtr

boost::shared_ptr<double> OpAdJointObjective::globObjectivePtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1950 of file adjoint.cpp.

◆ jacPtr

boost::shared_ptr<MatrixDouble> OpAdJointObjective::jacPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1944 of file adjoint.cpp.

◆ pythonPtr

boost::shared_ptr<ObjectiveFunctionData> OpAdJointObjective::pythonPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1942 of file adjoint.cpp.

◆ uPtr

boost::shared_ptr<MatrixDouble> OpAdJointObjective::uPtr
private
Examples
mofem/tutorials/vec-7/adjoint.cpp.

Definition at line 1948 of file adjoint.cpp.


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