MOAB: Mesh Oriented datABase  (version 5.4.1)
TShapeSizeOrientB1.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2006 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     (2006) [email protected]
00024 
00025   ***************************************************************** */
00026 
00027 /** \file TRel2DShapeSizeOrientBarrier.cpp
00028  *  \brief
00029  *  \author Jason Kraftcheck
00030  */
00031 
00032 #include "Mesquite.hpp"
00033 #include "TShapeSizeOrientB1.hpp"
00034 #include "MsqMatrix.hpp"
00035 #include "MsqError.hpp"
00036 #include "TMPDerivs.hpp"
00037 #include "TMPCommon.hpp"
00038 
00039 namespace MBMesquite
00040 {
00041 
00042 std::string TShapeSizeOrientB1::get_name() const
00043 {
00044     return "TShapeSizeOrientB1";
00045 }
00046 
00047 TShapeSizeOrientB1::~TShapeSizeOrientB1() {}
00048 
00049 bool TShapeSizeOrientB1::evaluate( const MsqMatrix< 2, 2 >& T, double& result, MsqError& err )
00050 {
00051     double tau = det( T );
00052     if( TMetric::invalid_determinant( tau ) )
00053     {
00054         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00055         return false;
00056     }
00057 
00058     MsqMatrix< 2, 2 > T_I( T );
00059     pluseq_scaled_I( T_I, -1 );
00060     result = sqr_Frobenius( T_I ) / ( 2 * tau );
00061     return true;
00062 }
00063 
00064 bool TShapeSizeOrientB1::evaluate_with_grad( const MsqMatrix< 2, 2 >& T,
00065                                              double& result,
00066                                              MsqMatrix< 2, 2 >& deriv_wrt_T,
00067                                              MsqError& err )
00068 {
00069     const double d = det( T );
00070     if( TMetric::invalid_determinant( d ) )
00071     {  // barrier
00072         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00073         return false;
00074     }
00075 
00076     deriv_wrt_T = T;
00077     pluseq_scaled_I( deriv_wrt_T, -1 );
00078     double inv_d = 1.0 / d;
00079     result       = 0.5 * sqr_Frobenius( deriv_wrt_T ) * inv_d;
00080 
00081     deriv_wrt_T -= result * transpose_adj( T );
00082     deriv_wrt_T *= inv_d;
00083 
00084     return true;
00085 }
00086 
00087 bool TShapeSizeOrientB1::evaluate_with_hess( const MsqMatrix< 2, 2 >& T,
00088                                              double& result,
00089                                              MsqMatrix< 2, 2 >& deriv_wrt_T,
00090                                              MsqMatrix< 2, 2 > second_wrt_T[3],
00091                                              MsqError& err )
00092 {
00093     const double d = det( T );
00094     if( TMetric::invalid_determinant( d ) )
00095     {  // barrier
00096         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00097         return false;
00098     }
00099 
00100     deriv_wrt_T = T;
00101     pluseq_scaled_I( deriv_wrt_T, -1.0 );
00102     double inv_d = 1.0 / d;
00103     result       = 0.5 * sqr_Frobenius( deriv_wrt_T ) * inv_d;
00104 
00105     MsqMatrix< 2, 2 > adjt = transpose_adj( T );
00106     set_scaled_outer_product( second_wrt_T, 2 * result * inv_d * inv_d, adjt );
00107     pluseq_scaled_sum_outer_product( second_wrt_T, -inv_d * inv_d, deriv_wrt_T, adjt );
00108     pluseq_scaled_2nd_deriv_of_det( second_wrt_T, -result * inv_d, T );
00109     pluseq_scaled_I( second_wrt_T, inv_d );
00110 
00111     deriv_wrt_T -= result * adjt;
00112     deriv_wrt_T *= inv_d;
00113 
00114     return true;
00115 }
00116 
00117 bool TShapeSizeOrientB1::evaluate( const MsqMatrix< 3, 3 >& T, double& result, MsqError& err )
00118 {
00119     double tau = det( T );
00120     if( TMetric::invalid_determinant( tau ) )
00121     {
00122         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00123         return false;
00124     }
00125 
00126     MsqMatrix< 3, 3 > T_I( T );
00127     pluseq_scaled_I( T_I, -1 );
00128     result = sqr_Frobenius( T_I ) / ( 2 * tau );
00129     return true;
00130 }
00131 
00132 bool TShapeSizeOrientB1::evaluate_with_grad( const MsqMatrix< 3, 3 >& T,
00133                                              double& result,
00134                                              MsqMatrix< 3, 3 >& deriv_wrt_T,
00135                                              MsqError& err )
00136 {
00137     const double d = det( T );
00138     if( TMetric::invalid_determinant( d ) )
00139     {  // barrier
00140         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00141         return false;
00142     }
00143 
00144     deriv_wrt_T = T;
00145     pluseq_scaled_I( deriv_wrt_T, -1 );
00146     double inv_d = 1.0 / d;
00147     result       = 0.5 * sqr_Frobenius( deriv_wrt_T ) * inv_d;
00148 
00149     deriv_wrt_T -= result * transpose_adj( T );
00150     deriv_wrt_T *= inv_d;
00151 
00152     return true;
00153 }
00154 
00155 bool TShapeSizeOrientB1::evaluate_with_hess( const MsqMatrix< 3, 3 >& T,
00156                                              double& result,
00157                                              MsqMatrix< 3, 3 >& deriv_wrt_T,
00158                                              MsqMatrix< 3, 3 > second_wrt_T[6],
00159                                              MsqError& err )
00160 {
00161     const double d = det( T );
00162     if( TMetric::invalid_determinant( d ) )
00163     {  // barrier
00164         MSQ_SETERR( err )( barrier_violated_msg, MsqError::BARRIER_VIOLATED );
00165         return false;
00166     }
00167 
00168     deriv_wrt_T = T;
00169     pluseq_scaled_I( deriv_wrt_T, -1.0 );
00170     double inv_d = 1.0 / d;
00171     result       = 0.5 * sqr_Frobenius( deriv_wrt_T ) * inv_d;
00172 
00173     MsqMatrix< 3, 3 > adjt = transpose_adj( T );
00174     set_scaled_outer_product( second_wrt_T, 2 * result * inv_d * inv_d, adjt );
00175     pluseq_scaled_sum_outer_product( second_wrt_T, -inv_d * inv_d, deriv_wrt_T, adjt );
00176     pluseq_scaled_2nd_deriv_of_det( second_wrt_T, -result * inv_d, T );
00177     pluseq_scaled_I( second_wrt_T, inv_d );
00178 
00179     deriv_wrt_T -= result * adjt;
00180     deriv_wrt_T *= inv_d;
00181 
00182     return true;
00183 }
00184 
00185 // TMP_T_TEMPL_IMPL_COMMON(TShapeSizeOrientB1)
00186 
00187 }  // namespace MBMesquite
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines