MOAB: Mesh Oriented datABase  (version 5.2.1)
AWMetricTest.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2008 Sandia National Laboratories.  Developed at the
00005     University of Wisconsin--Madison under SNL contract number
00006     624796.  The U.S. Government and the University of Wisconsin
00007     retain certain rights to this software.
00008 
00009     This library is free software; you can redistribute it and/or
00010     modify it under the terms of the GNU Lesser General Public
00011     License as published by the Free Software Foundation; either
00012     version 2.1 of the License, or (at your option) any later version.
00013 
00014     This library is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017     Lesser General Public License for more details.
00018 
00019     You should have received a copy of the GNU Lesser General Public License
00020     (lgpl.txt) along with this library; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 
00023     (2010) kraftche@cae.wisc.edu
00024 
00025   ***************************************************************** */
00026 
00027 /** \file AWMetricTest.cpp
00028  *  \brief Unit tests AWMetric base class functionality
00029  *  \author Jason Kraftcheck
00030  */
00031 
00032 #include "Mesquite.hpp"
00033 #include "AWMetric.hpp"
00034 #include "UnitUtil.hpp"
00035 #include "MsqError.hpp"
00036 
00037 using namespace MBMesquite;
00038 
00039 // Test functions implemented in class AWMetric
00040 class AWMetricTest : public CppUnit::TestFixture
00041 {
00042     CPPUNIT_TEST_SUITE( AWMetricTest );
00043     CPPUNIT_TEST( test_numerical_gradient_2D );
00044     CPPUNIT_TEST( test_numerical_hessian_2D );
00045     CPPUNIT_TEST( test_numerical_gradient_3D );
00046     CPPUNIT_TEST( test_numerical_hessian_3D );
00047     CPPUNIT_TEST_SUITE_END();
00048 
00049   public:
00050     void test_numerical_gradient_2D();
00051     void test_numerical_hessian_2D();
00052     void test_numerical_gradient_3D();
00053     void test_numerical_hessian_3D();
00054 };
00055 
00056 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( AWMetricTest, "Unit" );
00057 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( AWMetricTest, "AWMetricTest" );
00058 
00059 // implement metric such that the expected derivatives
00060 // at each location r,c in dm/dA are 3r+c+1
00061 class GradTestMetricAbs : public AWMetric
00062 {
00063   public:
00064     std::string get_name() const
00065     {
00066         return "GradTest";
00067     }
00068 
00069     static double grad( int r, int c )
00070     {
00071         return 3 * r + c + 1;
00072     }
00073 
00074     bool evaluate( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >& W, double& result, MsqError& )
00075     {
00076         result = 0;
00077         for( int r = 0; r < 2; ++r )
00078             for( int c = 0; c < 2; ++c )
00079                 result += grad( r, c ) * ( A( r, c ) - W( r, c ) );
00080         return true;
00081     }
00082 
00083     bool evaluate( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >& W, double& result, MsqError& )
00084     {
00085         result = 0;
00086         for( int r = 0; r < 3; ++r )
00087             for( int c = 0; c < 3; ++c )
00088                 result += grad( r, c ) * ( A( r, c ) - W( r, c ) );
00089         return true;
00090     }
00091 };
00092 
00093 // implement metric that is the |2A - A^T - W|^2
00094 // such that the Hessian is the constant:
00095 //  _          _
00096 // | 2  0  0  0 |
00097 // | 0 10 -8  0 |
00098 // | 0 -8 10  0 |
00099 // |_0  0  0  2_|
00100 //  _                         _
00101 // | 2  0  0  0  0  0  0  0  0 |
00102 // | 0 10  0 -8  0  0  0  0  0 |
00103 // | 0  0 10  0  0  0 -8  0  0 |
00104 // | 0 -8  0 10  0  0  0  0  0 |
00105 // | 0  0  0  0  2  0  0  0  0 |
00106 // | 0  0  0  0  0 10  0 -8  0 |
00107 // | 0  0 -8  0  0  0 10  0  0 |
00108 // | 0  0  0  0  0 -8  0 10  0 |
00109 // |_0  0  0  0  0  0  0  0  2_|
00110 class HessTestMetricAbs : public AWMetric
00111 {
00112   public:
00113     std::string get_name() const
00114     {
00115         return "HessTest";
00116     }
00117 
00118     bool evaluate( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >& W, double& result, MsqError& )
00119     {
00120         result = sqr_Frobenius( 2 * A - transpose( A ) - W );
00121         return true;
00122     }
00123     bool evaluate_with_grad( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >& W, double& result,
00124                              MsqMatrix< 2, 2 >& wrt_A, MsqError& )
00125     {
00126         result = sqr_Frobenius( 2 * A - transpose( A ) - W );
00127         wrt_A  = 10 * A - 8 * transpose( A ) - 4 * W + 2 * transpose( W );
00128         return true;
00129     }
00130 
00131     bool evaluate( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >& W, double& result, MsqError& )
00132     {
00133         result = sqr_Frobenius( 2 * A - transpose( A ) - W );
00134         return true;
00135     }
00136     bool evaluate_with_grad( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >& W, double& result,
00137                              MsqMatrix< 3, 3 >& wrt_A, MsqError& )
00138     {
00139         result = sqr_Frobenius( 2 * A - transpose( A ) - W );
00140         wrt_A  = 10 * A - 8 * transpose( A ) - 4 * W + 2 * transpose( W );
00141         return true;
00142     }
00143 };
00144 
00145 /** Simple target metric for testing second partial derivatives.
00146  *  \f$\mu(T) = |T|\f$
00147  *  \f$\frac{\partial\mu}{\partial T} = \frac{1}{|T|} T \f$
00148  *  \f$\frac{\partial^{2}\mu}{\partial t_{i,i}^2} = \frac{1}{|T|} - \frac{t_{i,i}^2}{|T|^3}\f$
00149  *  \f$\frac{\partial^{2}\mu}{\partial t_{i,j} \partial t_{k,l} (i \ne k or j \ne l)} =
00150  * -\frac{t_{i,j} a_{k,l}}{|T|^3}\f$
00151  */
00152 /** Simple target metric for testing second partial derivatives.
00153  *  \f$\mu(A,W) = |A|\f$
00154  *  \f$\frac{\partial\mu}{\partial A} = \frac{1}{|A|} A \f$
00155  *  \f$\frac{\partial^{2}\mu}{\partial a_{i,i}^2} = \frac{1}{|A|} - \frac{a_{i,i}^2}{|A|^3}\f$
00156  *  \f$\frac{\partial^{2}\mu}{\partial a_{i,j} \partial a_{k,l} (i \ne k or j \ne l)} =
00157  * -\frac{a_{i,j} a_{k,l}}{|A|^3}\f$
00158  */
00159 class HessTestMetricAbs_2 : public AWMetric
00160 {
00161   public:
00162     std::string get_name() const
00163     {
00164         return "HessTest2";
00165     }
00166 
00167     bool evaluate( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >&, double& result, MsqError& err )
00168     {
00169         result = Frobenius( A );
00170         return true;
00171     }
00172 
00173     bool evaluate_with_grad( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >&, double& result, MsqMatrix< 2, 2 >& d,
00174                              MsqError& err )
00175     {
00176         result = Frobenius( A );
00177         d      = A / result;
00178         return true;
00179     }
00180 
00181     bool evaluate_with_hess( const MsqMatrix< 2, 2 >& A, const MsqMatrix< 2, 2 >&, double& result, MsqMatrix< 2, 2 >& d,
00182                              MsqMatrix< 2, 2 > d2[3], MsqError& err )
00183     {
00184         result = Frobenius( A );
00185         d      = A / result;
00186         int h  = 0;
00187         for( int r = 0; r < 2; ++r )
00188         {
00189             int i = h;
00190             for( int c = r; c < 2; ++c )
00191                 d2[h++] = transpose( A.row( r ) ) * A.row( c ) / -( result * result * result );
00192             d2[i] += MsqMatrix< 2, 2 >( 1.0 / result );
00193         }
00194         return true;
00195     }
00196 
00197     bool evaluate( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >&, double& result, MsqError& err )
00198     {
00199         result = Frobenius( A );
00200         return true;
00201     }
00202 
00203     bool evaluate_with_grad( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >&, double& result, MsqMatrix< 3, 3 >& d,
00204                              MsqError& err )
00205     {
00206         result = Frobenius( A );
00207         d      = A / result;
00208         return true;
00209     }
00210 
00211     bool evaluate_with_hess( const MsqMatrix< 3, 3 >& A, const MsqMatrix< 3, 3 >&, double& result, MsqMatrix< 3, 3 >& d,
00212                              MsqMatrix< 3, 3 > d2[6], MsqError& err )
00213     {
00214         result = Frobenius( A );
00215         d      = A / result;
00216         int h  = 0;
00217         for( int r = 0; r < 3; ++r )
00218         {
00219             int i = h;
00220             for( int c = r; c < 3; ++c )
00221                 d2[h++] = transpose( A.row( r ) ) * A.row( c ) / -( result * result * result );
00222             d2[i] += MsqMatrix< 3, 3 >( 1.0 / result );
00223         }
00224         return true;
00225     }
00226 };
00227 
00228 void AWMetricTest::test_numerical_gradient_2D()
00229 {
00230     GradTestMetricAbs metric;
00231     HessTestMetricAbs_2 metric2;
00232     const double Avals[] = { 1, 2, 2, 5 };
00233     const double Bvals[] = { -0.1, -0.15, -0.25, -0.8 };
00234     const MsqMatrix< 2, 2 > I( 1.0 );
00235     const MsqMatrix< 2, 2 > A( Avals );
00236     const MsqMatrix< 2, 2 > B( Bvals );
00237 
00238     MsqError err;
00239     MsqMatrix< 2, 2 > d;
00240     bool valid;
00241     double val, gval;
00242 
00243     MsqMatrix< 2, 2 > expected;
00244     for( int r = 0; r < 2; ++r )
00245         for( int c = 0; c < 2; ++c )
00246             expected( r, c ) = metric.grad( r, c );
00247 
00248     valid = metric.evaluate( I, A, val, err );
00249     ASSERT_NO_ERROR( err );
00250     CPPUNIT_ASSERT( valid );
00251     valid = metric.evaluate_with_grad( I, A, gval, d, err );
00252     ASSERT_NO_ERROR( err );
00253     CPPUNIT_ASSERT( valid );
00254     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00255     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00256 
00257     valid = metric.evaluate( A, I, val, err );
00258     ASSERT_NO_ERROR( err );
00259     CPPUNIT_ASSERT( valid );
00260     valid = metric.evaluate_with_grad( A, I, gval, d, err );
00261     ASSERT_NO_ERROR( err );
00262     CPPUNIT_ASSERT( valid );
00263     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00264     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00265 
00266     valid = metric.evaluate( I, B, val, err );
00267     ASSERT_NO_ERROR( err );
00268     CPPUNIT_ASSERT( valid );
00269     valid = metric.evaluate_with_grad( I, B, gval, d, err );
00270     ASSERT_NO_ERROR( err );
00271     CPPUNIT_ASSERT( valid );
00272     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00273     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00274 
00275     valid = metric.evaluate( B, I, val, err );
00276     ASSERT_NO_ERROR( err );
00277     CPPUNIT_ASSERT( valid );
00278     valid = metric.evaluate_with_grad( B, I, gval, d, err );
00279     ASSERT_NO_ERROR( err );
00280     CPPUNIT_ASSERT( valid );
00281     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00282     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00283 
00284     valid = metric.evaluate( A, B, val, err );
00285     ASSERT_NO_ERROR( err );
00286     CPPUNIT_ASSERT( valid );
00287     valid = metric.evaluate_with_grad( A, B, gval, d, err );
00288     ASSERT_NO_ERROR( err );
00289     CPPUNIT_ASSERT( valid );
00290     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00291     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00292 
00293     valid = metric.evaluate( B, A, val, err );
00294     ASSERT_NO_ERROR( err );
00295     CPPUNIT_ASSERT( valid );
00296     valid = metric.evaluate_with_grad( B, A, gval, d, err );
00297     ASSERT_NO_ERROR( err );
00298     CPPUNIT_ASSERT( valid );
00299     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00300 
00301     MsqMatrix< 2, 2 > da;
00302     valid = metric2.evaluate_with_grad( A, I, val, da, err );
00303     ASSERT_NO_ERROR( err );
00304     CPPUNIT_ASSERT( valid );
00305     valid = metric2.AWMetric::evaluate_with_grad( A, I, gval, d, err );
00306     ASSERT_NO_ERROR( err );
00307     CPPUNIT_ASSERT( valid );
00308     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00309     ASSERT_MATRICES_EQUAL( da, d, 1e-6 );
00310 
00311     valid = metric2.evaluate_with_grad( B, I, val, da, err );
00312     ASSERT_NO_ERROR( err );
00313     CPPUNIT_ASSERT( valid );
00314     valid = metric2.AWMetric::evaluate_with_grad( B, I, gval, d, err );
00315     ASSERT_NO_ERROR( err );
00316     CPPUNIT_ASSERT( valid );
00317     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00318     ASSERT_MATRICES_EQUAL( da, d, 1e-6 );
00319 }
00320 
00321 void AWMetricTest::test_numerical_gradient_3D()
00322 {
00323     GradTestMetricAbs metric;
00324     HessTestMetricAbs_2 metric2;
00325     const double Avals[] = { 1, 2, 3, 4, 1, 4, 3, 2, 1 };
00326     const double Bvals[] = { 0.1, 0.15, 0.05, 0.2, -0.1, -0.15, -0.05, -0.2, 2 };
00327     const MsqMatrix< 3, 3 > I( 1.0 );
00328     const MsqMatrix< 3, 3 > A( Avals );
00329     const MsqMatrix< 3, 3 > B( Bvals );
00330 
00331     MsqError err;
00332     MsqMatrix< 3, 3 > d;
00333     bool valid;
00334     double val, gval;
00335 
00336     MsqMatrix< 3, 3 > expected;
00337     for( int r = 0; r < 3; ++r )
00338         for( int c = 0; c < 3; ++c )
00339             expected( r, c ) = metric.grad( r, c );
00340 
00341     valid = metric.evaluate( I, A, val, err );
00342     ASSERT_NO_ERROR( err );
00343     CPPUNIT_ASSERT( valid );
00344     valid = metric.evaluate_with_grad( I, A, gval, d, err );
00345     ASSERT_NO_ERROR( err );
00346     CPPUNIT_ASSERT( valid );
00347     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00348     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00349 
00350     valid = metric.evaluate( A, I, val, err );
00351     ASSERT_NO_ERROR( err );
00352     CPPUNIT_ASSERT( valid );
00353     valid = metric.evaluate_with_grad( A, I, gval, d, err );
00354     ASSERT_NO_ERROR( err );
00355     CPPUNIT_ASSERT( valid );
00356     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00357     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00358 
00359     valid = metric.evaluate( I, B, val, err );
00360     ASSERT_NO_ERROR( err );
00361     CPPUNIT_ASSERT( valid );
00362     valid = metric.evaluate_with_grad( I, B, gval, d, err );
00363     ASSERT_NO_ERROR( err );
00364     CPPUNIT_ASSERT( valid );
00365     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00366     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00367 
00368     valid = metric.evaluate( B, I, val, err );
00369     ASSERT_NO_ERROR( err );
00370     CPPUNIT_ASSERT( valid );
00371     valid = metric.evaluate_with_grad( B, I, gval, d, err );
00372     ASSERT_NO_ERROR( err );
00373     CPPUNIT_ASSERT( valid );
00374     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00375     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00376 
00377     valid = metric.evaluate( A, B, val, err );
00378     ASSERT_NO_ERROR( err );
00379     CPPUNIT_ASSERT( valid );
00380     valid = metric.evaluate_with_grad( A, B, gval, d, err );
00381     ASSERT_NO_ERROR( err );
00382     CPPUNIT_ASSERT( valid );
00383     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00384     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00385 
00386     valid = metric.evaluate( B, A, val, err );
00387     ASSERT_NO_ERROR( err );
00388     CPPUNIT_ASSERT( valid );
00389     valid = metric.evaluate_with_grad( B, A, gval, d, err );
00390     ASSERT_NO_ERROR( err );
00391     CPPUNIT_ASSERT( valid );
00392     ASSERT_MATRICES_EQUAL( expected, d, 1e-6 );
00393 
00394     MsqMatrix< 3, 3 > da;
00395     valid = metric2.evaluate_with_grad( A, I, val, da, err );
00396     ASSERT_NO_ERROR( err );
00397     CPPUNIT_ASSERT( valid );
00398     valid = metric2.AWMetric::evaluate_with_grad( A, I, gval, d, err );
00399     ASSERT_NO_ERROR( err );
00400     CPPUNIT_ASSERT( valid );
00401     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00402     ASSERT_MATRICES_EQUAL( da, d, 1e-6 );
00403 
00404     valid = metric2.evaluate_with_grad( B, I, val, da, err );
00405     ASSERT_NO_ERROR( err );
00406     CPPUNIT_ASSERT( valid );
00407     valid = metric2.AWMetric::evaluate_with_grad( B, I, gval, d, err );
00408     ASSERT_NO_ERROR( err );
00409     CPPUNIT_ASSERT( valid );
00410     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, gval, 1e-6 );
00411     ASSERT_MATRICES_EQUAL( da, d, 1e-6 );
00412 }
00413 
00414 void AWMetricTest::test_numerical_hessian_2D()
00415 {
00416     HessTestMetricAbs metric;
00417     HessTestMetricAbs_2 metric2;
00418     const double Avals[] = { 1, 2, 2, 5 };
00419     const double Bvals[] = { -0.1, -0.15, -0.25, -0.8 };
00420     const MsqMatrix< 2, 2 > I( 1.0 );
00421     const MsqMatrix< 2, 2 > A( Avals );
00422     const MsqMatrix< 2, 2 > B( Bvals );
00423 
00424     MsqError err;
00425     MsqMatrix< 2, 2 > g, gh;
00426     MsqMatrix< 2, 2 > h[6];
00427     bool valid;
00428     double val, hval;
00429 
00430     const double h_00[] = { 2, 0, 0, 10 };
00431     const double h_01[] = { 0, 0, -8, 0 };
00432     const double h_11[] = { 10, 0, 0, 2 };
00433     MsqMatrix< 2, 2 > h00( h_00 ), h01( h_01 ), h11( h_11 );
00434 
00435     valid = metric.evaluate_with_grad( I, A, val, g, err );
00436     ASSERT_NO_ERROR( err );
00437     CPPUNIT_ASSERT( valid );
00438     valid = metric.evaluate_with_hess( I, A, hval, gh, h, err );
00439     ASSERT_NO_ERROR( err );
00440     CPPUNIT_ASSERT( valid );
00441     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00442     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00443     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00444     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00445     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00446 
00447     valid = metric.evaluate_with_grad( A, I, val, g, err );
00448     ASSERT_NO_ERROR( err );
00449     CPPUNIT_ASSERT( valid );
00450     valid = metric.evaluate_with_hess( A, I, hval, gh, h, err );
00451     ASSERT_NO_ERROR( err );
00452     CPPUNIT_ASSERT( valid );
00453     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00454     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00455     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00456     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00457     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00458 
00459     valid = metric.evaluate_with_grad( I, B, val, g, err );
00460     ASSERT_NO_ERROR( err );
00461     CPPUNIT_ASSERT( valid );
00462     valid = metric.evaluate_with_hess( I, B, hval, gh, h, err );
00463     ASSERT_NO_ERROR( err );
00464     CPPUNIT_ASSERT( valid );
00465     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00466     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00467     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00468     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00469     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00470 
00471     valid = metric.evaluate_with_grad( B, I, val, g, err );
00472     ASSERT_NO_ERROR( err );
00473     CPPUNIT_ASSERT( valid );
00474     valid = metric.evaluate_with_hess( B, I, hval, gh, h, err );
00475     ASSERT_NO_ERROR( err );
00476     CPPUNIT_ASSERT( valid );
00477     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00478     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00479     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00480     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00481     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00482 
00483     valid = metric.evaluate_with_grad( A, B, val, g, err );
00484     ASSERT_NO_ERROR( err );
00485     CPPUNIT_ASSERT( valid );
00486     valid = metric.evaluate_with_hess( A, B, hval, gh, h, err );
00487     ASSERT_NO_ERROR( err );
00488     CPPUNIT_ASSERT( valid );
00489     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00490     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00491     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00492     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00493     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00494 
00495     valid = metric.evaluate_with_grad( B, A, val, g, err );
00496     ASSERT_NO_ERROR( err );
00497     CPPUNIT_ASSERT( valid );
00498     valid = metric.evaluate_with_hess( B, A, hval, gh, h, err );
00499     ASSERT_NO_ERROR( err );
00500     CPPUNIT_ASSERT( valid );
00501     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00502     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00503     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00504     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00505     ASSERT_MATRICES_EQUAL( h11, h[2], 1e-6 );
00506 
00507     MsqMatrix< 2, 2 > ah[6];
00508     valid = metric2.evaluate_with_hess( A, I, val, g, ah, err );
00509     ASSERT_NO_ERROR( err );
00510     CPPUNIT_ASSERT( valid );
00511     valid = metric2.AWMetric::evaluate_with_hess( A, I, hval, gh, h, err );
00512     ASSERT_NO_ERROR( err );
00513     CPPUNIT_ASSERT( valid );
00514     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00515     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00516     ASSERT_MATRICES_EQUAL( ah[0], h[0], 1e-6 );
00517     ASSERT_MATRICES_EQUAL( ah[1], h[1], 1e-6 );
00518     ASSERT_MATRICES_EQUAL( ah[2], h[2], 1e-6 );
00519 
00520     valid = metric2.evaluate_with_hess( B, I, val, g, ah, err );
00521     ASSERT_NO_ERROR( err );
00522     CPPUNIT_ASSERT( valid );
00523     valid = metric2.AWMetric::evaluate_with_hess( B, I, hval, gh, h, err );
00524     ASSERT_NO_ERROR( err );
00525     CPPUNIT_ASSERT( valid );
00526     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00527     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00528     ASSERT_MATRICES_EQUAL( ah[0], h[0], 1e-6 );
00529     ASSERT_MATRICES_EQUAL( ah[1], h[1], 1e-6 );
00530     ASSERT_MATRICES_EQUAL( ah[2], h[2], 1e-6 );
00531 }
00532 
00533 void AWMetricTest::test_numerical_hessian_3D()
00534 {
00535     HessTestMetricAbs metric;
00536     HessTestMetricAbs_2 metric2;
00537     const double Avals[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9 };
00538     const double Bvals[] = { 0.1, 0.15, 0.05, 0.2, -0.1, -0.15, -0.05, -0.2, 2 };
00539     const MsqMatrix< 3, 3 > I( 1.0 );
00540     const MsqMatrix< 3, 3 > A( Avals );
00541     const MsqMatrix< 3, 3 > B( Bvals );
00542 
00543     MsqError err;
00544     MsqMatrix< 3, 3 > g, gh;
00545     MsqMatrix< 3, 3 > h[6];
00546     bool valid;
00547     double val, hval;
00548 
00549     const double h_00[] = { 2, 0, 0, 0, 10, 0, 0, 0, 10 };
00550     const double h_01[] = { 0, 0, 0, -8, 0, 0, 0, 0, 0 };
00551     const double h_02[] = { 0, 0, 0, 0, 0, 0, -8, 0, 0 };
00552     const double h_11[] = { 10, 0, 0, 0, 2, 0, 0, 0, 10 };
00553     const double h_12[] = { 0, 0, 0, 0, 0, 0, 0, -8, 0 };
00554     const double h_22[] = { 10, 0, 0, 0, 10, 0, 0, 0, 2 };
00555     MsqMatrix< 3, 3 > h00( h_00 ), h01( h_01 ), h02( h_02 ), h11( h_11 ), h12( h_12 ), h22( h_22 );
00556 
00557     valid = metric.evaluate_with_grad( I, A, val, g, err );
00558     ASSERT_NO_ERROR( err );
00559     CPPUNIT_ASSERT( valid );
00560     valid = metric.evaluate_with_hess( I, A, hval, gh, h, err );
00561     ASSERT_NO_ERROR( err );
00562     CPPUNIT_ASSERT( valid );
00563     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00564     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00565     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00566     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00567     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00568     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00569     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00570     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00571 
00572     valid = metric.evaluate_with_grad( A, I, val, g, err );
00573     ASSERT_NO_ERROR( err );
00574     CPPUNIT_ASSERT( valid );
00575     valid = metric.evaluate_with_hess( A, I, hval, gh, h, err );
00576     ASSERT_NO_ERROR( err );
00577     CPPUNIT_ASSERT( valid );
00578     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00579     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00580     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00581     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00582     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00583     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00584     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00585     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00586 
00587     valid = metric.evaluate_with_grad( I, B, val, g, err );
00588     ASSERT_NO_ERROR( err );
00589     CPPUNIT_ASSERT( valid );
00590     valid = metric.evaluate_with_hess( I, B, hval, gh, h, err );
00591     ASSERT_NO_ERROR( err );
00592     CPPUNIT_ASSERT( valid );
00593     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00594     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00595     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00596     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00597     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00598     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00599     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00600     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00601 
00602     valid = metric.evaluate_with_grad( B, I, val, g, err );
00603     ASSERT_NO_ERROR( err );
00604     CPPUNIT_ASSERT( valid );
00605     valid = metric.evaluate_with_hess( B, I, hval, gh, h, err );
00606     ASSERT_NO_ERROR( err );
00607     CPPUNIT_ASSERT( valid );
00608     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00609     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00610     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00611     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00612     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00613     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00614     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00615     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00616 
00617     valid = metric.evaluate_with_grad( A, B, val, g, err );
00618     ASSERT_NO_ERROR( err );
00619     CPPUNIT_ASSERT( valid );
00620     valid = metric.evaluate_with_hess( A, B, hval, gh, h, err );
00621     ASSERT_NO_ERROR( err );
00622     CPPUNIT_ASSERT( valid );
00623     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00624     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00625     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00626     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00627     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00628     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00629     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00630     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00631 
00632     valid = metric.evaluate_with_grad( B, A, val, g, err );
00633     ASSERT_NO_ERROR( err );
00634     CPPUNIT_ASSERT( valid );
00635     valid = metric.evaluate_with_hess( B, A, hval, gh, h, err );
00636     ASSERT_NO_ERROR( err );
00637     CPPUNIT_ASSERT( valid );
00638     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00639     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00640     ASSERT_MATRICES_EQUAL( h00, h[0], 1e-6 );
00641     ASSERT_MATRICES_EQUAL( h01, h[1], 1e-6 );
00642     ASSERT_MATRICES_EQUAL( h02, h[2], 1e-6 );
00643     ASSERT_MATRICES_EQUAL( h11, h[3], 1e-6 );
00644     ASSERT_MATRICES_EQUAL( h12, h[4], 1e-6 );
00645     ASSERT_MATRICES_EQUAL( h22, h[5], 1e-6 );
00646 
00647     MsqMatrix< 3, 3 > ah[6];
00648     valid = metric2.evaluate_with_hess( A, I, val, g, ah, err );
00649     ASSERT_NO_ERROR( err );
00650     CPPUNIT_ASSERT( valid );
00651     valid = metric2.AWMetric::evaluate_with_hess( A, I, hval, gh, h, err );
00652     ASSERT_NO_ERROR( err );
00653     CPPUNIT_ASSERT( valid );
00654     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00655     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00656     ASSERT_MATRICES_EQUAL( ah[0], h[0], 1e-6 );
00657     ASSERT_MATRICES_EQUAL( ah[1], h[1], 1e-6 );
00658     ASSERT_MATRICES_EQUAL( ah[2], h[2], 1e-6 );
00659     ASSERT_MATRICES_EQUAL( ah[3], h[3], 1e-6 );
00660     ASSERT_MATRICES_EQUAL( ah[4], h[4], 1e-6 );
00661     ASSERT_MATRICES_EQUAL( ah[5], h[5], 1e-6 );
00662 
00663     valid = metric2.evaluate_with_hess( B, I, val, g, ah, err );
00664     ASSERT_NO_ERROR( err );
00665     CPPUNIT_ASSERT( valid );
00666     valid = metric2.AWMetric::evaluate_with_hess( B, I, hval, gh, h, err );
00667     ASSERT_NO_ERROR( err );
00668     CPPUNIT_ASSERT( valid );
00669     CPPUNIT_ASSERT_DOUBLES_EQUAL( val, hval, 1e-6 );
00670     ASSERT_MATRICES_EQUAL( g, gh, 1e-6 );
00671     ASSERT_MATRICES_EQUAL( ah[0], h[0], 1e-6 );
00672     ASSERT_MATRICES_EQUAL( ah[1], h[1], 1e-6 );
00673     ASSERT_MATRICES_EQUAL( ah[2], h[2], 1e-6 );
00674     ASSERT_MATRICES_EQUAL( ah[3], h[3], 1e-6 );
00675     ASSERT_MATRICES_EQUAL( ah[4], h[4], 1e-6 );
00676     ASSERT_MATRICES_EQUAL( ah[5], h[5], 1e-6 );
00677 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines