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. More...
void	assemble_A_and_F ()
	The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F. More...
virtual Number	objective (const NumericVector< Number > &soln, OptimizationSystem &)
	Evaluate the objective function. More...
virtual void	gradient (const NumericVector< Number > &soln, NumericVector< Number > &grad_f, OptimizationSystem &)
	Evaluate the gradient. More...
virtual void	hessian (const NumericVector< Number > &soln, SparseMatrix< Number > &H_f, OptimizationSystem &)
	Evaluate the Hessian. More...
	AssembleOptimization (OptimizationSystem &sys_in)
	Constructor. More...
void	assemble_A_and_F ()
	The optimization problem we consider here is: min_U 0.5*U^T A U - U^T F. More...
virtual Number	objective (const NumericVector< Number > &soln, OptimizationSystem &)
	Evaluate the objective function. More...
virtual void	gradient (const NumericVector< Number > &soln, NumericVector< Number > &grad_f, OptimizationSystem &)
	Evaluate the gradient. More...
virtual void	hessian (const NumericVector< Number > &soln, SparseMatrix< Number > &H_f, OptimizationSystem &)
	Evaluate the Hessian. More...
virtual void	equality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_eq, OptimizationSystem &)
	Evaluate the equality constraints. More...
virtual void	equality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_eq_jac, OptimizationSystem &)
	Evaluate the equality constraints Jacobian. More...
virtual void	inequality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_ineq, OptimizationSystem &)
	Evaluate the inequality constraints. More...
virtual void	inequality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, OptimizationSystem &)
	Evaluate the inequality constraints Jacobian. More...
virtual void	lower_and_upper_bounds (OptimizationSystem &sys)
	Evaluate the lower and upper bounds vectors. More...
virtual Number	objective (const NumericVector< Number > &X, sys_type &S)=0
	This function will be called to compute the objective function to be minimized, and must be implemented by the user in a derived class. More...
virtual void	gradient (const NumericVector< Number > &X, NumericVector< Number > &grad_f, sys_type &S)=0
	This function will be called to compute the gradient of the objective function, and must be implemented by the user in a derived class. More...
virtual void	hessian (const NumericVector< Number > &X, SparseMatrix< Number > &H_f, sys_type &S)=0
	This function will be called to compute the Hessian of the objective function, and must be implemented by the user in a derived class. More...
virtual void	equality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_eq, sys_type &S)=0
	This function will be called to evaluate the equality constraints vector C_eq(X). More...
virtual void	equality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_eq_jac, sys_type &S)=0
	This function will be called to evaluate the Jacobian of C_eq(X). More...
virtual void	inequality_constraints (const NumericVector< Number > &X, NumericVector< Number > &C_ineq, sys_type &S)=0
	This function will be called to evaluate the equality constraints vector C_ineq(X). More...
virtual void	inequality_constraints_jacobian (const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, sys_type &S)=0
	This function will be called to evaluate the Jacobian of C_ineq(X). More...
virtual void	lower_and_upper_bounds (sys_type &S)=0
	This function should update the following two vectors: this->get_vector("lower_bounds"), this->get_vector("upper_bounds"). More...

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

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

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 66 of file optimization_ex1.C.

Constructor & Destructor Documentation

AssembleOptimization() [1/2]

AssembleOptimization::AssembleOptimization

(

OptimizationSystem &

sys_in

)

Constructor.

Definition at line 126 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 130 of file optimization_ex1.C.

References _sys, A_matrix, libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), 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::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().

{
  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();
}

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() [1/2]

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

virtual

Evaluate the equality constraints.

Definition at line 276 of file optimization_ex2.C.

References libMesh::ParallelObject::comm(), libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::init(), 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().

{
  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]);
}

equality_constraints() [2/2]

virtual void libMesh::OptimizationSystem::ComputeEqualityConstraints::equality_constraints ( const NumericVector< Number > &  X, NumericVector< Number > &  C_eq, sys_type &  S  )

pure virtualinherited

This function will be called to evaluate the equality constraints vector C_eq(X).

This will impose the constraints C_eq(X) = 0.

Referenced by libMesh::__libmesh_nlopt_equality_constraints(), and libMesh::__libmesh_tao_equality_constraints().

equality_constraints_jacobian() [1/2]

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

virtual

Evaluate the equality constraints Jacobian.

Definition at line 298 of file optimization_ex2.C.

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

{
  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);
        }
}

equality_constraints_jacobian() [2/2]

virtual void libMesh::OptimizationSystem::ComputeEqualityConstraintsJacobian::equality_constraints_jacobian ( const NumericVector< Number > &  X, SparseMatrix< Number > &  C_eq_jac, sys_type &  S  )

pure virtualinherited

This function will be called to evaluate the Jacobian of C_eq(X).

Referenced by libMesh::__libmesh_nlopt_equality_constraints(), and libMesh::__libmesh_tao_equality_constraints_jacobian().

gradient() [1/3]

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

virtual

Evaluate the gradient.

Definition at line 216 of file optimization_ex1.C.

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

{
  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);
}

gradient() [2/3]

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

virtual

Evaluate the gradient.

gradient() [3/3]

virtual void libMesh::OptimizationSystem::ComputeGradient::gradient ( const NumericVector< Number > &  X, NumericVector< Number > &  grad_f, sys_type &  S  )

pure virtualinherited

This function will be called to compute the gradient of the objective function, and must be implemented by the user in a derived class.

Set grad_f to be the gradient at the iterate X.

Referenced by libMesh::__libmesh_nlopt_objective(), and libMesh::__libmesh_tao_gradient().

hessian() [1/3]

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

virtual

Evaluate the Hessian.

hessian() [2/3]

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

virtual

Evaluate the Hessian.

Definition at line 230 of file optimization_ex1.C.

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

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

hessian() [3/3]

virtual void libMesh::OptimizationSystem::ComputeHessian::hessian ( const NumericVector< Number > &  X, SparseMatrix< Number > &  H_f, sys_type &  S  )

pure virtualinherited

This function will be called to compute the Hessian of the objective function, and must be implemented by the user in a derived class.

Set H_f to be the gradient at the iterate X.

Referenced by libMesh::__libmesh_tao_hessian().

inequality_constraints() [1/2]

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

virtual

Evaluate the inequality constraints.

Definition at line 335 of file optimization_ex2.C.

References libMesh::ParallelObject::comm(), libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::init(), 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().

{
  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]);
}

inequality_constraints() [2/2]

virtual void libMesh::OptimizationSystem::ComputeInequalityConstraints::inequality_constraints ( const NumericVector< Number > &  X, NumericVector< Number > &  C_ineq, sys_type &  S  )

pure virtualinherited

This function will be called to evaluate the equality constraints vector C_ineq(X).

This will impose the constraints C_ineq(X) >= 0.

Referenced by libMesh::__libmesh_nlopt_inequality_constraints(), and libMesh::__libmesh_tao_inequality_constraints().

inequality_constraints_jacobian() [1/2]

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

virtual

Evaluate the inequality constraints Jacobian.

Definition at line 354 of file optimization_ex2.C.

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

{
  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);
        }

}

inequality_constraints_jacobian() [2/2]

virtual void libMesh::OptimizationSystem::ComputeInequalityConstraintsJacobian::inequality_constraints_jacobian ( const NumericVector< Number > &  X, SparseMatrix< Number > &  C_ineq_jac, sys_type &  S  )

pure virtualinherited

This function will be called to evaluate the Jacobian of C_ineq(X).

Referenced by libMesh::__libmesh_nlopt_inequality_constraints(), and libMesh::__libmesh_tao_inequality_constraints_jacobian().

lower_and_upper_bounds() [1/2]

void AssembleOptimization::lower_and_upper_bounds ( OptimizationSystem &  sys )

virtual

Evaluate the lower and upper bounds vectors.

Definition at line 387 of file optimization_ex2.C.

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

{
  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.);
    }
}

lower_and_upper_bounds() [2/2]

virtual void libMesh::OptimizationSystem::ComputeLowerAndUpperBounds::lower_and_upper_bounds ( sys_type &  S )

pure virtualinherited

This function should update the following two vectors: this->get_vector("lower_bounds"), this->get_vector("upper_bounds").

objective() [1/3]

virtual Number libMesh::OptimizationSystem::ComputeObjective::objective ( const NumericVector< Number > &  X, sys_type &  S  )

pure virtualinherited

This function will be called to compute the objective function to be minimized, and must be implemented by the user in a derived class.

Returns

The value of the objective function at iterate X.

Referenced by libMesh::__libmesh_nlopt_objective(), and libMesh::__libmesh_tao_objective().

objective() [2/3]

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

virtual

Evaluate the objective function.

objective() [3/3]

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

virtual

Evaluate the objective function.

Definition at line 202 of file optimization_ex1.C.

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

{
  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;
}

Member Data Documentation

_sys

OptimizationSystem & AssembleOptimization::_sys

private

Keep a reference to the OptimizationSystem.

Definition at line 76 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 117 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 123 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