libMesh
Public Member Functions | Public Attributes | Private Attributes | List of all members
AssembleOptimization Class Referenceabstract

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

Inheritance diagram for AssembleOptimization:
[legend]

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::AssembleOptimization ( OptimizationSystem sys_in)

Constructor.

Definition at line 126 of file optimization_ex1.C.

126  :
127  _sys(sys_in)
128 {}
OptimizationSystem & _sys
Keep a reference to the OptimizationSystem.
AssembleOptimization::AssembleOptimization ( OptimizationSystem sys_in)

Constructor.

Member Function Documentation

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.

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::MeshBase::active_local_element_ptr_range(), 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::SparseMatrix< T >::zero(), and libMesh::NumericVector< T >::zero().

Referenced by main().

131 {
132  A_matrix->zero();
133  F_vector->zero();
134 
135  const MeshBase & mesh = _sys.get_mesh();
136 
137  const unsigned int dim = mesh.mesh_dimension();
138  const unsigned int u_var = _sys.variable_number ("u");
139 
140  const DofMap & dof_map = _sys.get_dof_map();
141  FEType fe_type = dof_map.variable_type(u_var);
142  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
143  QGauss qrule (dim, fe_type.default_quadrature_order());
144  fe->attach_quadrature_rule (&qrule);
145 
146  const std::vector<Real> & JxW = fe->get_JxW();
147  const std::vector<std::vector<Real>> & phi = fe->get_phi();
148  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
149 
150  std::vector<dof_id_type> dof_indices;
151 
154 
155  for (const auto & elem : mesh.active_local_element_ptr_range())
156  {
157  dof_map.dof_indices (elem, dof_indices);
158 
159  const unsigned int n_dofs = dof_indices.size();
160 
161  fe->reinit (elem);
162 
163  Ke.resize (n_dofs, n_dofs);
164  Fe.resize (n_dofs);
165 
166  for (unsigned int qp=0; qp<qrule.n_points(); qp++)
167  {
168  for (unsigned int dof_i=0; dof_i<n_dofs; dof_i++)
169  {
170  for (unsigned int dof_j=0; dof_j<n_dofs; dof_j++)
171  {
172  Ke(dof_i, dof_j) += JxW[qp] * (dphi[dof_j][qp]* dphi[dof_i][qp]);
173  }
174  Fe(dof_i) += JxW[qp] * phi[dof_i][qp];
175  }
176  }
177 
178  // This will zero off-diagonal entries of Ke corresponding to
179  // Dirichlet dofs.
180  dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
181 
182  // We want the diagonal of constrained dofs to be zero too
183  for (unsigned int local_dof_index=0; local_dof_index<n_dofs; local_dof_index++)
184  {
185  dof_id_type global_dof_index = dof_indices[local_dof_index];
186  if (dof_map.is_constrained_dof(global_dof_index))
187  {
188  Ke(local_dof_index, local_dof_index) = 0.;
189  }
190  }
191 
192  A_matrix->add_matrix (Ke, dof_indices);
193  F_vector->add_vector (Fe, dof_indices);
194  }
195 
196  A_matrix->close();
197  F_vector->close();
198 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
OptimizationSystem & _sys
Keep a reference to the OptimizationSystem.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
unsigned int dim
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
MeshBase & mesh
NumericVector< Number > * F_vector
Vector for storing F.
This is the MeshBase class.
Definition: mesh_base.h:68
virtual void zero()=0
Set all entries to zero.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
SparseMatrix< Number > * A_matrix
Sparse matrix for storing the matrix A.
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1263
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
const MeshBase & get_mesh() const
Definition: system.h:2014
bool is_constrained_dof(const dof_id_type dof) const
Definition: dof_map.h:1734
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual void zero()=0
Set all entries to 0.
Order default_quadrature_order() const
Definition: fe_type.h:332
void constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix and vector.
Definition: dof_map.h:1806
virtual void close()=0
Calls the NumericVector&#39;s internal assembly routines, ensuring that the values are consistent across ...
virtual void close()=0
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
This class implements specific orders of Gauss quadrature.
uint8_t dof_id_type
Definition: id_types.h:64
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
void AssembleOptimization::equality_constraints ( const NumericVector< Number > &  X,
NumericVector< Number > &  C_eq,
OptimizationSystem  
)
virtual

Evaluate the equality constraints.

Definition at line 275 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().

278 {
279  C_eq.zero();
280 
281  UniquePtr<NumericVector<Number>> X_localized =
283  X_localized->init(X.size(), false, SERIAL);
284  X.localize(*X_localized);
285 
286  std::vector<Number> constraint_values(3);
287  constraint_values[0] = (*X_localized)(17);
288  constraint_values[1] = (*X_localized)(23);
289  constraint_values[2] = (*X_localized)(98) + (*X_localized)(185);
290 
291  for (std::size_t i=0; i<constraint_values.size(); i++)
292  if ((C_eq.first_local_index() <= i) &&
293  (i < C_eq.last_local_index()))
294  C_eq.set(i, constraint_values[i]);
295 }
virtual numeric_index_type size() const =0
virtual numeric_index_type last_local_index() const =0
virtual void init(const numeric_index_type, const numeric_index_type, const bool=false, const ParallelType=AUTOMATIC)=0
Change the dimension of the vector to N.
Numeric vector.
Definition: dof_map.h:66
virtual void zero()=0
Set all entries to zero.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const Parallel::Communicator & comm() const
virtual numeric_index_type first_local_index() const =0
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.
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().

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 297 of file optimization_ex2.C.

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

300 {
301  C_eq_jac.zero();
302 
303  std::vector<std::vector<Number>> constraint_jac_values(3);
304  std::vector<std::vector<dof_id_type>> constraint_jac_indices(3);
305 
306  constraint_jac_values[0].resize(1);
307  constraint_jac_indices[0].resize(1);
308  constraint_jac_values[0][0] = 1.;
309  constraint_jac_indices[0][0] = 17;
310 
311  constraint_jac_values[1].resize(1);
312  constraint_jac_indices[1].resize(1);
313  constraint_jac_values[1][0] = 1.;
314  constraint_jac_indices[1][0] = 23;
315 
316  constraint_jac_values[2].resize(2);
317  constraint_jac_indices[2].resize(2);
318  constraint_jac_values[2][0] = 1.;
319  constraint_jac_values[2][1] = 1.;
320  constraint_jac_indices[2][0] = 98;
321  constraint_jac_indices[2][1] = 185;
322 
323  for (std::size_t i=0; i<constraint_jac_values.size(); i++)
324  for (std::size_t j=0; j<constraint_jac_values[i].size(); j++)
325  if ((sys.C_eq->first_local_index() <= i) &&
326  (i < sys.C_eq->last_local_index()))
327  {
328  dof_id_type col_index = constraint_jac_indices[i][j];
329  Number value = constraint_jac_values[i][j];
330  C_eq_jac.set(i, col_index, value);
331  }
332 }
virtual void set(const numeric_index_type i, const numeric_index_type j, const T value)=0
Set the element (i,j) to value.
virtual void zero()=0
Set all entries to 0.
static const bool value
Definition: xdr_io.C:108
UniquePtr< NumericVector< Number > > C_eq
The vector that stores equality constraints.
uint8_t dof_id_type
Definition: id_types.h:64
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().

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

Evaluate the gradient.

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

Evaluate the gradient.

Definition at line 214 of file optimization_ex1.C.

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

217 {
218  grad_f.zero();
219 
220  // Since we've enforced constraints on soln, A and F,
221  // this automatically sets grad_f to zero for constrained
222  // dofs.
223  A_matrix->vector_mult(grad_f, soln);
224  grad_f.add(-1, *F_vector);
225 }
NumericVector< Number > * F_vector
Vector for storing F.
void vector_mult(NumericVector< T > &dest, const NumericVector< T > &arg) const
Multiplies the matrix by the NumericVector arg and stores the result in NumericVector dest...
virtual void zero()=0
Set all entries to zero.
SparseMatrix< Number > * A_matrix
Sparse matrix for storing the matrix A.
virtual void add(const numeric_index_type i, const T value)=0
Adds value to each entry of the vector.
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().

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

Evaluate the Hessian.

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

Evaluate the Hessian.

Definition at line 228 of file optimization_ex1.C.

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

231 {
232  H_f.zero();
233  H_f.add(1., *A_matrix);
234 }
virtual void add(const numeric_index_type i, const numeric_index_type j, const T value)=0
Add value to the element (i,j).
SparseMatrix< Number > * A_matrix
Sparse matrix for storing the matrix A.
virtual void zero()=0
Set all entries to 0.
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().

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

Evaluate the inequality constraints.

Definition at line 334 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().

337 {
338  C_ineq.zero();
339 
340  UniquePtr<NumericVector<Number>> X_localized =
342  X_localized->init(X.size(), false, SERIAL);
343  X.localize(*X_localized);
344 
345  std::vector<Number> constraint_values(1);
346  constraint_values[0] = (*X_localized)(200)*(*X_localized)(200) + (*X_localized)(201) - 5.;
347 
348  for (std::size_t i=0; i<constraint_values.size(); i++)
349  if ((C_ineq.first_local_index() <= i) && (i < C_ineq.last_local_index()))
350  C_ineq.set(i, constraint_values[i]);
351 }
virtual numeric_index_type size() const =0
virtual numeric_index_type last_local_index() const =0
virtual void init(const numeric_index_type, const numeric_index_type, const bool=false, const ParallelType=AUTOMATIC)=0
Change the dimension of the vector to N.
Numeric vector.
Definition: dof_map.h:66
virtual void zero()=0
Set all entries to zero.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const Parallel::Communicator & comm() const
virtual numeric_index_type first_local_index() const =0
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.
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().

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 353 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().

356 {
357  C_ineq_jac.zero();
358 
359  UniquePtr<NumericVector<Number>> X_localized =
361  X_localized->init(X.size(), false, SERIAL);
362  X.localize(*X_localized);
363 
364  std::vector<std::vector<Number>> constraint_jac_values(1);
365  std::vector<std::vector<dof_id_type>> constraint_jac_indices(1);
366 
367  constraint_jac_values[0].resize(2);
368  constraint_jac_indices[0].resize(2);
369  constraint_jac_values[0][0] = 2.* (*X_localized)(200);
370  constraint_jac_values[0][1] = 1.;
371  constraint_jac_indices[0][0] = 200;
372  constraint_jac_indices[0][1] = 201;
373 
374  for (std::size_t i=0; i<constraint_jac_values.size(); i++)
375  for (std::size_t j=0; j<constraint_jac_values[i].size(); j++)
376  if ((sys.C_ineq->first_local_index() <= i) &&
377  (i < sys.C_ineq->last_local_index()))
378  {
379  dof_id_type col_index = constraint_jac_indices[i][j];
380  Number value = constraint_jac_values[i][j];
381  C_ineq_jac.set(i, col_index, value);
382  }
383 
384 }
virtual numeric_index_type size() const =0
virtual void init(const numeric_index_type, const numeric_index_type, const bool=false, const ParallelType=AUTOMATIC)=0
Change the dimension of the vector to N.
Numeric vector.
Definition: dof_map.h:66
virtual void set(const numeric_index_type i, const numeric_index_type j, const T value)=0
Set the element (i,j) to value.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual void zero()=0
Set all entries to 0.
const Parallel::Communicator & comm() const
static const bool value
Definition: xdr_io.C:108
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.
uint8_t dof_id_type
Definition: id_types.h:64
UniquePtr< NumericVector< Number > > C_ineq
The vector that stores inequality constraints.
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().

void AssembleOptimization::lower_and_upper_bounds ( OptimizationSystem sys)
virtual

Evaluate the lower and upper bounds vectors.

Definition at line 386 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().

387 {
388  for (unsigned int i=sys.get_dof_map().first_dof(); i<sys.get_dof_map().end_dof(); i++)
389  {
390  sys.get_vector("lower_bounds").set(i, -2.);
391  sys.get_vector("upper_bounds").set(i, 2.);
392  }
393 }
dof_id_type first_dof(const processor_id_type proc) const
Definition: dof_map.h:535
dof_id_type end_dof(const processor_id_type proc) const
Definition: dof_map.h:577
const DofMap & get_dof_map() const
Definition: system.h:2030
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:794
virtual void set(const numeric_index_type i, const T value)=0
Sets v(i) = value.
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").

Referenced by libMesh::TaoOptimizationSolver< T >::solve(), and libMesh::NloptOptimizationSolver< T >::solve().

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(), libMesh::__libmesh_tao_objective(), and libMesh::OptimizationSystem::ComputeObjective::~ComputeObjective().

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

Evaluate the objective function.

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

Evaluate the objective function.

Definition at line 200 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().

202 {
204 
205  A_matrix->vector_mult(*AxU, soln);
206 
207  Number
208  UTxAxU = AxU->dot(soln),
209  UTxF = F_vector->dot(soln);
210 
211  return 0.5 * UTxAxU - UTxF;
212 }
NumericVector< Number > * F_vector
Vector for storing F.
void vector_mult(NumericVector< T > &dest, const NumericVector< T > &arg) const
Multiplies the matrix by the NumericVector arg and stores the result in NumericVector dest...
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
SparseMatrix< Number > * A_matrix
Sparse matrix for storing the matrix A.
virtual UniquePtr< NumericVector< T > > zero_clone() const =0
virtual T dot(const NumericVector< T > &v) const =0

Member Data Documentation

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().

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().

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: