AssembleOptimization Class Reference

This class encapsulate all functionality required for assembling the objective function, gradient, and hessian. More...

Inheritance diagram for AssembleOptimization:

Public Member Functions
	AssembleOptimization (OptimizationSystem &sys_in)
	Constructor.
void	assemble_A_and_F ()
	The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F.
virtual Number	objective (const NumericVector< Number > &soln, OptimizationSystem &)
	Evaluate the objective function.
virtual void	gradient (const NumericVector< Number > &soln, NumericVector< Number > &grad_f, OptimizationSystem &)
	Evaluate the gradient.
virtual void	hessian (const NumericVector< Number > &soln, SparseMatrix< Number > &H_f, OptimizationSystem &)
	Evaluate the Hessian.
	AssembleOptimization (OptimizationSystem &sys_in)
	Constructor.
void	assemble_A_and_F ()
	The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F.
virtual Number	objective (const NumericVector< Number > &soln, OptimizationSystem &)
	Evaluate the objective function.
virtual void	gradient (const NumericVector< Number > &soln, NumericVector< Number > &grad_f, OptimizationSystem &)
	Evaluate the gradient.
virtual void	hessian (const NumericVector< Number > &soln, SparseMatrix< Number > &H_f, OptimizationSystem &)
	Evaluate the Hessian.
virtual void	equality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_eq, OptimizationSystem &)
	Evaluate the equality constraints.
virtual void	equality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_eq_jac, OptimizationSystem &)
	Evaluate the equality constraints Jacobian.
virtual void	inequality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_ineq, OptimizationSystem &)
	Evaluate the inequality constraints.
virtual void	inequality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, OptimizationSystem &)
	Evaluate the inequality constraints Jacobian.
virtual void	lower_and_upper_bounds (OptimizationSystem &sys)
	Evaluate the lower and upper bounds vectors.

Public Attributes
SparseMatrix< Number > *	A_matrix
	Sparse matrix for storing the matrix A.
NumericVector< Number > *	F_vector
	Vector for storing F.

Private Attributes
OptimizationSystem &	_sys
	Keep a reference to the OptimizationSystem.

Detailed Description

This class encapsulate all functionality required for assembling the objective function, gradient, and hessian.

This class encapsulate all functionality required for assembling the objective function, gradient, hessian, and constraints.

Definition at line 67 of file optimization_ex1.C.

Constructor & Destructor Documentation

AssembleOptimization() [1/2]

AssembleOptimization::AssembleOptimization

(

OptimizationSystem &

sys_in

)

Constructor.

Definition at line 127 of file optimization_ex1.C.

                                                                      :
  _sys(sys_in)
{}

AssembleOptimization() [2/2]

AssembleOptimization::AssembleOptimization

(

OptimizationSystem &

sys_in

)

Constructor.

Member Function Documentation

assemble_A_and_F() [1/2]

void AssembleOptimization::assemble_A_and_F

(

)

The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F.

In this method, we assemble A and F.

Definition at line 131 of file optimization_ex1.C.

{
  A_matrix->zero();
  F_vector->zero();
 
  const MeshBase & mesh = _sys.get_mesh();
 
  const unsigned int dim = mesh.mesh_dimension();
  const unsigned int u_var = _sys.variable_number ("u");
 
  const DofMap & dof_map = _sys.get_dof_map();
  FEType fe_type = dof_map.variable_type(u_var);
  std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
  QGauss qrule (dim, fe_type.default_quadrature_order());
  fe->attach_quadrature_rule (&qrule);
 
  const std::vector<Real> & JxW = fe->get_JxW();
  const std::vector<std::vector<Real>> & phi = fe->get_phi();
  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
 
  std::vector<dof_id_type> dof_indices;
 
  DenseMatrix<Number> Ke;
  DenseVector<Number> Fe;
 
  for (const auto & elem : mesh.active_local_element_ptr_range())
    {
      dof_map.dof_indices (elem, dof_indices);
 
      const unsigned int n_dofs = dof_indices.size();
 
      fe->reinit (elem);
 
      Ke.resize (n_dofs, n_dofs);
      Fe.resize (n_dofs);
 
      for (unsigned int qp=0; qp<qrule.n_points(); qp++)
        {
          for (unsigned int dof_i=0; dof_i<n_dofs; dof_i++)
            {
              for (unsigned int dof_j=0; dof_j<n_dofs; dof_j++)
                {
                  Ke(dof_i, dof_j) += JxW[qp] * (dphi[dof_j][qp]* dphi[dof_i][qp]);
                }
              Fe(dof_i) += JxW[qp] * phi[dof_i][qp];
            }
        }
 
#ifdef LIBMESH_ENABLE_CONSTRAINTS
      // This will zero off-diagonal entries of Ke corresponding to
      // Dirichlet dofs.
      dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
 
      // We want the diagonal of constrained dofs to be zero too
      for (unsigned int local_dof_index=0; local_dof_index<n_dofs; local_dof_index++)
        {
          dof_id_type global_dof_index = dof_indices[local_dof_index];
          if (dof_map.is_constrained_dof(global_dof_index))
            {
              Ke(local_dof_index, local_dof_index) = 0.;
            }
        }
#endif // LIBMESH_ENABLE_CONSTRAINTS
 
      A_matrix->add_matrix (Ke, dof_indices);
      F_vector->add_vector (Fe, dof_indices);
    }
 
  A_matrix->close();
  F_vector->close();
}

References _sys, A_matrix, libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::NumericVector< T >::close(), libMesh::SparseMatrix< T >::close(), libMesh::DofMap::constrain_element_matrix_and_vector(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), F_vector, libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::DofMap::is_constrained_dof(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::System::variable_number(), libMesh::DofMap::variable_type(), libMesh::NumericVector< T >::zero(), and libMesh::SparseMatrix< T >::zero().

Referenced by main().

assemble_A_and_F() [2/2]

void AssembleOptimization::assemble_A_and_F

(

)

The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F.

In this method, we assemble A and F.

equality_constraints()

void AssembleOptimization::equality_constraints ( const NumericVector< Number > &  X, NumericVector< Number > &  C_eq, OptimizationSystem &    )

virtual

Evaluate the equality constraints.

Implements libMesh::OptimizationSystem::ComputeEqualityConstraints.

Definition at line 276 of file optimization_ex2.C.

{
  C_eq.zero();
 
  std::unique_ptr<NumericVector<Number>> X_localized =
    NumericVector<Number>::build(X.comm());
  X_localized->init(X.size(), false, SERIAL);
  X.localize(*X_localized);
 
  std::vector<Number> constraint_values(3);
  constraint_values[0] = (*X_localized)(17);
  constraint_values[1] = (*X_localized)(23);
  constraint_values[2] = (*X_localized)(98) + (*X_localized)(185);
 
  for (std::size_t i=0; i<constraint_values.size(); i++)
    if ((C_eq.first_local_index() <= i) &&
        (i < C_eq.last_local_index()))
      C_eq.set(i, constraint_values[i]);
}

References libMesh::NumericVector< T >::build(), libMesh::ParallelObject::comm(), libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::last_local_index(), libMesh::NumericVector< T >::localize(), libMesh::SERIAL, libMesh::NumericVector< T >::set(), libMesh::NumericVector< T >::size(), and libMesh::NumericVector< T >::zero().

equality_constraints_jacobian()

void AssembleOptimization::equality_constraints_jacobian ( const NumericVector< Number > &  X, SparseMatrix< Number > &  C_eq_jac, OptimizationSystem &  sys  )

virtual

Evaluate the equality constraints Jacobian.

Implements libMesh::OptimizationSystem::ComputeEqualityConstraintsJacobian.

Definition at line 298 of file optimization_ex2.C.

{
  C_eq_jac.zero();
 
  std::vector<std::vector<Number>> constraint_jac_values(3);
  std::vector<std::vector<dof_id_type>> constraint_jac_indices(3);
 
  constraint_jac_values[0].resize(1);
  constraint_jac_indices[0].resize(1);
  constraint_jac_values[0][0] = 1.;
  constraint_jac_indices[0][0] = 17;
 
  constraint_jac_values[1].resize(1);
  constraint_jac_indices[1].resize(1);
  constraint_jac_values[1][0] = 1.;
  constraint_jac_indices[1][0] = 23;
 
  constraint_jac_values[2].resize(2);
  constraint_jac_indices[2].resize(2);
  constraint_jac_values[2][0] = 1.;
  constraint_jac_values[2][1] = 1.;
  constraint_jac_indices[2][0] = 98;
  constraint_jac_indices[2][1] = 185;
 
  for (std::size_t i=0; i<constraint_jac_values.size(); i++)
    for (std::size_t j=0; j<constraint_jac_values[i].size(); j++)
      if ((sys.C_eq->first_local_index() <= i) &&
          (i < sys.C_eq->last_local_index()))
        {
          dof_id_type col_index = constraint_jac_indices[i][j];
          Number value = constraint_jac_values[i][j];
          C_eq_jac.set(i, col_index, value);
        }
}

References libMesh::OptimizationSystem::C_eq, libMesh::SparseMatrix< T >::set(), value, and libMesh::SparseMatrix< T >::zero().

gradient() [1/2]

void AssembleOptimization::gradient ( const NumericVector< Number > &  soln, NumericVector< Number > &  grad_f, OptimizationSystem &    )

virtual

Evaluate the gradient.

Implements libMesh::OptimizationSystem::ComputeGradient.

Definition at line 217 of file optimization_ex1.C.

{
  grad_f.zero();
 
  // Since we've enforced constraints on soln, A and F,
  // this automatically sets grad_f to zero for constrained
  // dofs.
  A_matrix->vector_mult(grad_f, soln);
  grad_f.add(-1, *F_vector);
}

References A_matrix, libMesh::NumericVector< T >::add(), F_vector, libMesh::SparseMatrix< T >::vector_mult(), and libMesh::NumericVector< T >::zero().

gradient() [2/2]

virtual void AssembleOptimization::gradient ( const NumericVector< Number > &  soln, NumericVector< Number > &  grad_f, OptimizationSystem &    )

virtual

Evaluate the gradient.

Implements libMesh::OptimizationSystem::ComputeGradient.

hessian() [1/2]

void AssembleOptimization::hessian ( const NumericVector< Number > &  soln, SparseMatrix< Number > &  H_f, OptimizationSystem &  sys  )

virtual

Evaluate the Hessian.

Implements libMesh::OptimizationSystem::ComputeHessian.

Definition at line 231 of file optimization_ex1.C.

{
  LOG_SCOPE("hessian()", "AssembleOptimization");
  H_f.zero();
  H_f.add(1., *A_matrix);
}

References A_matrix, libMesh::SparseMatrix< T >::add(), and libMesh::SparseMatrix< T >::zero().

hessian() [2/2]

virtual void AssembleOptimization::hessian ( const NumericVector< Number > &  soln, SparseMatrix< Number > &  H_f, OptimizationSystem &    )

virtual

Evaluate the Hessian.

Implements libMesh::OptimizationSystem::ComputeHessian.

inequality_constraints()

void AssembleOptimization::inequality_constraints ( const NumericVector< Number > &  X, NumericVector< Number > &  C_ineq, OptimizationSystem &    )

virtual

Evaluate the inequality constraints.

Implements libMesh::OptimizationSystem::ComputeInequalityConstraints.

Definition at line 335 of file optimization_ex2.C.

{
  C_ineq.zero();
 
  std::unique_ptr<NumericVector<Number>> X_localized =
    NumericVector<Number>::build(X.comm());
  X_localized->init(X.size(), false, SERIAL);
  X.localize(*X_localized);
 
  std::vector<Number> constraint_values(1);
  constraint_values[0] = (*X_localized)(200)*(*X_localized)(200) + (*X_localized)(201) - 5.;
 
  for (std::size_t i=0; i<constraint_values.size(); i++)
    if ((C_ineq.first_local_index() <= i) && (i < C_ineq.last_local_index()))
      C_ineq.set(i, constraint_values[i]);
}

References libMesh::NumericVector< T >::build(), libMesh::ParallelObject::comm(), libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::last_local_index(), libMesh::NumericVector< T >::localize(), libMesh::SERIAL, libMesh::NumericVector< T >::set(), libMesh::NumericVector< T >::size(), and libMesh::NumericVector< T >::zero().

inequality_constraints_jacobian()

void AssembleOptimization::inequality_constraints_jacobian ( const NumericVector< Number > &  X, SparseMatrix< Number > &  C_ineq_jac, OptimizationSystem &  sys  )

virtual

Evaluate the inequality constraints Jacobian.

Implements libMesh::OptimizationSystem::ComputeInequalityConstraintsJacobian.

Definition at line 354 of file optimization_ex2.C.

{
  C_ineq_jac.zero();
 
  std::unique_ptr<NumericVector<Number>> X_localized =
    NumericVector<Number>::build(X.comm());
  X_localized->init(X.size(), false, SERIAL);
  X.localize(*X_localized);
 
  std::vector<std::vector<Number>> constraint_jac_values(1);
  std::vector<std::vector<dof_id_type>> constraint_jac_indices(1);
 
  constraint_jac_values[0].resize(2);
  constraint_jac_indices[0].resize(2);
  constraint_jac_values[0][0] = 2.* (*X_localized)(200);
  constraint_jac_values[0][1] = 1.;
  constraint_jac_indices[0][0] = 200;
  constraint_jac_indices[0][1] = 201;
 
  for (std::size_t i=0; i<constraint_jac_values.size(); i++)
    for (std::size_t j=0; j<constraint_jac_values[i].size(); j++)
      if ((sys.C_ineq->first_local_index() <= i) &&
          (i < sys.C_ineq->last_local_index()))
        {
          dof_id_type col_index = constraint_jac_indices[i][j];
          Number value = constraint_jac_values[i][j];
          C_ineq_jac.set(i, col_index, value);
        }
 
}

References libMesh::NumericVector< T >::build(), libMesh::OptimizationSystem::C_ineq, libMesh::ParallelObject::comm(), libMesh::NumericVector< T >::localize(), libMesh::SERIAL, libMesh::SparseMatrix< T >::set(), libMesh::NumericVector< T >::size(), value, and libMesh::SparseMatrix< T >::zero().

lower_and_upper_bounds()

void AssembleOptimization::lower_and_upper_bounds ( OptimizationSystem &  sys )

virtual

Evaluate the lower and upper bounds vectors.

Implements libMesh::OptimizationSystem::ComputeLowerAndUpperBounds.

Definition at line 387 of file optimization_ex2.C.

{
  for (unsigned int i=sys.get_dof_map().first_dof(); i<sys.get_dof_map().end_dof(); i++)
    {
      sys.get_vector("lower_bounds").set(i, -2.);
      sys.get_vector("upper_bounds").set(i, 2.);
    }
}

References libMesh::DofMapBase::end_dof(), libMesh::DofMapBase::first_dof(), libMesh::System::get_dof_map(), libMesh::System::get_vector(), and libMesh::NumericVector< T >::set().

objective() [1/2]

Number AssembleOptimization::objective ( const NumericVector< Number > &  soln, OptimizationSystem &    )

virtual

Evaluate the objective function.

Implements libMesh::OptimizationSystem::ComputeObjective.

Definition at line 203 of file optimization_ex1.C.

{
  std::unique_ptr<NumericVector<Number>> AxU = soln.zero_clone();
 
  A_matrix->vector_mult(*AxU, soln);
 
  Number
    UTxAxU = AxU->dot(soln),
    UTxF   = F_vector->dot(soln);
 
  return 0.5 * UTxAxU - UTxF;
}

References A_matrix, libMesh::NumericVector< T >::dot(), F_vector, libMesh::SparseMatrix< T >::vector_mult(), and libMesh::NumericVector< T >::zero_clone().

objective() [2/2]

virtual Number AssembleOptimization::objective ( const NumericVector< Number > &  soln, OptimizationSystem &    )

virtual

Evaluate the objective function.

Implements libMesh::OptimizationSystem::ComputeObjective.

Member Data Documentation

_sys

OptimizationSystem & AssembleOptimization::_sys

private

Keep a reference to the OptimizationSystem.

Definition at line 77 of file optimization_ex1.C.

Referenced by assemble_A_and_F().

A_matrix

SparseMatrix< Number > * AssembleOptimization::A_matrix

Sparse matrix for storing the matrix A.

We use this to facilitate computation of objective, gradient and hessian.

Definition at line 118 of file optimization_ex1.C.

Referenced by assemble_A_and_F(), gradient(), hessian(), main(), and objective().

F_vector

NumericVector< Number > * AssembleOptimization::F_vector

Vector for storing F.

We use this to facilitate computation of objective, gradient and hessian.

Definition at line 124 of file optimization_ex1.C.

Referenced by assemble_A_and_F(), gradient(), main(), and objective().

---

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

• examples/optimization/optimization_ex1/optimization_ex1.C
• examples/optimization/optimization_ex2/optimization_ex2.C