/usr/include/trilinos/QualityMetric.hpp is in libtrilinos-dev 10.4.0.dfsg-1ubuntu2.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 | /* *****************************************************************
MESQUITE -- The Mesh Quality Improvement Toolkit
Copyright 2004 Sandia Corporation and Argonne National
Laboratory. Under the terms of Contract DE-AC04-94AL85000
with Sandia Corporation, the U.S. Government retains certain
rights in this software.
This library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public License
(lgpl.txt) along with this library; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov,
pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov
(2006) kraftche@cae.wisc.edu
***************************************************************** */
/*! \file QualityMetric.hpp
\brief
Header file for the Mesquite::QualityMetric class
\author Thomas Leurent
\author Michael Brewer
\date 2002-05-01
*/
#ifndef QualityMetric_hpp
#define QualityMetric_hpp
#include <cmath>
#include <vector>
#include <algorithm>
#include "Mesquite.hpp"
#include "Vector3D.hpp"
#include "Matrix3D.hpp"
#ifdef _MSC_VER
typedef unsigned uint32_t;
#elif defined(MSQ_HAVE_STDINT_H)
# include <stdint.h>
#elif defined(MSQ_HAVE_INTTYPES_H)
# include <inttypes.h>
#endif
namespace MESQUITE_NS
{
/*! \class QualityMetric
\brief Base class for concrete quality metrics.
*/
class PatchData;
class MsqMeshEntity;
class Mesh;
class MeshDomain;
class Settings;
class QualityMetric
{
protected:
QualityMetric( ) {}
public:
enum MetricType
{
VERTEX_BASED, /**< Iterate over vertices to evaluate metric. */
ELEMENT_BASED /**< Iterate over elements to evaluate metric. */
};
MESQUITE_EXPORT virtual ~QualityMetric()
{}
MESQUITE_EXPORT virtual MetricType get_metric_type() const = 0;
MESQUITE_EXPORT virtual std::string get_name() const = 0;
//! 1 if metric should be minimized, -1 if metric should be maximized.
MESQUITE_EXPORT virtual int get_negate_flag() const = 0;
/**\brief Get locations at which metric can be evaluated
*
* Different metrics are evaluated for different things within
* a patch. For example, an element-based metric will be evaluated
* once for each element in patch, a vertex-based metric once for
* each veretx, and a target/sample-point based metric will be
* evaluated once for each samle point in each element. This method
* returns a list of handles, one for each location in the patch
* at which the metric can be evaluated. The handle values are used
* as input to the evaluate methods.
*\param pd The patch
*\param handles Output list of handles
*\param free_vertices_only If true, only pass back evaluation points
* that depend on at least one free vertex.
*/
MESQUITE_EXPORT virtual
void get_evaluations( PatchData& pd,
std::vector<size_t>& handles,
bool free_vertices_only,
MsqError& err ) = 0;
/**\brief Get locations at which metric can be evaluated for
* use in BCD intialization and QualityAssessor.
*
* For element-based, sample-based, and vertex-based metrics,
* this function is the same as get_evaluations. For edge-based
* metrics it returns only a subset of the results for get_evaluations
* such that each edge in the mesh is visited only once even though
* it would normally be visited twice when iterating over patches
* of the mesh. This assumes that no vertex occurs in more than one
* patch with its MSQ_PATCH_VERTEX set. This assumption is true for
* both element-on-vertex and global patches.
*\param pd The patch
*\param handles Output list of handles
*\param free_vertices_only If true, only pass back evaluation points
* that depend on at least one free vertex.
*/
MESQUITE_EXPORT virtual
void get_single_pass( PatchData& pd,
std::vector<size_t>& handles,
bool free_vertices_only,
MsqError& err );
/**\brief Get metric value at a logical location in the patch.
*
* Evaluate the metric at one location in the PatchData.
*\param pd The patch.
*\param handle The location in the patch (as passed back from get_evaluations).
*\param value The output metric value.
*/
MESQUITE_EXPORT virtual
bool evaluate( PatchData& pd,
size_t handle,
double& value,
MsqError& err ) = 0;
/**\brief Get metric value at a logical location in the patch.
*
* Evaluate the metric at one location in the PatchData.
*\param pd The patch.
*\param handle The location in the patch (as passed back from get_evaluations).
*\param value The output metric value.
*\param indices The free vertices that the evaluation is a function
* of, specified as vertex indices in the PatchData.
*/
MESQUITE_EXPORT virtual
bool evaluate_with_indices( PatchData& pd,
size_t handle,
double& value,
std::vector<size_t>& indices,
MsqError& err ) = 0;
/**\brief Get metric value and gradient at a logical location in the patch.
*
* Evaluate the metric at one location in the PatchData.
*\param pd The patch.
*\param handle The location in the patch (as passed back from get_evaluations).
*\param value The output metric value.
*\param indices The free vertices that the evaluation is a function
* of, specified as vertex indices in the PatchData.
*\param gradient The gradient of the metric as a function of the
* coordinates of the free vertices passed back in
* the indices list.
*/
MESQUITE_EXPORT virtual
bool evaluate_with_gradient( PatchData& pd,
size_t handle,
double& value,
std::vector<size_t>& indices,
std::vector<Vector3D>& gradient,
MsqError& err );
/**\brief Get metric value and gradient at a logical location in the patch.
*
* Evaluate the metric at one location in the PatchData.
*\param pd The patch.
*\param handle The location in the patch (as passed back from get_evaluations).
*\param value The output metric value.
*\param indices The free vertices that the evaluation is a function
* of, specified as vertex indices in the PatchData.
*\param gradient The gradient of the metric as a function of the
* coordinates of the free vertices passed back in
* the indices list.
*\param Hessian_diagonal The 3x3 blocks along the diagonal of
* the Hessian matrix.
*/
MESQUITE_EXPORT virtual
bool evaluate_with_Hessian_diagonal( PatchData& pd,
size_t handle,
double& value,
std::vector<size_t>& indices,
std::vector<Vector3D>& gradient,
std::vector<SymMatrix3D>& Hessian_diagonal,
MsqError& err );
/**\brief Get metric value and deravitives at a logical location in the patch.
*
* Evaluate the metric at one location in the PatchData.
*\param pd The patch.
*\param handle The location in the patch (as passed back from get_evaluations).
*\param value The output metric value.
*\param indices The free vertices that the evaluation is a function
* of, specified as vertex indices in the PatchData.
*\param gradient The gradient of the metric as a function of the
* coordinates of the free vertices passed back in
* the indices list.
*\param Hessian The Hessian of the metric as a function of the
* coordinates. The Hessian is passed back as the
* upper-triangular portion of the matrix in row-major
* order, where each Matrix3D is the portion of the
* Hessian with respect to the vertices at the
* corresponding positions in the indices list.
*/
MESQUITE_EXPORT virtual
bool evaluate_with_Hessian( PatchData& pd,
size_t handle,
double& value,
std::vector<size_t>& indices,
std::vector<Vector3D>& gradient,
std::vector<Matrix3D>& Hessian,
MsqError& err );
//!Escobar Barrier Function for Shape and Other Metrics
// det = signed determinant of Jacobian Matrix at a Vertex
// delta = scaling parameter
static inline double vertex_barrier_function(double det, double delta)
{ return 0.5*(det+sqrt(det*det+4*delta*delta)); }
//protected:
/** \brief Remove from vector any gradient terms corresponding
* to a fixed vertex.
*
* Remove terms from vector that correspond to fixed vertices.
*\param type Element type
*\param fixed_vertices Bit flags, one per vertex, 1 if
* vertex is fixed.
*\param gradients Array of gradients
*/
MESQUITE_EXPORT
static void remove_fixed_gradients( EntityTopology type,
uint32_t fixed_vertices,
std::vector<Vector3D>& gradients );
/** \brief Remove from vectors any gradient terms and hessian
* diagonal blcoks corresponding to a fixed vertex.
*
* Remove terms from vector that correspond to fixed vertices.
*\param type Element type
*\param fixed_vertices Bit flags, one per vertex, 1 if
* vertex is fixed.
*\param gradients Array of gradients
*\param hess_diagonal_blocks Array of diagonal blocks of Hessian matrix.
*/
MESQUITE_EXPORT
static void remove_fixed_diagonals( EntityTopology type,
uint32_t fixed_vertices,
std::vector<Vector3D>& gradients,
std::vector<SymMatrix3D>& hess_diagonal_blocks );
/** \brief Remove from vector any Hessian blocks corresponding
* to a fixed vertex.
*
* Remove blocks from vector that correspond to fixed vertices.
*\param type Element type
*\param fixed_vertices Bit flags, one per vertex, 1 if
* vertex is fixed.
*\param hessians Array of Hessian blocks (upper trianguler, row-major)
*/
MESQUITE_EXPORT
static void remove_fixed_hessians ( EntityTopology type,
uint32_t fixed_vertices,
std::vector<Matrix3D>& hessians );
/** \brief Convert fixed vertex format from list to bit flags
*
* Given list of pointers to fixed vertices as passed to
* evaluation functions, convert to bit flag format used
* for many utility functions in this class. Bits correspond
* to vertices in the canonical vertex order, beginning with
* the least-significant bit. The bit is cleared for free
* vertices and set (1) for fixed vertices.
*/
MESQUITE_EXPORT
static uint32_t fixed_vertex_bitmap( PatchData& pd,
const MsqMeshEntity* elem,
std::vector<size_t>& free_indices );
//! takes an array of coefficients and an array of metrics (both of length num_value)
//! and averages the contents using averaging method 'method'.
MESQUITE_EXPORT
double weighted_average_metrics(const double coef[],
const double metric_values[],
const int& num_values, MsqError &err);
/*!AveragingMethod allows you to set how the quality metric values
attained at each sample point will be averaged together to produce
a single metric value for an element.
*/
enum AveragingMethod
{
LINEAR, //!< the linear average
RMS, //!< the root-mean-squared average
HMS, //!< the harmonic-mean-squared average
SUM, //!< the sum of the values
SUM_SQUARED, //!< the sum of the squares of the values
HARMONIC, //!< the harmonic average
LAST_WITH_HESSIAN=HARMONIC,
MINIMUM, //!< the minimum value
MAXIMUM, //!< the maximum value
GEOMETRIC, //!< the geometric average
LAST_WITH_GRADIENT=GEOMETRIC,
STANDARD_DEVIATION, //!< the standard deviation squared of the values
MAX_OVER_MIN, //!< the maximum value minus the minum value
MAX_MINUS_MIN, //!< the maximum value divided by the minimum value
SUM_OF_RATIOS_SQUARED //!< (1/(N^2))*(SUM (SUM (v_i/v_j)^2))
};
//!\brief Called at start of instruction queue processing
//!
//! Do any preliminary global initialization, consistency checking,
//! etc. Default implementation does nothing.
virtual void initialize_queue( Mesh* mesh,
MeshDomain* domain,
const Settings* settings,
MsqError& err );
private:
int feasible;
std::vector<Matrix3D> tmpHess;
};
} //namespace
#endif // QualityMetric_hpp
|