Branch data Line data Source code
1 : : /* *****************************************************************
2 : : MESQUITE -- The Mesh Quality Improvement Toolkit
3 : :
4 : : Copyright 2006 Sandia National Laboratories. Developed at the
5 : : University of Wisconsin--Madison under SNL contract number
6 : : 624796. The U.S. Government and the University of Wisconsin
7 : : retain certain rights to this software.
8 : :
9 : : This library is free software; you can redistribute it and/or
10 : : modify it under the terms of the GNU Lesser General Public
11 : : License as published by the Free Software Foundation; either
12 : : version 2.1 of the License, or (at your option) any later version.
13 : :
14 : : This library is distributed in the hope that it will be useful,
15 : : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 : : Lesser General Public License for more details.
18 : :
19 : : You should have received a copy of the GNU Lesser General Public License
20 : : (lgpl.txt) along with this library; if not, write to the Free Software
21 : : Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22 : :
23 : : (2006) [email protected]
24 : :
25 : : ***************************************************************** */
26 : :
27 : : /** \file TShapeSize3DNB1.cpp
28 : : * \brief
29 : : * \author Jason Kraftcheck
30 : : */
31 : :
32 : : #include "Mesquite.hpp"
33 : : #include "TShapeSize3DNB1.hpp"
34 : : #include "MsqMatrix.hpp"
35 : : #include "TMPDerivs.hpp"
36 : :
37 : : namespace MBMesquite
38 : : {
39 : :
40 : 6 : std::string TShapeSize3DNB1::get_name() const
41 : : {
42 [ + - ]: 6 : return "TShapeSize3DNB1";
43 : : }
44 : :
45 [ - + ]: 46 : TShapeSize3DNB1::~TShapeSize3DNB1() {}
46 : :
47 : 0 : bool TShapeSize3DNB1::evaluate( const MsqMatrix< 3, 3 >& T, double& result, MsqError& /*err*/ )
48 : : {
49 : 0 : const double nT = Frobenius( T );
50 : 0 : const double tau = det( T );
51 : 0 : const double tau1 = tau - 1;
52 : 0 : result = nT * nT * nT - 3 * MSQ_SQRT_THREE * tau + mGamma * tau1 * tau1;
53 : 0 : return true;
54 : : }
55 : :
56 : 0 : bool TShapeSize3DNB1::evaluate_with_grad( const MsqMatrix< 3, 3 >& T, double& result, MsqMatrix< 3, 3 >& wrt_T,
57 : : MsqError& /*err*/ )
58 : : {
59 : 0 : const double nT = Frobenius( T );
60 : 0 : const double tau = det( T );
61 : 0 : const double tau1 = tau - 1;
62 : 0 : result = nT * nT * nT - 3 * MSQ_SQRT_THREE * tau + mGamma * tau1 * tau1;
63 : :
64 : 0 : wrt_T = T;
65 : 0 : wrt_T *= 3 * nT;
66 [ # # ][ # # ]: 0 : wrt_T -= ( 3 * MSQ_SQRT_THREE - 2 * mGamma * tau1 ) * transpose_adj( T );
67 : :
68 : 0 : return true;
69 : : }
70 : :
71 : 0 : bool TShapeSize3DNB1::evaluate_with_hess( const MsqMatrix< 3, 3 >& T, double& result, MsqMatrix< 3, 3 >& wrt_T,
72 : : MsqMatrix< 3, 3 > second[6], MsqError& /*err*/ )
73 : : {
74 [ # # ]: 0 : const double nT = Frobenius( T );
75 [ # # ]: 0 : const double tau = det( T );
76 : 0 : const double tau1 = tau - 1;
77 : 0 : result = nT * nT * nT - 3 * MSQ_SQRT_THREE * tau + mGamma * tau1 * tau1;
78 : :
79 : 0 : const double f = ( 3 * MSQ_SQRT_THREE - 2 * mGamma * tau1 );
80 [ # # ]: 0 : const MsqMatrix< 3, 3 > adjt = transpose_adj( T );
81 : 0 : wrt_T = T;
82 [ # # ]: 0 : wrt_T *= 3 * nT;
83 [ # # ][ # # ]: 0 : wrt_T -= f * adjt;
84 : :
85 [ # # ]: 0 : set_scaled_outer_product( second, 2 * mGamma, adjt );
86 [ # # ]: 0 : pluseq_scaled_2nd_deriv_of_det( second, -f, T );
87 [ # # ]: 0 : pluseq_scaled_I( second, 3 * nT );
88 : : // Could perturb T a bit if the norm is zero, but that would just
89 : : // result in the coefficent of the outer product being practically
90 : : // zero, so just skip the outer product in that case.
91 : : // Anyway nT approaches zero as T does, so the limit of this term
92 : : // as nT approaches zero is zero.
93 [ # # ]: 0 : if( nT > 1e-100 ) // NOTE: nT is always positive
94 [ # # ]: 0 : pluseq_scaled_outer_product( second, 3 / nT, T );
95 : :
96 : 0 : return true;
97 : : }
98 : :
99 : : } // namespace MBMesquite
|