MoochoPack : Framework for Large-Scale Optimization Algorithms  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
MoochoPack_ReducedHessianSerialization_Step.cpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
5 // Copyright (2003) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #include <fstream>
43 
44 #include "MoochoPack_ReducedHessianSerialization_Step.hpp"
45 #include "MoochoPack_moocho_algo_conversion.hpp"
46 #include "MoochoPack_Exceptions.hpp"
47 #include "IterationPack_print_algorithm_step.hpp"
48 #include "AbstractLinAlgPack_VectorSpace.hpp"
49 #include "AbstractLinAlgPack_VectorStdOps.hpp"
50 #include "AbstractLinAlgPack_VectorMutable.hpp"
51 #include "AbstractLinAlgPack_VectorOut.hpp"
52 #include "AbstractLinAlgPack_MatrixSymOpNonsing.hpp"
53 #include "AbstractLinAlgPack_MatrixSymInitDiag.hpp"
54 #include "AbstractLinAlgPack_MatrixOpOut.hpp"
55 #include "AbstractLinAlgPack_LinAlgOpPack.hpp"
56 #include "SerializationPack_Serializable.hpp"
57 #include "Teuchos_dyn_cast.hpp"
58 
59 namespace MoochoPack {
60 
62  const std::string &reduced_hessian_input_file_name
63  ,const std::string &reduced_hessian_output_file_name
64  )
65  :reduced_hessian_input_file_name_(reduced_hessian_input_file_name)
66  ,reduced_hessian_output_file_name_(reduced_hessian_output_file_name)
67 {}
68 
69 // Overridden from AlgorithmStep
70 
72  Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
73  ,poss_type assoc_step_poss
74  )
75 {
76  using Teuchos::dyn_cast;
78 
79  NLPAlgo &algo = rsqp_algo(_algo);
80  NLPAlgoState &s = algo.rsqp_state();
81 
82  EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
83  EJournalOutputLevel ns_olevel = algo.algo_cntr().null_space_journal_output_level();
84  std::ostream& out = algo.track().journal_out();
85 
86  // print step header.
87  if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) {
88  using IterationPack::print_algorithm_step;
89  print_algorithm_step( _algo, step_poss, type, assoc_step_poss, out );
90  }
91 
92  IterQuantityAccess<MatrixSymOp> &rHL_iq = s.rHL();
93 
94  if( !rHL_iq.updated_k(0) && reduced_hessian_input_file_name().length() ) {
95  int k_last_offset = rHL_iq.last_updated();
96  if( k_last_offset == IterQuantity::NONE_UPDATED ) {
97  if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
98  out
99  << "\nNo previous matrix rHL was found!\n"
100  << "\nReading in the matrix rHL_k from the file \""<<reduced_hessian_input_file_name()<<"\" ...\n";
101  }
102  MatrixSymOp &rHL_k = rHL_iq.set_k(0);
103  Serializable &rHL_serializable = dyn_cast<Serializable>(rHL_k);
104  std::ifstream reduced_hessian_input_file(reduced_hessian_input_file_name().c_str());
106  !reduced_hessian_input_file, std::logic_error
107  ,"ReducedHessianSerialization_Step::do_step(...): Error, the file \""<<reduced_hessian_input_file_name()<<"\""
108  " could not be opened or contains no input!"
109  );
110  rHL_serializable.unserialize(reduced_hessian_input_file);
111  if( (int)ns_olevel >= (int)PRINT_ALGORITHM_STEPS ) {
112  out << "\nrHL_k.rows() = " << rHL_k.rows() << std::endl;
113  out << "\nrHL_k.cols() = " << rHL_k.cols() << std::endl;
114  if(algo.algo_cntr().calc_matrix_norms())
115  out << "\n||rHL_k||inf = " << rHL_k.calc_norm(MatrixOp::MAT_NORM_INF).value << std::endl;
116  if(algo.algo_cntr().calc_conditioning()) {
117  const MatrixSymOpNonsing
118  *rHL_ns_k = dynamic_cast<const MatrixSymOpNonsing*>(&rHL_k);
119  if(rHL_ns_k)
120  out << "\ncond_inf(rHL_k) = " << rHL_ns_k->calc_cond_num(MatrixOp::MAT_NORM_INF).value << std::endl;
121  }
122  }
123  if( (int)ns_olevel >= (int)PRINT_ITERATION_QUANTITIES )
124  out << "\nrHL_k = \n" << rHL_k;
125  // Validate the space
126  const MatrixOp &Z_k = s.Z().get_k(0);
127  const VectorSpace &null_space = Z_k.space_rows();
129  !null_space.is_compatible(rHL_k.space_cols()) || !null_space.is_compatible(rHL_k.space_rows())
130  ,std::runtime_error
131  ,"ReducedHessianSerialization_Step::do_step(...): Error, the read-in reduced Hessian of dimension "
132  << rHL_k.rows() << " x " << rHL_k.cols() << " is not compatible with the null space of dimension "
133  "Z_k.cols() = " << Z_k.cols() << "!"
134  );
135  }
136  }
137 
138  return true;
139 
140 }
141 
143  Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
144  ,poss_type assoc_step_poss
145  )
146 {
147  using Teuchos::dyn_cast;
149 
150  const NLPAlgo &algo = rsqp_algo(_algo);
151  const NLPAlgoState &s = algo.rsqp_state();
152 
153  EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
154  std::ostream& out = algo.track().journal_out();
155 
156  const IterQuantityAccess<MatrixSymOp> &rHL_iq = s.rHL();
157 
158  int k_last_offset = rHL_iq.last_updated();
159 
160  if( k_last_offset != IterQuantity::NONE_UPDATED && reduced_hessian_output_file_name().length() ) {
161  if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
162  out
163  << "\nSerializing the matrix rHL_k("<<k_last_offset<<") to the file "
164  << "\""<<reduced_hessian_output_file_name()<<"\" ...\n";
165  }
166  const MatrixSymOp &rHL_k = rHL_iq.get_k(k_last_offset);
167  const Serializable &rHL_serializable = dyn_cast<const Serializable>(rHL_k);
168  std::ofstream reduced_hessian_output_file(reduced_hessian_output_file_name().c_str());
170  !reduced_hessian_output_file, std::logic_error
171  ,"ReducedHessianSerialization_Step::finalize_step(...): Error, the file \""<<reduced_hessian_output_file_name()<<"\""
172  " could not be opened!"
173  );
174  rHL_serializable.serialize(reduced_hessian_output_file);
175  }
176 }
177 
179  const Algorithm& algo, poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss
180  ,std::ostream& out, const std::string& L
181  ) const
182 {
183  out
184  << L << "*** Read in the reduced Hessian of the Lagrangian rHL from a file.\n"
185  << L << "if (rHL_k is not updated and reduced_hessian_input_file_name != \"\") then\n"
186  << L << " k_last_offset = last iteration rHL was updated for\n"
187  << L << " if k_last_offset==NONE_UPDATED then\n"
188  << L << " rHL_serializable = dyn_cast<Serializable>(rHL_k) *** Throws exception if fails!\n"
189  << L << " Unserialize into rHL_serializable from the file \""<<reduced_hessian_input_file_name()<<"\"\n"
190  << L << " else\n"
191  << L << " *** There is some reduced Hessian that exists for some past iteration so\n"
192  << L << " *** we will let some other step object initialize itQ\n"
193  << L << " end\n"
194  << L << "end\n"
195  << L << "*** Note: On finalization, this step object will serialize rHL_k to the file:\n"
196  << L << "*** \""<<reduced_hessian_output_file_name()<<"\""
197  ;
198 }
199 
200 } // end namespace MoochoPack
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
T_To & dyn_cast(T_From &from)
rSQP Algorithm control class.
void finalize_step(Algorithm &algo, poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss)
ReducedHessianSerialization_Step(const std::string &reduced_hessian_input_file_name="reduced_hessian.in", const std::string &reduced_hessian_output_file_name="reduced_hessian.out")
bool do_step(Algorithm &algo, poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss)
virtual std::ostream & journal_out() const
Reduced space SQP state encapsulation interface.
AlgorithmTracker & track()
void print_step(const Algorithm &algo, poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss, std::ostream &out, const std::string &leading_str) const
NLPAlgoState & rsqp_state()
<<std aggr>="">> members for algo_cntr