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 retian 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) kraftche@cae.wisc.edu
24
25 ***************************************************************** */
26
27
28 /** \file IdealTargetTest.cpp
29 * \brief Test the IdealShapeTarget
30 * \author Jason Kraftcheck
31 */
32
33 #include "Mesquite.hpp"
34 #include "IdealShapeTarget.hpp"
35 #include "MsqError.hpp"
36 #include "Vector3D.hpp"
37 #include "PatchData.hpp"
38 #include "PlanarDomain.hpp"
39 #include "MappingFunction.hpp"
40 #include "Settings.hpp"
41 #include "IdealElements.hpp"
42 #include "ElemSampleQM.hpp"
43 #include "cppunit/extensions/HelperMacros.h"
44 #include "UnitUtil.hpp"
45
46 using namespace MBMesquite;
47
48 class IdealTargetTest : public CppUnit::TestFixture
49 {
50 private:
51 CPPUNIT_TEST_SUITE(IdealTargetTest);
52 CPPUNIT_TEST (test_tri_corner);
53 CPPUNIT_TEST (test_tri_edge);
54 CPPUNIT_TEST (test_tri_center);
55 CPPUNIT_TEST (test_hex_corner);
56 CPPUNIT_TEST (test_hex_edge);
57 CPPUNIT_TEST (test_hex_face);
58 CPPUNIT_TEST (test_hex_center);
59 CPPUNIT_TEST_SUITE_END();
60
61 public:
62
63 void test_tri_corner();
64 void test_tri_edge();
65 void test_tri_center();
66 void test_hex_corner();
67 void test_hex_edge();
68 void test_hex_face();
69 void test_hex_center();
70
71 private:
72
73 void get_calc_target( EntityTopology type,
74 Sample sample,
75 MsqMatrix<3,3>&,
76 MsqMatrix<2,2>& );
77
78 void get_ideal_target( EntityTopology type,
79 Sample sample,
80 MsqMatrix<3,3>&,
81 MsqMatrix<2,2>& );
82
83 void do_test( EntityTopology type, Sample location );
84
85 Settings settings;
86 };
87
88 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION(IdealTargetTest, "IdealTargetTest");
89 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION(IdealTargetTest, "Unit");
90
test_tri_corner()91 void IdealTargetTest::test_tri_corner()
92 {
93 do_test( TRIANGLE, Sample(0, 0) );
94 do_test( TRIANGLE, Sample(0, 1) );
95 do_test( TRIANGLE, Sample(0, 2) );
96 }
97
test_tri_edge()98 void IdealTargetTest::test_tri_edge()
99 {
100 do_test( TRIANGLE, Sample(1, 0) );
101 do_test( TRIANGLE, Sample(1, 1) );
102 do_test( TRIANGLE, Sample(1, 2) );
103 }
104
test_tri_center()105 void IdealTargetTest::test_tri_center()
106 {
107 do_test( TRIANGLE, Sample(2, 0) );
108 }
109
test_hex_corner()110 void IdealTargetTest::test_hex_corner()
111 {
112 do_test( HEXAHEDRON, Sample(0, 0) );
113 do_test( HEXAHEDRON, Sample(0, 1) );
114 do_test( HEXAHEDRON, Sample(0, 6) );
115 do_test( HEXAHEDRON, Sample(0, 7) );
116 }
117
test_hex_edge()118 void IdealTargetTest::test_hex_edge()
119 {
120 do_test( HEXAHEDRON, Sample(1, 1) );
121 do_test( HEXAHEDRON, Sample(1, 4) );
122 do_test( HEXAHEDRON, Sample(1, 11) );
123 }
124
test_hex_face()125 void IdealTargetTest::test_hex_face()
126 {
127 do_test( HEXAHEDRON, Sample(2, 0) );
128 do_test( HEXAHEDRON, Sample(2, 1) );
129 do_test( HEXAHEDRON, Sample(2, 2) );
130 do_test( HEXAHEDRON, Sample(2, 3) );
131 do_test( HEXAHEDRON, Sample(2, 4) );
132 do_test( HEXAHEDRON, Sample(2, 5) );
133 }
134
test_hex_center()135 void IdealTargetTest::test_hex_center()
136 {
137 do_test( HEXAHEDRON, Sample(3, 0) );
138 }
139
get_calc_target(EntityTopology type,Sample location,MsqMatrix<3,3> & w3,MsqMatrix<2,2> & w2)140 void IdealTargetTest::get_calc_target( EntityTopology type,
141 Sample location,
142 MsqMatrix<3,3>& w3,
143 MsqMatrix<2,2>& w2 )
144 {
145 MsqPrintError err( std::cout );
146 const int elem_dim = TopologyInfo::dimension(type);
147
148 // create a patch -- actual coords and such don't really matter
149 std::vector<double> coords( 24, 0.0 );
150 const size_t conn[] = { 0, 1, 2, 3, 4, 5, 6, 7 };
151 PatchData pd;
152 pd.fill( 8, arrptr(coords), 1, type, conn, 0, err );
153 CPPUNIT_ASSERT(!MSQ_CHKERR(err));
154 pd.attach_settings( &settings );
155
156 IdealShapeTarget tc;
157 if (elem_dim == 2)
158 tc.get_2D_target( pd, 0, location, w2, err );
159 else
160 tc.get_3D_target( pd, 0, location, w3, err );
161 CPPUNIT_ASSERT(!MSQ_CHKERR(err));
162 }
163
get_ideal_target(EntityTopology type,Sample location,MsqMatrix<3,3> & w3,MsqMatrix<2,2> & w2)164 void IdealTargetTest::get_ideal_target( EntityTopology type,
165 Sample location,
166 MsqMatrix<3,3>& w3,
167 MsqMatrix<2,2>& w2 )
168 {
169 MsqPrintError err( std::cout );
170 const unsigned elem_dim = TopologyInfo::dimension(type);
171
172 // get the target matrix for an ideal element
173 size_t indices[100];
174 size_t num_vtx;
175 const Vector3D* coords = unit_element( type );
176 Vector3D c[3];
177 if (elem_dim == 2) {
178 MsqVector<2> derivs[100];
179 const MappingFunction2D* func = settings.get_mapping_function_2D( type );
180 func->derivatives( location, NodeSet(), indices, derivs, num_vtx, err );
181 CPPUNIT_ASSERT(!MSQ_CHKERR(err));
182
183 MsqMatrix<3,2> J;
184 for (size_t i = 0; i < num_vtx; ++i)
185 for (unsigned j = 0; j < 2; ++j)
186 c[j] += derivs[i][j] * coords[indices[i]];
187
188 for (unsigned i = 0; i < 3; ++i)
189 for (unsigned j = 0; j < 2; ++j)
190 J(i,j) = c[j][i];
191
192 w2 = TargetCalculator::skew(J);
193 }
194 else {
195 MsqVector<3> derivs[100];
196 const MappingFunction3D* func = settings.get_mapping_function_3D( type );
197 func->derivatives( location, NodeSet(), indices, derivs, num_vtx, err );
198 CPPUNIT_ASSERT(!MSQ_CHKERR(err));
199
200 for (size_t i = 0; i < num_vtx; ++i)
201 for (unsigned j = 0; j < 3; ++j)
202 c[j] += derivs[i][j] * coords[indices[i]];
203
204 for (unsigned i = 0; i < 3; ++i)
205 for (unsigned j = 0; j < 3; ++j)
206 w3(i,j) = c[j][i];
207
208 w3 = TargetCalculator::skew(w3);
209 }
210
211 }
212
do_test(EntityTopology type,Sample location)213 void IdealTargetTest::do_test( EntityTopology type, Sample location )
214 {
215 MsqMatrix<3,3> w3_calc, w3_exp;
216 MsqMatrix<2,2> w2_calc, w2_exp;
217 get_calc_target( type, location, w3_calc, w2_calc );
218 get_ideal_target( type, location, w3_exp, w2_exp );
219 if (TopologyInfo::dimension(type) == 2) {
220 CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det(w2_calc), 1e-6 );
221 CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det(w2_exp), 1e-6 );
222 // a rotation of the expected matrix is acceptable.
223 MsqMatrix<2,2> R = inverse(w2_calc) * w2_exp;
224 ASSERT_MATRICES_EQUAL( transpose(R), inverse(R), 1e-6 );
225 }
226 else {
227 MsqMatrix<3,3> R = inverse(w3_calc) *w3_exp;
228 ASSERT_MATRICES_EQUAL( transpose(R), inverse(R), 1e-6 );
229 }
230 }
231
232
233