NLPInterfacePack_NLPSerialPreprocessExplJac.hpp

// @HEADER
// ***********************************************************************
// 
// Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
//                  Copyright (2003) Sandia Corporation
// 
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
// 
// This library is free software; you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; either version 2.1 of the
// License, or (at your option) any later version.
//  
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//  
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
// USA
// Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
// 
// ***********************************************************************
// @HEADER

#ifndef NLP_SERIAL_PREPROCESS_EXPL_JAC_H
#define NLP_SERIAL_PREPROCESS_EXPL_JAC_H

#include <valarray>

#include "NLPInterfacePack_NLPSerialPreprocess.hpp"
#include "NLPInterfacePack_NLPFirstOrder.hpp"
#include "AbstractLinAlgPack_BasisSystemFactoryStd.hpp"
#include "DenseLinAlgPack_DVectorClass.hpp"
#include "Teuchos_AbstractFactory.hpp"
#include "Teuchos_StandardCompositionMacros.hpp"

namespace NLPInterfacePack {

class NLPSerialPreprocessExplJac
  : virtual public NLPSerialPreprocess
  , virtual public NLPFirstOrder
{
public:

  
  typedef Teuchos::RCP<
    const Teuchos::AbstractFactory<MatrixOp> >    factory_mat_ptr_t;



  STANDARD_COMPOSITION_MEMBERS( BasisSystemFactory, basis_sys_fcty );

  NLPSerialPreprocessExplJac(
    const basis_sys_fcty_ptr_t  &basis_sys_fcty  = Teuchos::rcp(new BasisSystemFactoryStd())
    ,const factory_mat_ptr_t    &factory_Gc_full = Teuchos::null
    );

  void set_factory_Gc_full( const factory_mat_ptr_t &factory_Gc_full );



  void set_options( const options_ptr_t& options );
  const options_ptr_t& get_options() const;
  void initialize(bool test_setup); 
  bool is_initialized() const;


  
  const mat_fcty_ptr_t factory_Gc() const;
  const basis_sys_ptr_t basis_sys() const;
  void set_Gc(MatrixOp* Gc);



  bool get_next_basis(
    Permutation*  P_var,   Range1D* var_dep
    ,Permutation* P_equ,   Range1D* equ_decomp
    );
  void set_basis(
    const Permutation   &P_var,   const Range1D  &var_dep
    ,const Permutation  *P_equ,   const Range1D  *equ_decomp
    );



protected:


  void imp_calc_Gc(
    const Vector& x, bool newx
    ,const FirstOrderInfo& first_order_info
    ) const;
  


  struct FirstOrderExplInfo {
    typedef std::valarray<value_type>    val_t;
    typedef std::valarray<index_type>    ivect_t;
    //
    typedef std::valarray<index_type>    jvect_t;
    FirstOrderExplInfo()
      :Gc_val(NULL), Gc_ivect(NULL), Gc_jvect(NULL)
      ,Gh_val(NULL), Gh_ivect(NULL), Gh_jvect(NULL)
      ,f(NULL)
    {}
    FirstOrderExplInfo(
      index_type* Gc_nz_in, val_t* Gc_val_in, ivect_t* Gc_ivect_in, jvect_t* Gc_jvect_in
      ,index_type* Gh_nz_in, val_t* Gh_val_in, ivect_t* Gh_ivect_in, jvect_t* Gh_jvect_in
      ,const ObjGradInfoSerial& obj_grad
      )
      :Gc_nz(Gc_nz_in), Gc_val(Gc_val_in), Gc_ivect(Gc_ivect_in), Gc_jvect(Gc_jvect_in)
      ,Gh_nz(Gh_nz_in), Gh_val(Gh_val_in), Gh_ivect(Gh_ivect_in), Gh_jvect(Gh_jvect_in)
      ,Gf(obj_grad.Gf), f(obj_grad.f), c(obj_grad.c), h(obj_grad.h)
    {}
    size_type*    Gc_nz;
    val_t*        Gc_val;
    ivect_t*      Gc_ivect;
    jvect_t*      Gc_jvect;
    size_type*    Gh_nz;
    val_t*        Gh_val;
    ivect_t*      Gh_ivect;
    jvect_t*      Gh_jvect;
    DVector*       Gf;
    value_type*   f;
    DVector*       c;
    DVector*       h;
  }; // end struct FirstOrderExplInfo



  virtual size_type imp_Gc_nz_orig() const = 0;

  virtual size_type imp_Gh_nz_orig() const = 0;

  virtual void imp_calc_Gc_orig(
    const DVectorSlice& x_full, bool newx
    , const FirstOrderExplInfo& first_order_expl_info
    ) const = 0;

  virtual void imp_calc_Gh_orig(
    const DVectorSlice& x_full, bool newx
    , const FirstOrderExplInfo& first_order_expl_info
    ) const = 0;



  void assert_initialized() const;

  const FirstOrderExplInfo first_order_expl_info() const;


private:

  // ////////////////////////////////////////
  // Private data members
  
  bool initialized_;              // Flag for if the NLP has has been properly initialized
  bool test_setup_;               // Flag for if to test the setup of things or not
  options_ptr_t options_;         // The options being used

  factory_mat_ptr_t   factory_Gc_full_;
  mat_fcty_ptr_t      factory_Gc_;

  mutable size_type   Gc_nz_orig_;    // Number of nonzeros in the original NLP Gc
  mutable size_type   Gh_nz_orig_;    // Number of nonzeros in the original NLP Gh
  mutable size_type   Gc_nz_full_;    // Number of nonzeros in the full NLP Gc
  mutable size_type   Gh_nz_full_;    // Number of nonzeros in the full NLP Gh
  mutable FirstOrderExplInfo::val_t    Gc_val_orig_;   // Storage for explicit nonzeros of full Gc
  mutable FirstOrderExplInfo::ivect_t  Gc_ivect_orig_;
  mutable FirstOrderExplInfo::jvect_t  Gc_jvect_orig_;
  mutable FirstOrderExplInfo::val_t    Gh_val_orig_;   // Storage for explicit nonzeros of orig Gh
  mutable FirstOrderExplInfo::ivect_t  Gh_ivect_orig_;
  mutable FirstOrderExplInfo::jvect_t  Gh_jvect_orig_;

  mutable bool                         Gc_perm_new_basis_updated_;  // Flag for if a new basis was set!

  // ////////////////////////////
  // Private member functions

  //
  void imp_calc_Gc_or_Gh(
    bool calc_Gc
    ,const Vector& x, bool newx
    ,const FirstOrderInfo& first_order_info
    ) const;

  //
  void imp_fill_jacobian_entries(
    size_type           n             // [in]
    ,size_type          n_full        // [in]
    ,bool               load_struct   // [in] If true, then the structure is loaded also
    ,const index_type   col_offset    // [in] Offset for filled column indexes
    ,const value_type   *val_full     // [in] Values (!=NULL)
    ,const value_type   *val_full_end // [in] Values end (!=NULL)
    ,const index_type   *ivect_full   // [in] Row indexes (!=NULL)
    ,const index_type   *jvect_full   // [in] Column indexes (!=NULL)
    ,index_type         *nz           // [in/out] Number of nonzeros added (!=NULL)            
    ,value_type         *val_itr      // [out] Values to fill (!=NULL)
    ,index_type         *ivect_itr    // [out] Row indexes (can be NULL if load_struct == false)
    ,index_type         *jvect_itr    // [out] Column indexes  (can be NULL if load_struct == false)
    ) const;

};  // end class NLPSerialPreprocessExplJac

// ///////////////////////////
// inline members

inline
const NLPSerialPreprocessExplJac::FirstOrderExplInfo
NLPSerialPreprocessExplJac::first_order_expl_info() const
{
  return FirstOrderExplInfo(
    &Gc_nz_orig_
    ,&Gc_val_orig_
    ,&Gc_ivect_orig_
    ,&Gc_jvect_orig_
    ,&Gh_nz_orig_
    ,&Gh_val_orig_
    ,&Gh_ivect_orig_
    ,&Gh_jvect_orig_
    ,obj_grad_orig_info()
    );
}

} // end namespace NLPInterfacePack 

#endif // NLP_SERIAL_PREPROCESS_EXPL_JAC_H