NLPInterfacePack: C++ Interfaces and Implementation for Non-Linear Programs  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
NLPInterfacePack_NLPDirectTester.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 <assert.h>
43 #include <math.h>
44 
45 #include <ostream>
46 #include <iomanip>
47 #include <sstream>
48 #include <limits>
49 
50 #include "NLPInterfacePack_NLPDirectTester.hpp"
51 #include "NLPInterfacePack_NLPDirect.hpp"
52 #include "AbstractLinAlgPack_VectorMutable.hpp"
53 #include "AbstractLinAlgPack_VectorOut.hpp"
54 #include "AbstractLinAlgPack_VectorSpace.hpp"
55 #include "AbstractLinAlgPack_VectorStdOps.hpp"
56 #include "AbstractLinAlgPack_LinAlgOpPack.hpp"
57 #include "AbstractLinAlgPack_LinAlgOpPack.hpp"
58 #include "AbstractLinAlgPack_assert_print_nan_inf.hpp"
59 #include "TestingHelperPack_update_success.hpp"
60 #include "Teuchos_Assert.hpp"
61 
62 namespace {
63 template< class T >
64 inline
65 T my_max( const T& v1, const T& v2 ) { return v1 > v2 ? v1 : v2; }
66 } // end namespace
67 
68 namespace NLPInterfacePack {
69 
71  const calc_fd_prod_ptr_t &calc_fd_prod
72  ,ETestingMethod Gf_testing_method
73  ,ETestingMethod Gc_testing_method
74  ,value_type Gf_warning_tol
75  ,value_type Gf_error_tol
76  ,value_type Gc_warning_tol
77  ,value_type Gc_error_tol
78  ,size_type num_fd_directions
79  ,bool dump_all
80  )
81  :calc_fd_prod_(calc_fd_prod)
82  ,Gf_testing_method_(Gf_testing_method)
83  ,Gc_testing_method_(Gc_testing_method)
84  ,Gf_warning_tol_(Gf_warning_tol)
85  ,Gf_error_tol_(Gf_error_tol)
86  ,Gc_warning_tol_(Gc_warning_tol)
87  ,Gc_error_tol_(Gc_error_tol)
88  ,num_fd_directions_(num_fd_directions)
89  ,dump_all_(dump_all)
90 {
91  if(calc_fd_prod_ == Teuchos::null)
92  calc_fd_prod_ = Teuchos::rcp(new CalcFiniteDiffProd());
93 }
94 
96  NLPDirect *nlp
97  ,const Vector &xo
98  ,const Vector *xl
99  ,const Vector *xu
100  ,const Vector *c
101  ,const Vector *Gf
102  ,const Vector *py
103  ,const Vector *rGf
104  ,const MatrixOp *GcU
105  ,const MatrixOp *D
106  ,const MatrixOp *Uz
107  ,bool print_all_warnings
108  ,std::ostream *out
109  ) const
110 {
111 
112  using std::setw;
113  using std::endl;
114  using std::right;
115 
120  using AbstractLinAlgPack::assert_print_nan_inf;
121  using LinAlgOpPack::V_StV;
122  using LinAlgOpPack::V_StMtV;
123  using LinAlgOpPack::Vp_MtV;
124  using LinAlgOpPack::M_StM;
125  using LinAlgOpPack::M_StMtM;
126 
127  typedef VectorSpace::vec_mut_ptr_t vec_mut_ptr_t;
128 
129 // using AbstractLinAlgPack::TestingPack::CompareDenseVectors;
130 // using AbstractLinAlgPack::TestingPack::CompareDenseSparseMatrices;
131 
132  using TestingHelperPack::update_success;
133 
134  bool success = true, preformed_fd;
135  if(out) {
136  *out << std::boolalpha
137  << std::endl
138  << "*********************************************************\n"
139  << "*** NLPDirectTester::finite_diff_check(...) ***\n"
140  << "*********************************************************\n";
141  }
142 
143  const Range1D
144  var_dep = nlp->var_dep(),
145  var_indep = nlp->var_indep(),
146  con_decomp = nlp->con_decomp(),
147  con_undecomp = nlp->con_undecomp();
149  space_x = nlp->space_x(),
150  space_c = nlp->space_c();
151 
152  // //////////////////////////////////////////////
153  // Validate the input
154 
156  py && !c, std::invalid_argument
157  ,"NLPDirectTester::finite_diff_check(...) : "
158  "Error, if py!=NULL then c!=NULL must also be true!" );
159 
160  const CalcFiniteDiffProd
161  &fd_deriv_prod = this->calc_fd_prod();
162 
163  const value_type
164  rand_y_l = -1.0, rand_y_u = 1.0,
165  small_num = ::sqrt(std::numeric_limits<value_type>::epsilon());
166 
167  try {
168 
169  // ///////////////////////////////////////////////
170  // (1) Check Gf
171 
172  if(Gf) {
173  switch( Gf_testing_method() ) {
174  case FD_COMPUTE_ALL: {
175  // Compute FDGf outright
176  TEUCHOS_TEST_FOR_EXCEPT(true); // ToDo: update above!
177  break;
178  }
179  case FD_DIRECTIONAL: {
180  // Compute FDGF'*y using random y's
181  if(out)
182  *out
183  << "\nComparing products Gf'*y with finite difference values FDGf'*y for "
184  << "random y's ...\n";
185  vec_mut_ptr_t y = space_x->create_member();
186  value_type max_warning_viol = 0.0;
187  int num_warning_viol = 0;
188  const int num_fd_directions_used = ( num_fd_directions() > 0 ? num_fd_directions() : 1 );
189  for( int direc_i = 1; direc_i <= num_fd_directions_used; ++direc_i ) {
190  if( num_fd_directions() > 0 ) {
191  random_vector( rand_y_l, rand_y_u, y.get() );
192  if(out)
193  *out
194  << "\n****"
195  << "\n**** Random directional vector " << direc_i << " ( ||y||_1 / n = "
196  << (y->norm_1() / y->dim()) << " )"
197  << "\n***\n";
198  }
199  else {
200  *y = 1.0;
201  if(out)
202  *out
203  << "\n****"
204  << "\n**** Ones vector y ( ||y||_1 / n = "<<(y->norm_1()/y->dim())<<" )"
205  << "\n***\n";
206  }
207  value_type
208  Gf_y = dot( *Gf, *y ),
209  FDGf_y;
210  preformed_fd = fd_deriv_prod.calc_deriv_product(
211  xo,xl,xu
212  ,*y,NULL,NULL,true,nlp,&FDGf_y,NULL,out,dump_all(),dump_all()
213  );
214  if( !preformed_fd )
215  goto FD_NOT_PREFORMED;
216  assert_print_nan_inf(FDGf_y, "FDGf'*y",true,out);
217  const value_type
218  calc_err = ::fabs( ( Gf_y - FDGf_y )/( ::fabs(Gf_y) + ::fabs(FDGf_y) + small_num ) );
219  if( calc_err >= Gf_warning_tol() ) {
220  max_warning_viol = my_max( max_warning_viol, calc_err );
221  ++num_warning_viol;
222  }
223  if(out)
224  *out
225  << "\nrel_err(Gf'*y,FDGf'*y) = "
226  << "rel_err(" << Gf_y << "," << FDGf_y << ") = "
227  << calc_err << endl;
228  if( calc_err >= Gf_error_tol() ) {
229  if(out) {
230  *out
231  << "Error, above relative error exceeded Gf_error_tol = " << Gf_error_tol() << endl;
232  if(dump_all()) {
233  *out << "\ny =\n" << *y;
234  }
235  }
236  }
237  }
238  if(out && num_warning_viol)
239  *out
240  << "\nThere were " << num_warning_viol << " warning tolerance "
241  << "violations out of num_fd_directions = " << num_fd_directions()
242  << " computations of FDGf'*y\n"
243  << "and the maximum violation was " << max_warning_viol
244  << " > Gf_waring_tol = " << Gf_warning_tol() << endl;
245  break;
246  }
247  default:
248  TEUCHOS_TEST_FOR_EXCEPT(true); // Invalid value
249  }
250  }
251 
252  // /////////////////////////////////////////////
253  // (2) Check py = -inv(C)*c
254  //
255  // We want to check;
256  //
257  // FDC * (inv(C)*c) \approx c
258  // \_________/
259  // -py
260  //
261  // We can compute this as:
262  //
263  // FDC * py = [ FDC, FDN ] * [ -py ; 0 ]
264  // \__________/
265  // FDA'
266  //
267  // t1 = [ -py ; 0 ]
268  //
269  // t2 = FDA'*t1
270  //
271  // Compare t2 \approx c
272  //
273  if(py) {
274  if(out)
275  *out
276  << "\nComparing c with finite difference product FDA'*[ -py; 0 ] = -FDC*py ...\n";
277  // t1 = [ -py ; 0 ]
278  VectorSpace::vec_mut_ptr_t
279  t1 = space_x->create_member();
280  V_StV( t1->sub_view(var_dep).get(), -1.0, *py );
281  *t1->sub_view(var_indep) = 0.0;
282  // t2 = FDA'*t1
283  VectorSpace::vec_mut_ptr_t
284  t2 = nlp->space_c()->create_member();
285  preformed_fd = fd_deriv_prod.calc_deriv_product(
286  xo,xl,xu
287  ,*t1,NULL,NULL,true,nlp,NULL,t2.get(),out,dump_all(),dump_all()
288  );
289  if( !preformed_fd )
290  goto FD_NOT_PREFORMED;
291  const value_type
292  sum_c = sum(*c),
293  sum_t2 = sum(*t2);
294  assert_print_nan_inf(sum_t2, "sum(-FDC*py)",true,out);
295  const value_type
296  calc_err = ::fabs( ( sum_c - sum_t2 )/( ::fabs(sum_c) + ::fabs(sum_t2) + small_num ) );
297  if(out)
298  *out
299  << "\nrel_err(sum(c),sum(-FDC*py)) = "
300  << "rel_err(" << sum_c << "," << sum_t2 << ") = "
301  << calc_err << endl;
302  if( calc_err >= Gc_error_tol() ) {
303  if(out)
304  *out
305  << "Error, above relative error exceeded Gc_error_tol = " << Gc_error_tol() << endl;
306  if(print_all_warnings)
307  *out << "\nt1 = [ -py; 0 ] =\n" << *t1
308  << "\nt2 = FDA'*t1 = -FDC*py =\n" << *t2;
309  update_success( false, &success );
310  }
311  if( calc_err >= Gc_warning_tol() ) {
312  if(out)
313  *out
314  << "\nWarning, above relative error exceeded Gc_warning_tol = " << Gc_warning_tol() << endl;
315  }
316  }
317 
318  // /////////////////////////////////////////////
319  // (3) Check D = -inv(C)*N
320 
321  if(D) {
322  switch( Gc_testing_method() ) {
323  case FD_COMPUTE_ALL: {
324  //
325  // Compute FDN outright and check
326  // -FDC * D \aprox FDN
327  //
328  // We want to compute:
329  //
330  // FDC * -D = [ FDC, FDN ] * [ -D; 0 ]
331  // \__________/
332  // FDA'
333  //
334  // To compute the above we perform:
335  //
336  // T = FDA' * [ -D; 0 ] (one column at a time)
337  //
338  // Compare T \approx FDN
339  //
340 /*
341  // FDN
342  DMatrix FDN(m,n-m);
343  fd_deriv_computer.calc_deriv( xo, xl, xu
344  , Range1D(m+1,n), nlp, NULL
345  , &FDN() ,BLAS_Cpp::trans, out );
346 
347  // T = FDA' * [ -D; 0 ] (one column at a time)
348  DMatrix T(m,n-m);
349  DVector t(n);
350  t(m+1,n) = 0.0;
351  for( int s = 1; s <= n-m; ++s ) {
352  // t = [ -D(:,s); 0 ]
353  V_StV( &t(1,m), -1.0, D->col(s) );
354  // T(:,s) = FDA' * t
355  fd_deriv_prod.calc_deriv_product(
356  xo,xl,xu,t(),NULL,NULL,nlp,NULL,&T.col(s),out);
357  }
358 
359  // Compare T \approx FDN
360  if(out)
361  *out
362  << "\nChecking the computed D = -inv(C)*N\n"
363  << "where D(i,j) = (-FDC*D)(i,j), dM(i,j) = FDN(i,j) ...\n";
364  result = comp_M.comp(
365  T(), FDN, BLAS_Cpp::no_trans
366  , CompareDenseSparseMatrices::FULL_MATRIX
367  , CompareDenseSparseMatrices::REL_ERR_BY_COL
368  , Gc_warning_tol(), Gc_error_tol()
369  , print_all_warnings, out );
370  update_success( result, &success );
371  if(!result) return false;
372 */
373  TEUCHOS_TEST_FOR_EXCEPT(true); // Todo: Implement above!
374  break;
375  }
376  case FD_DIRECTIONAL: {
377  //
378  // Compute -FDC * D * v \aprox FDN * v
379  // for random v's
380  //
381  // We will compute this as:
382  //
383  // t1 = [ 0; y ] <: R^(n)
384  //
385  // t2 = FDA' * t1 ( FDN * y ) <: R^(m)
386  //
387  // t1 = [ -D * y ; 0 ] <: R^(n)
388  //
389  // t3 = FDA' * t1 ( -FDC * D * y ) <: R^(m)
390  //
391  // Compare t2 \approx t3
392  //
393  if(out)
394  *out
395  << "\nComparing finite difference products -FDC*D*y with FDN*y for "
396  "random vectors y ...\n";
397  VectorSpace::vec_mut_ptr_t
398  y = space_x->sub_space(var_indep)->create_member(),
399  t1 = space_x->create_member(),
400  t2 = space_c->create_member(),
401  t3 = space_c->create_member();
402  value_type max_warning_viol = 0.0;
403  int num_warning_viol = 0;
404  const int num_fd_directions_used = ( num_fd_directions() > 0 ? num_fd_directions() : 1 );
405  for( int direc_i = 1; direc_i <= num_fd_directions_used; ++direc_i ) {
406  if( num_fd_directions() > 0 ) {
407  random_vector( rand_y_l, rand_y_u, y.get() );
408  if(out)
409  *out
410  << "\n****"
411  << "\n**** Random directional vector " << direc_i << " ( ||y||_1 / n = "
412  << (y->norm_1() / y->dim()) << " )"
413  << "\n***\n";
414  }
415  else {
416  *y = 1.0;
417  if(out)
418  *out
419  << "\n****"
420  << "\n**** Ones vector y ( ||y||_1 / n = "<<(y->norm_1()/y->dim())<<" )"
421  << "\n***\n";
422  }
423  // t1 = [ 0; y ] <: R^(n)
424  *t1->sub_view(var_dep) = 0.0;
425  *t1->sub_view(var_indep) = *y;
426  // t2 = FDA' * t1 ( FDN * y ) <: R^(m)
427  preformed_fd = fd_deriv_prod.calc_deriv_product(
428  xo,xl,xu
429  ,*t1,NULL,NULL,true,nlp,NULL,t2.get(),out,dump_all(),dump_all()
430  );
431  if( !preformed_fd )
432  goto FD_NOT_PREFORMED;
433  // t1 = [ -D * y ; 0 ] <: R^(n)
434  V_StMtV( t1->sub_view(var_dep).get(), -1.0, *D, BLAS_Cpp::no_trans, *y );
435  *t1->sub_view(var_indep) = 0.0;
436  // t3 = FDA' * t1 ( -FDC * D * y ) <: R^(m)
437  preformed_fd = fd_deriv_prod.calc_deriv_product(
438  xo,xl,xu
439  ,*t1,NULL,NULL,true,nlp,NULL,t3.get(),out,dump_all(),dump_all()
440  );
441  // Compare t2 \approx t3
442  const value_type
443  sum_t2 = sum(*t2),
444  sum_t3 = sum(*t3);
445  const value_type
446  calc_err = ::fabs( ( sum_t2 - sum_t3 )/( ::fabs(sum_t2) + ::fabs(sum_t3) + small_num ) );
447  if(out)
448  *out
449  << "\nrel_err(sum(-FDC*D*y),sum(FDN*y)) = "
450  << "rel_err(" << sum_t3 << "," << sum_t2 << ") = "
451  << calc_err << endl;
452  if( calc_err >= Gc_warning_tol() ) {
453  max_warning_viol = my_max( max_warning_viol, calc_err );
454  ++num_warning_viol;
455  }
456  if( calc_err >= Gc_error_tol() ) {
457  if(out)
458  *out
459  << "Error, above relative error exceeded Gc_error_tol = " << Gc_error_tol() << endl
460  << "Stoping the tests!\n";
461  if(print_all_warnings)
462  *out << "\ny =\n" << *y
463  << "\nt1 = [ -D*y; 0 ] =\n" << *t1
464  << "\nt2 = FDA' * [ 0; y ] = FDN * y =\n" << *t2
465  << "\nt3 = FDA' * t1 = -FDC * D * y =\n" << *t3;
466  update_success( false, &success );
467  }
468  }
469  if(out && num_warning_viol)
470  *out
471  << "\nThere were " << num_warning_viol << " warning tolerance "
472  << "violations out of num_fd_directions = " << num_fd_directions()
473  << " computations of sum(FDC*D*y) and sum(FDN*y)\n"
474  << "and the maximum relative iolation was " << max_warning_viol
475  << " > Gc_waring_tol = " << Gc_warning_tol() << endl;
476  break;
477  }
478  default:
480  }
481  }
482 
483  // ///////////////////////////////////////////////
484  // (4) Check rGf = Gf(var_indep) + D'*Gf(var_dep)
485 
486  if(rGf) {
487  if( Gf && D ) {
488  if(out)
489  *out
490  << "\nComparing rGf_tmp = Gf(var_indep) - D'*Gf(var_dep) with rGf ...\n";
491  VectorSpace::vec_mut_ptr_t
492  rGf_tmp = space_x->sub_space(var_indep)->create_member();
493  *rGf_tmp = *Gf->sub_view(var_indep);
494  Vp_MtV( rGf_tmp.get(), *D, BLAS_Cpp::trans, *Gf->sub_view(var_dep) );
495  const value_type
496  sum_rGf_tmp = sum(*rGf_tmp),
497  sum_rGf = sum(*rGf);
498  const value_type
499  calc_err = ::fabs( ( sum_rGf_tmp - sum_rGf )/( ::fabs(sum_rGf_tmp) + ::fabs(sum_rGf) + small_num ) );
500  if(out)
501  *out
502  << "\nrel_err(sum(rGf_tmp),sum(rGf)) = "
503  << "rel_err(" << sum_rGf_tmp << "," << sum_rGf << ") = "
504  << calc_err << endl;
505  if( calc_err >= Gc_error_tol() ) {
506  if(out)
507  *out
508  << "Error, above relative error exceeded Gc_error_tol = " << Gc_error_tol() << endl;
509  if(print_all_warnings)
510  *out << "\nrGf_tmp =\n" << *rGf_tmp
511  << "\nrGf =\n" << *rGf;
512  update_success( false, &success );
513  }
514  if( calc_err >= Gc_warning_tol() ) {
515  if(out)
516  *out
517  << "\nWarning, above relative error exceeded Gc_warning_tol = "
518  << Gc_warning_tol() << endl;
519  }
520  }
521  else if( D ) {
522  if(out)
523  *out
524  << "\nComparing rGf'*y with the finite difference product"
525  << " fd_prod(f,[D*y;y]) for random vectors y ...\n";
526  VectorSpace::vec_mut_ptr_t
527  y = space_x->sub_space(var_indep)->create_member(),
528  t = space_x->create_member();
529  value_type max_warning_viol = 0.0;
530  int num_warning_viol = 0;
531  const int num_fd_directions_used = ( num_fd_directions() > 0 ? num_fd_directions() : 1 );
532  for( int direc_i = 1; direc_i <= num_fd_directions_used; ++direc_i ) {
533  if( num_fd_directions() > 0 ) {
534  random_vector( rand_y_l, rand_y_u, y.get() );
535  if(out)
536  *out
537  << "\n****"
538  << "\n**** Random directional vector " << direc_i << " ( ||y||_1 / n = "
539  << (y->norm_1() / y->dim()) << " )"
540  << "\n***\n";
541  }
542  else {
543  *y = 1.0;
544  if(out)
545  *out
546  << "\n****"
547  << "\n**** Ones vector y ( ||y||_1 / n = "<<(y->norm_1()/y->dim())<<" )"
548  << "\n***\n";
549  }
550  // t = [ D*y; y ]
551  LinAlgOpPack::V_MtV(&*t->sub_view(var_dep),*D,BLAS_Cpp::no_trans,*y);
552  *t->sub_view(var_indep) = *y;
553  value_type fd_rGf_y = 0.0;
554  // fd_Gf_y
555  preformed_fd = fd_deriv_prod.calc_deriv_product(
556  xo,xl,xu
557  ,*t,NULL,NULL,true,nlp,&fd_rGf_y,NULL,out,dump_all(),dump_all()
558  );
559  if( !preformed_fd )
560  goto FD_NOT_PREFORMED;
561  if(out) *out << "fd_prod(f,[D*y;y]) = " << fd_rGf_y << "\n";
562  // rGf_y = rGf'*y
563  const value_type rGf_y = dot(*rGf,*y);
564  if(out) *out << "rGf'*y = " << rGf_y << "\n";
565  // Compare fd_rGf_y and rGf*y
566  const value_type
567  calc_err = ::fabs( ( rGf_y - fd_rGf_y )/( ::fabs(rGf_y) + ::fabs(fd_rGf_y) + small_num ) );
568  if( calc_err >= Gc_warning_tol() ) {
569  max_warning_viol = my_max( max_warning_viol, calc_err );
570  ++num_warning_viol;
571  }
572  if(out)
573  *out
574  << "\nrel_err(rGf'*y,fd_prod(f,[D*y;y])) = "
575  << "rel_err(" << fd_rGf_y << "," << rGf_y << ") = "
576  << calc_err << endl;
577  if( calc_err >= Gf_error_tol() ) {
578  if(out)
579  *out << "Error, above relative error exceeded Gc_error_tol = " << Gc_error_tol() << endl;
580  if(print_all_warnings)
581  *out << "\ny =\n" << *y
582  << "\nt = [ D*y; y ] =\n" << *t;
583  update_success( false, &success );
584  }
585  }
586  }
587  else {
588  TEUCHOS_TEST_FOR_EXCEPT(true); // ToDo: Test rGf without D? (This is not going to be easy!)
589  }
590  }
591 
592  // ///////////////////////////////////////////////////
593  // (5) Check GcU, and/or Uz (for undecomposed equalities)
594 
595  if(GcU || Uz) {
596  TEUCHOS_TEST_FOR_EXCEPT(true); // ToDo: Implement!
597  }
598 
599 FD_NOT_PREFORMED:
600 
601  if(!preformed_fd) {
602  if(out)
603  *out
604  << "\nError, the finite difference computation was not preformed due to cramped bounds\n"
605  << "Finite difference test failed!\n" << endl;
606  return false;
607  }
608 
609  } // end try
610  catch( const AbstractLinAlgPack::NaNInfException& except ) {
611  if(out)
612  *out
613  << "Error, found a NaN or Inf. Stoping tests\n";
614  success = false;
615  }
616 
617  if( out ) {
618  if( success )
619  *out
620  << "\nCongradulations, all the finite difference errors where within the\n"
621  "specified error tolerances!\n";
622  else
623  *out
624  << "\nOh no, at least one of the above finite difference tests failed!\n";
625  }
626 
627  return success;
628 
629 }
630 
631 } // end namespace NLPInterfacePack
NLPDirectTester(const calc_fd_prod_ptr_t &calc_fd_prod=Teuchos::null, ETestingMethod Gf_testing_method=FD_DIRECTIONAL, ETestingMethod Gc_testing_method=FD_DIRECTIONAL, value_type Gf_warning_tol=1e-6, value_type Gf_error_tol=1e-1, value_type Gc_warning_tol=1e-6, value_type Gc_error_tol=1e-1, size_type num_fd_directions=1, bool dump_all=false)
Constructor.
virtual Range1D con_decomp() const
Return the range of decomposed equality constraints.
virtual vec_space_ptr_t space_x() const =0
Vector space object for unknown variables x (dimension n).
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void Vp_StV(VectorMutable *v_lhs, const value_type &alpha, const Vector &v_rhs)
void V_StV(VectorMutable *v_lhs, value_type alpha, const V &V_rhs)
Interface providing only direct first order sensitivity information.
virtual Range1D con_undecomp() const
Return the range of undecomposed equality constraints.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
bool finite_diff_check(NLPDirect *nlp, const Vector &xo, const Vector *xl, const Vector *xu, const Vector *c, const Vector *Gf, const Vector *py, const Vector *rGf, const MatrixOp *GcU, const MatrixOp *D, const MatrixOp *Uz, bool print_all_warnings, std::ostream *out) const
This function takes an NLP object and its computed derivatives and function values and validates the ...
Strategy interface for computing the product of the derivatives of the functions of an NLP along give...
void M_StM(MatrixOp *M_lhs, value_type alpha, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
virtual vec_space_ptr_t space_c() const =0
Vector space object for general equality constraints c(x) (dimension m).
size_t size_type
virtual Range1D var_dep() const
Return the range of dependent (i.e. basic) variables.
virtual Range1D var_indep() const
Return the range of independent (i.e. nonbasic) variables.
value_type dot(const Vector &v_rhs1, const Vector &v_rhs2)
virtual bool calc_deriv_product(const Vector &xo, const Vector *xl, const Vector *xu, const Vector &v, const value_type *fo, const Vector *co, bool check_nan_inf, NLP *nlp, value_type *Gf_prod, VectorMutable *Gc_prod, std::ostream *out, bool trace=false, bool dump_all=false) const
Compute the directional derivatives by finite differences.
value_type sum(const Vector &v_rhs)
void random_vector(value_type l, value_type u, VectorMutable *v)
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)