ROL
|
#include <example_01.hpp>
Public Member Functions | |
Objective_GrossPitaevskii (const Real &g, const Vector< Real > &V) | |
Real | value (const Vector< Real > &psi, Real &tol) |
Evaluate \(J[\psi]\). More... | |
void | gradient (Vector< Real > &g, const Vector< Real > &psi, Real &tol) |
Evaluate \(\nabla J[\psi]\). More... | |
void | hessVec (Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol) |
Evaluate \(\nabla^2 J[\psi] v\). More... | |
Objective_GrossPitaevskii (const Real &g, const Vector< Real > &V, Teuchos::RCP< FiniteDifference< Real > > fd) | |
Real | value (const Vector< Real > &psi, Real &tol) |
Evaluate \(J[\psi]\). More... | |
void | gradient (Vector< Real > &g, const Vector< Real > &psi, Real &tol) |
Evaluate \(\nabla J[\psi]\). More... | |
void | hessVec (Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol) |
Evaluate \(\nabla^2 J[\psi] v\). More... | |
Objective_GrossPitaevskii (const int ni, const Real gnl, Teuchos::RCP< NodalBasis< Real > > nb, Teuchos::RCP< InnerProductMatrix< Real > > kinetic, Teuchos::RCP< InnerProductMatrix< Real > > potential, Teuchos::RCP< InnerProductMatrix< Real > > nonlinear) | |
Real | value (const Vector< Real > &psi, Real &tol) |
Compute value. More... | |
void | gradient (Vector< Real > &g, const Vector< Real > &psi, Real &tol) |
Compute gradient. More... | |
void | hessVec (Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol) |
Apply Hessian approximation to vector. More... | |
![]() | |
virtual | ~Objective () |
virtual void | update (const Vector< Real > &x, bool flag=true, int iter=-1) |
Update objective function. More... | |
virtual Real | dirDeriv (const Vector< Real > &x, const Vector< Real > &d, Real &tol) |
Compute directional derivative. More... | |
virtual void | invHessVec (Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol) |
Apply inverse Hessian approximation to vector. More... | |
virtual void | precond (Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol) |
Apply preconditioner to vector. More... | |
virtual std::vector < std::vector< Real > > | checkGradient (const Vector< Real > &x, const Vector< Real > &d, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1) |
Finite-difference gradient check. More... | |
virtual std::vector < std::vector< Real > > | checkGradient (const Vector< Real > &x, const Vector< Real > &g, const Vector< Real > &d, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1) |
Finite-difference gradient check. More... | |
virtual std::vector < std::vector< Real > > | checkGradient (const Vector< Real > &x, const Vector< Real > &d, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1) |
Finite-difference gradient check with specified step sizes. More... | |
virtual std::vector < std::vector< Real > > | checkGradient (const Vector< Real > &x, const Vector< Real > &g, const Vector< Real > &d, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1) |
Finite-difference gradient check with specified step sizes. More... | |
virtual std::vector < std::vector< Real > > | checkHessVec (const Vector< Real > &x, const Vector< Real > &v, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1) |
Finite-difference Hessian-applied-to-vector check. More... | |
virtual std::vector < std::vector< Real > > | checkHessVec (const Vector< Real > &x, const Vector< Real > &hv, const Vector< Real > &v, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1) |
Finite-difference Hessian-applied-to-vector check. More... | |
virtual std::vector < std::vector< Real > > | checkHessVec (const Vector< Real > &x, const Vector< Real > &v, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1) |
Finite-difference Hessian-applied-to-vector check with specified step sizes. More... | |
virtual std::vector < std::vector< Real > > | checkHessVec (const Vector< Real > &x, const Vector< Real > &hv, const Vector< Real > &v, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1) |
Finite-difference Hessian-applied-to-vector check with specified step sizes. More... | |
virtual std::vector< Real > | checkHessSym (const Vector< Real > &x, const Vector< Real > &v, const Vector< Real > &w, const bool printToStream=true, std::ostream &outStream=std::cout) |
Hessian symmetry check. More... | |
virtual std::vector< Real > | checkHessSym (const Vector< Real > &x, const Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &w, const bool printToStream=true, std::ostream &outStream=std::cout) |
Hessian symmetry check. More... | |
Private Types | |
typedef std::vector< Real > | svec |
typedef StdVector< Real > | rvec |
typedef Teuchos::RCP< const svec > | pcsv |
typedef Teuchos::RCP< svec > | psv |
Private Member Functions | |
void | applyK (const Vector< Real > &v, Vector< Real > &kv) |
Apply finite difference operator. More... | |
void | applyK (const Vector< Real > &v, Vector< Real > &kv) |
Apply finite difference operator. More... | |
void | updateNonlinear (Teuchos::RCP< const std::vector< Real > > psip) |
Private Attributes | |
Real | g_ |
int | nx_ |
Real | dx_ |
pcsv | Vp_ |
Teuchos::RCP< const std::vector< Real > > | Vp_ |
Teuchos::RCP< FiniteDifference < Real > > | fd_ |
const int | ni_ |
const Real | gnl_ |
Teuchos::RCP< NodalBasis< Real > > | nb_ |
Teuchos::RCP < InnerProductMatrix< Real > > | kinetic_ |
Teuchos::RCP < InnerProductMatrix< Real > > | potential_ |
Teuchos::RCP < InnerProductMatrix< Real > > | nonlinear_ |
Objective Function Class
Definition at line 86 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 89 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 92 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 95 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 98 of file gross-pitaevskii/example_01.hpp.
|
inline |
Definition at line 140 of file gross-pitaevskii/example_01.hpp.
|
inline |
Definition at line 466 of file gross-pitaevskii/example_02.hpp.
|
inline |
Definition at line 455 of file gross-pitaevskii/example_03.hpp.
|
inlineprivate |
Apply finite difference operator.
Compute \(K\psi\), where \(K\) is the finite difference approximation of \(-D_x^2\)
Definition at line 118 of file gross-pitaevskii/example_01.hpp.
|
inlinevirtual |
Evaluate \(J[\psi]\).
\[ J[\psi]=\frac{1}{2} \int\limits_0^1 |\psi'|^2 + V(x)|\psi|^2+g|\psi|^4\,\mathrm{d}x \]
where the integral is approximated with the trapezoidal rule and the derivative is approximated using finite differences
Implements ROL::Objective< Real >.
Definition at line 152 of file gross-pitaevskii/example_01.hpp.
|
inlinevirtual |
Evaluate \(\nabla J[\psi]\).
\[ \nabla J[\psi] = -\psi'' + V(x)\psi+2g|\psi|^3 \]
Reimplemented from ROL::Objective< Real >.
Definition at line 177 of file gross-pitaevskii/example_01.hpp.
|
inlinevirtual |
Evaluate \(\nabla^2 J[\psi] v\).
\[ \nabla^2 J[\psi]v = -v'' + V(x)v+6g|\psi|^2 v \]
Reimplemented from ROL::Objective< Real >.
Definition at line 201 of file gross-pitaevskii/example_01.hpp.
|
inlineprivate |
Apply finite difference operator.
Compute \(K\psi\), where \(K\) is the finite difference approximation of \(-D_x^2\)
Definition at line 444 of file gross-pitaevskii/example_02.hpp.
|
inlinevirtual |
Evaluate \(J[\psi]\).
\[ J[\psi]=\frac{1}{2} \int\limits_0^1 |\psi'|^2 + V(x)|\psi|^2+g|\psi|^4\,\mathrm{d}x \]
where the integral is approximated with the trapezoidal rule and the derivative is approximated using finite differences
Implements ROL::Objective< Real >.
Definition at line 478 of file gross-pitaevskii/example_02.hpp.
|
inlinevirtual |
Evaluate \(\nabla J[\psi]\).
\[ \nabla J[\psi] = -\psi'' + V(x)\psi+2g|\psi|^3 \]
Reimplemented from ROL::Objective< Real >.
Definition at line 506 of file gross-pitaevskii/example_02.hpp.
|
inlinevirtual |
Evaluate \(\nabla^2 J[\psi] v\).
\[ \nabla^2 J[\psi]v = -v'' + V(x)v+6g|\psi|^2 v \]
Reimplemented from ROL::Objective< Real >.
Definition at line 532 of file gross-pitaevskii/example_02.hpp.
|
inlineprivate |
Definition at line 439 of file gross-pitaevskii/example_03.hpp.
|
inlinevirtual |
Compute value.
This function returns the objective function value.
[in] | x | is the current iterate. |
[in] | tol | is a tolerance for inexact objective function computation. |
Implements ROL::Objective< Real >.
Definition at line 464 of file gross-pitaevskii/example_03.hpp.
|
inlinevirtual |
Compute gradient.
This function returns the objective function gradient.
[out] | g | is the gradient. |
[in] | x | is the current iterate. |
[in] | tol | is a tolerance for inexact objective function computation. |
The default implementation is a finite-difference approximation based on the function value. This requires the definition of a basis \(\{\phi_i\}\) for the optimization vectors x and the definition of a basis \(\{\psi_j\}\) for the dual optimization vectors (gradient vectors g). The bases must be related through the Riesz map, i.e., \( R \{\phi_i\} = \{\psi_j\}\), and this must be reflected in the implementation of the ROL::Vector::dual() method.
Reimplemented from ROL::Objective< Real >.
Definition at line 481 of file gross-pitaevskii/example_03.hpp.
|
inlinevirtual |
Apply Hessian approximation to vector.
This function applies the Hessian of the objective function to the vector \(v\).
[out] | hv | is the the action of the Hessian on \(v\). |
[in] | v | is the direction vector. |
[in] | x | is the current iterate. |
[in] | tol | is a tolerance for inexact objective function computation. |
Reimplemented from ROL::Objective< Real >.
Definition at line 498 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 104 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 107 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 110 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 113 of file gross-pitaevskii/example_01.hpp.
|
private |
Definition at line 437 of file gross-pitaevskii/example_02.hpp.
|
private |
Definition at line 439 of file gross-pitaevskii/example_02.hpp.
|
private |
Definition at line 431 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 432 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 433 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 435 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 436 of file gross-pitaevskii/example_03.hpp.
|
private |
Definition at line 437 of file gross-pitaevskii/example_03.hpp.