Branch data Line data Source code
1 : : /* *****************************************************************
2 : : MESQUITE -- The Mesh Quality Improvement Toolkit
3 : :
4 : : Copyright 2004 Sandia Corporation and Argonne National
5 : : Laboratory. Under the terms of Contract DE-AC04-94AL85000
6 : : with Sandia Corporation, the U.S. Government retains certain
7 : : rights in 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 : : [email protected], [email protected], [email protected],
24 : : [email protected], [email protected], [email protected]
25 : :
26 : : ***************************************************************** */
27 : : /*!
28 : : \file MeshTransform.cpp
29 : : \brief
30 : :
31 : : The MeshTransform Class is the base class for all the smoothing algorythms
32 : :
33 : : \author Michael Brewer
34 : : \date 2004-11-06
35 : : */
36 : :
37 : : #include "MeshTransform.hpp"
38 : : #include "MeshInterface.hpp"
39 : : #include "MsqVertex.hpp"
40 : : #include "MsqError.hpp"
41 : :
42 : : namespace MBMesquite
43 : : {
44 : :
45 [ - + ]: 2 : MeshTransform::~MeshTransform() {}
46 : :
47 : : /*!
48 : : Actually apply the affine transformation
49 : : */
50 : 1 : double MeshTransform::loop_over_mesh( MeshDomainAssoc* mesh_and_domain, const Settings*, MsqError& err )
51 : : {
52 [ + - ]: 1 : Mesh* mesh = mesh_and_domain->get_mesh();
53 : :
54 [ + - ]: 1 : std::vector< Mesh::VertexHandle > handle_list;
55 [ + - ]: 1 : mesh->get_all_vertices( handle_list, err );
56 [ + - ][ - + ]: 1 : if( MSQ_CHKERR( err ) ) return 1.0;
[ # # ][ # # ]
[ - + ]
57 : :
58 [ + - ]: 2 : std::vector< bool > fixed( 1 );
59 [ + - ]: 1 : MsqVertex vertex;
60 : 1 : std::vector< Mesh::VertexHandle >::const_iterator iter;
61 [ + - ][ + - ]: 56 : for( iter = handle_list.begin(); iter != handle_list.end(); ++iter )
[ + - ][ + + ]
62 : : {
63 [ + - ][ + - ]: 55 : mesh->vertices_get_coordinates( &*iter, &vertex, 1, err );
64 [ + - ][ - + ]: 55 : if( MSQ_CHKERR( err ) ) return 1.0;
[ # # ][ # # ]
[ - + ]
65 : :
66 [ - + ]: 55 : if( skipFixed )
67 : : {
68 [ # # ][ # # ]: 0 : mesh->vertices_get_fixed_flag( &*iter, fixed, 1, err );
69 [ # # ][ # # ]: 0 : if( MSQ_CHKERR( err ) ) return 1.0;
[ # # ][ # # ]
[ # # ]
70 [ # # ][ # # ]: 0 : if( fixed.front() ) continue;
71 : : }
72 : :
73 [ + - ][ + - ]: 55 : vertex = mMat * vertex + mVec;
[ + - ]
74 : :
75 [ + - ][ + - ]: 55 : mesh->vertex_set_coordinates( *iter, vertex, err );
76 [ + - ][ - + ]: 55 : if( MSQ_CHKERR( err ) ) return 1.0;
[ # # ][ # # ]
[ - + ]
77 : : }
78 : :
79 : 2 : return 0.0;
80 : : }
81 : :
82 : 0 : void MeshTransform::add_translation( const Vector3D& offset )
83 : : {
84 : 0 : mVec += offset;
85 : 0 : }
86 : :
87 : 0 : void MeshTransform::add_rotation( const Vector3D& axis, double radians )
88 : : {
89 : 0 : const double c = cos( radians );
90 : 0 : const double s = sin( radians );
91 [ # # ][ # # ]: 0 : const Vector3D a = axis / axis.length();
92 [ # # ][ # # ]: 0 : const Matrix3D m1( c, -a[2] * s, a[1] * s, a[2] * s, c, -a[0] * s, -a[1] * s, a[0] * s, c );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
93 [ # # ]: 0 : Matrix3D m2;
94 [ # # ]: 0 : m2.outer_product( a, a );
95 [ # # ][ # # ]: 0 : Matrix3D rot = m1 + ( 1.0 - c ) * m2;
96 [ # # ][ # # ]: 0 : mMat = rot * mMat;
97 [ # # ][ # # ]: 0 : mVec = rot * mVec;
98 : 0 : }
99 : :
100 : 0 : void MeshTransform::add_scale( double factor )
101 : : {
102 [ # # ]: 0 : add_scale( Vector3D( factor ) );
103 : 0 : }
104 : :
105 : 0 : void MeshTransform::add_scale( const Vector3D& f )
106 : : {
107 [ # # ]: 0 : for( int i = 0; i < 3; ++i )
108 : : {
109 : 0 : mVec[i] *= f[i];
110 : 0 : mMat[i][0] *= f[i];
111 : 0 : mMat[i][1] *= f[i];
112 : 0 : mMat[i][2] *= f[i];
113 : : }
114 : 0 : }
115 : :
116 : 0 : void MeshTransform::initialize_queue( MeshDomainAssoc*, const Settings*, MsqError& ) {}
117 : :
118 [ + - ][ + - ]: 4 : } // namespace MBMesquite
|