ROL_SchurComplement.hpp

Go to the documentation of this file.

// @HEADER
// *****************************************************************************
//               Rapid Optimization Library (ROL) Package
//
// Copyright 2014 NTESS and the ROL contributors.
// SPDX-License-Identifier: BSD-3-Clause
// *****************************************************************************
// @HEADER
 
#ifndef ROL_SCHURCOMPLEMENT_H
#define ROL_SCHURCOMPLEMENT_H
 
#include "ROL_BlockOperator2UnitUpper.hpp"
#include "ROL_BlockOperator2UnitLower.hpp"
#include "ROL_BlockOperator2Diagonal.hpp"
 
namespace ROL {

 
template class<Real> 
class SchurComplement {
 
  typedef Vector<Real>               V;
  typedef PartitionedVector<Real>    PV;
  typedef LinearOperator<Real>       OP;
  typedef BlockOperator2UnitUpper    UPPER;
  typedef BlockOperator2UnitLower    LOWER;
   
private:
 
  ROL::Ptr<OP> A_, B_, C_, D_;
 
  ROL::Ptr<OP> L_,U_;   
  ROL::Ptr<V>  scratch1_;
 
  
 
public:
 
  SchurComplement( ROL::Ptr<OP> &A, ROL::Ptr<OP> &B,  
                   ROL::Ptr<OP> &C, ROL::Ptr<OP> &D, 
                   ROL::Ptr<V> &scratch1 ) : 
                     A_(A), B_(B), C_(C), D_(D), scratch1_(scratch1) {
 
    U_ = ROL::makePtr<UPPER>(B_);
    L_ = ROL::makePtr<LOWER>(C_);
 
  }
 
 
  SchurComplement( BlockOperator2<Real> &op, ROL::Ptr<Vector<Real> > &scratch1 ) :
    scratch1_(scratch1) {}
 
  
 
  A_ = op.getOperator(0,0);
  B_ = op.getOperator(0,1);
  C_ = op.getOperator(1,0);
  D_ = op.getOperator(1,1);
 
  U_ = ROL::makePtr<UPPER>(B_);
  L_ = ROL::makePtr<LOWER>(C_);
 
  void applyLower( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) { 
    L_->apply(Hv,v,tol);
  }
 
  void applyLowerInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) { 
    L_->applyInverse(Hv,v,tol);
  }
 
  void applyUpper( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) { 
    U_->apply(Hv,v,tol);
  }
 
  void applyUpperInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) { 
    U_->applyInverse(Hv,v,tol);
  }
 
  ROL::Ptr<OP> getS11( void ) { 
    return ROL::makePtr<BlockOperator2Determinant<Real>>(A_,B_,C_,D_,scratch1_);
  }
 
  void solve2( Vector<Real> &Hv2, const Vector<Real> &v2, Real &tol ) {
    D_->applyInverse(Hv2,v2,tol);
  }
 
}; // class SchurComplement
 
 
} // namesapce ROL
 
#endif // ROL_SCHURCOMPLEMENT_H