Unverified Commit fb582cf4 authored by corefmark's avatar corefmark Committed by GitHub
Browse files

Merge pull request #1 from spyridon97/merge-kitware-verdict

Merge Kitware Verdict
No related merge requests found
Showing with 6400 additions and 6829 deletions
+6400 -6829
# Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
# Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
# Under the terms of Contract DE-NA0003525 with NTESS,
# the U.S. Government retains certain rights in this software.
cmake_minimum_required( VERSION 3.16 )
cmake_minimum_required(VERSION 3.16)
project(verdict)
if(POLICY CMP0063)
project(verdict VERSION 1.4.0)
# Includes
include(GNUInstallDirs)
include(CMakePackageConfigHelpers)
if (POLICY CMP0063)
cmake_policy(SET CMP0063 NEW)
endif()
endif ()
set( verdict_MAJOR_VERSION "1")
set( verdict_MINOR_VERSION "4")
set( verdict_BUILD_VERSION "0")
set( verdict_VERSION_FLAT "${verdict_MAJOR_VERSION}${verdict_MINOR_VERSION}${verdict_BUILD_VERSION}" )
set( verdict_VERSION "${verdict_MAJOR_VERSION}.${verdict_MINOR_VERSION}.${verdict_BUILD_VERSION}" )
set(verdict_VERSION_FLAT "${CMAKE_PROJECT_VERSION_MAJOR}${CMAKE_PROJECT_VERSION_MINOR}${CMAKE_PROJECT_VERSION_PATCH}")
set(verdict_VERSION "${CMAKE_PROJECT_VERSION}")
option(VERDICT_BUILD_DOC "Build the 2007 Verdict User Manual" OFF)
option(VERDICT_MANGLE "Mangle verdict names for inclusion in a larger library?" OFF)
if ( VERDICT_MANGLE )
set( VERDICT_MANGLE_PREFIX "verdict" CACHE STRING "The namespace enclosing verdict function names and classes." )
mark_as_advanced( VERDICT_MANGLE_PREFIX )
if (VERDICT_MANGLE)
set(VERDICT_MANGLE_PREFIX "verdict" CACHE STRING "The namespace enclosing verdict function names and classes.")
mark_as_advanced(VERDICT_MANGLE_PREFIX)
endif ()
mark_as_advanced( VERDICT_MANGLE )
mark_as_advanced(VERDICT_MANGLE)
option( VERDICT_ENABLE_TESTING "Should tests of the VERDICT library be built?" ON )
option(VERDICT_ENABLE_TESTING "Should tests of the VERDICT library be built?" ON)
set(verdict_HEADERS
verdict.h
VerdictVector.hpp
verdict_defines.hpp)
set( verdict_SRCS
set(verdict_SOURCES
V_EdgeMetric.cpp
V_GaussIntegration.cpp
V_GaussIntegration.hpp
......@@ -35,101 +42,100 @@ set( verdict_SRCS
V_QuadMetric.cpp
V_TetMetric.cpp
V_TriMetric.cpp
v_vector.h
V_WedgeMetric.cpp
verdict.h
VerdictVector.cpp
VerdictVector.hpp
verdict_defines.hpp
)
VerdictVector.cpp)
configure_file(
${verdict_SOURCE_DIR}/verdict_config.h.in
${verdict_BINARY_DIR}/verdict_config.h
@ONLY
)
${CMAKE_CURRENT_SOURCE_DIR}/verdict_config.h.in
${CMAKE_CURRENT_BINARY_DIR}/verdict_config.h
@ONLY)
add_library( verdict ${verdict_SRCS} )
target_include_directories(verdict PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}> $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}>)
add_library(verdict ${verdict_SOURCES})
target_include_directories(verdict PUBLIC
$<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}> $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}>)
# Setting the VERSION and SOVERSION of a library will include
# version information either in the library, or in the library
# name (depending on the platform). You may choose to exclude
# this information.
if ( NOT VERDICT_NO_LIBRARY_VERSION )
set_target_properties( verdict PROPERTIES
if (NOT VERDICT_NO_LIBRARY_VERSION)
set_target_properties(verdict PROPERTIES
VERSION "${verdict_VERSION}"
SOVERSION "${verdict_MAJOR_VERSION}.${verdict_MINOR_VERSION}"
)
SOVERSION "${CMAKE_PROJECT_VERSION_MAJOR}.${CMAKE_PROJECT_VERSION_MINOR}")
endif ()
if ( NOT VERDICT_EXPORT_GROUP )
if (NOT VERDICT_EXPORT_GROUP)
set(VERDICT_EXPORT_GROUP VerdictExport)
endif ()
if ( VERDICT_ENABLE_TESTING )
if (VERDICT_ENABLE_TESTING)
enable_testing()
if(NOT TARGET GTest::GTest)
if (NOT TARGET GTest::GTest)
find_package(GTest REQUIRED)
endif()
function(ADD_VERDICT_UNITTESTS DIR)
ADD_SUBDIRECTORY(${verdict_SOURCE_DIR}/${DIR} ${verdict_BINARY_DIR}/${DIR})
FOREACH(TEST ${ARGN})
ADD_TEST(NAME UT-${TEST}
WORKING_DIRECTORY ${verdict_BINARY_DIR}/${DIR}
COMMAND ${TEST})
ENDFOREACH(TEST ${ARGN})
endfunction()
endif ()
ADD_VERDICT_UNITTESTS(unittests unittests_verdict)
endif ()
function(add_verdict_unittests DIR)
add_subdirectory(${CMAKE_CURRENT_SOURCE_DIR}/${DIR} ${CMAKE_CURRENT_BINARY_DIR}/${DIR})
if ( NOT verdict_INSTALL_DOC_DIR )
set (verdict_INSTALL_DOC_DIR doc)
endif ()
if ( NOT verdict_INSTALL_INCLUDE_DIR)
set (verdict_INSTALL_INCLUDE_DIR include)
endif ()
if ( NOT verdict_INSTALL_INCLUDE_SUBDIR)
set (verdict_INSTALL_INCLUDE_SUBDIR verdict)
endif ()
if ( NOT verdict_INSTALL_BIN_DIR)
set (verdict_INSTALL_BIN_DIR bin)
endif ()
if ( NOT verdict_INSTALL_LIB_DIR)
set (verdict_INSTALL_LIB_DIR lib)
foreach (TEST ${ARGN})
add_test(NAME UT-${TEST}
WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/${DIR}
COMMAND ${TEST})
endforeach (TEST ${ARGN})
endfunction()
add_verdict_unittests(unittests unittests_verdict)
endif ()
if ( VERDICT_BUILD_DOC )
add_subdirectory( docs/VerdictUserManual2007 )
if (VERDICT_BUILD_DOC)
add_subdirectory(docs/VerdictUserManual2007)
endif ()
#
# Installation stuff
#
install(FILES
README
DESTINATION ${verdict_INSTALL_DOC_DIR}/verdict/${verdict_VERSION}/
COMPONENT VerdictDevelopment
)
install(
FILES ${verdict_BINARY_DIR}/verdict_config.h
DESTINATION ${verdict_INSTALL_INCLUDE_DIR}
COMPONENT VerdictDevelopment
)
#########################
# Installation commands #
#########################
install(TARGETS verdict EXPORT ${VERDICT_EXPORT_GROUP}
RUNTIME DESTINATION ${verdict_INSTALL_BIN_DIR} COMPONENT Runtime # .exe, .dll
LIBRARY DESTINATION ${verdict_INSTALL_LIB_DIR} COMPONENT Runtime # .so, .dll
ARCHIVE DESTINATION ${verdict_INSTALL_LIB_DIR} COMPONENT VerdictDevelopment # .a, .lib
)
# Install documentation
install(FILES
README.md
DESTINATION ${CMAKE_INSTALL_DOCDIR}/verdict/${verdict_VERSION}/ COMPONENT VerdictDevelopment)
export(EXPORT ${VERDICT_EXPORT_GROUP} FILE VerdictConfig.cmake)
# Install required header files
install(FILES
${CMAKE_CURRENT_BINARY_DIR}/verdict_config.h
verdict.h
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR} COMPONENT VerdictDevelopment)
install(EXPORT ${VERDICT_EXPORT_GROUP} FILE VerdictConfig.cmake
DESTINATION ${verdict_INSTALL_LIB_DIR} COMPONENT VerdictDevelopment
)
# Install library
install(TARGETS verdict EXPORT ${VERDICT_EXPORT_GROUP}
RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR} COMPONENT Runtime # .exe, .dll
LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} COMPONENT Runtime # .so, .dll
ARCHIVE DESTINATION ${CMAKE_INSTALL_LIBDIR} COMPONENT VerdictDevelopment #[[.a, .lib]])
# Export Targets file
export(EXPORT ${VERDICT_EXPORT_GROUP}
FILE "${CMAKE_CURRENT_BINARY_DIR}/cmake/VerdictTargets.cmake"
NAMESPACE Verdict::)
# Install Targets file
install(EXPORT ${VERDICT_EXPORT_GROUP}
FILE "VerdictTargets.cmake"
NAMESPACE Verdict::
DESTINATION ${CMAKE_INSTALL_LIBDIR}/cmake/verdict COMPONENT VerdictDevelopment)
# Create Config file
configure_package_config_file(${CMAKE_CURRENT_SOURCE_DIR}/VerdictConfig.cmake.in
"${CMAKE_CURRENT_BINARY_DIR}/VerdictConfig.cmake"
INSTALL_DESTINATION ${CMAKE_INSTALL_LIBDIR}/cmake/verdict)
# Generate the version file for the Config file
write_basic_package_version_file(
"${CMAKE_CURRENT_BINARY_DIR}/VerdictConfigVersion.cmake"
VERSION "${verdict_VERSION}"
COMPATIBILITY AnyNewerVersion)
# Install Config and ConfigVersion files
install(FILES
"${CMAKE_CURRENT_BINARY_DIR}/VerdictConfig.cmake"
"${CMAKE_CURRENT_BINARY_DIR}/VerdictConfigVersion.cmake"
DESTINATION ${CMAKE_INSTALL_LIBDIR}/cmake/verdict COMPONENT VerdictDevelopment)
VERDICT
Computes quality functions of 2 and 3-dimensional Finite Elements.
What's new for 1.4:
- Add inradius metrics for several element types
- Bug fixes
What's new for 1.3:
- Thread safety for all functions and simplified API
- Some existing metrics implemented for pyramid and wedge elements
- New inradius and timestep metrics
- New implementation for some higher order elements
# Verdict
## Compute quality functions of 2 and 3-dimensional regions.
## Changelog
What's new for 1.4:
- Add inradius metrics for several element types
- Bug fixes
What's new for 1.3:
- Thread safety for all functions and simplified API
- Some existing metrics implemented for pyramid and wedge elements
- New inradius and timestep metrics
- New implementation for some higher order elements
## Building Verdict
To build with CMake
1. Set up a build directory
2. Change to this directory
3. Type `ccmake /path/to/verdict` but replace
`/path/to/verdict` with the path to the directory
containing this read-me file.
4. Fill the required fields and press the 'c' key
NB: This process is iterative;
you may need to change values and reconfigure before continuing.
5. When the 'g' option becomes available, press the 'g' key to generate
the Makefile and exit CMake.
6. Build with `make`
7. Install with `make install`
## Using Verdict
If you do not already have verdict installed on your system,
see below for instructions to compile it.
Once you have verdict installed, you may use it in your
own project by adding
```cmake
find_package(Verdict)
```
to your `CMakeLists.txt` file.
Now you can link your executable to the verdict library.
Assume you have an executable named `example`;
verdict uses CMake's _namespace_ feature, so you should
link to the verdict library like so:
```cmake
target_link_libraries(example
PUBLIC
Verdict::verdict
)
```
Note that executables do not need to specify `PUBLIC`, `PRIVATE`,
or `INTERFACE` linkage, but libraries must so that transitive
dependencies can be resolved by cmake.
`verdict.h` is the only header file you need to include.
......@@ -3,36 +3,36 @@
Module: V_EdgeMetric.cpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* V_EdgeMetric.cpp contains quality calcultions for edges
* V_EdgeMetric.cpp contains quality calculations for edges
*
* This file is part of VERDICT
*
*/
#include "verdict.h"
#include <math.h>
#include <cmath>
namespace VERDICT_NAMESPACE
{
/*!\brief Length of and edge.
* Length is calculated by taking the distance between the end nodes.
*/
double edge_length( int /*num_nodes*/, double coordinates[][3] )
double edge_length(int /*num_nodes*/, const double coordinates[][3])
{
double x = coordinates[1][0] - coordinates[0][0];
double y = coordinates[1][1] - coordinates[0][1];
double z = coordinates[1][2] - coordinates[0][2];
return (double)( sqrt (x*x + y*y + z*z) );
return (double)(std::sqrt(x * x + y * y + z * z));
}
} // namespace verdict
This diff is collapsed.
......@@ -3,13 +3,13 @@
Module: V_GaussIntegration.hpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* GaussIntegration.hpp declaration of gauss integration functions
......@@ -18,110 +18,99 @@
*
*/
// .SECTION Thanks
// Prior to its inclusion within VTK, this code was developed by the CUBIT
// project at Sandia National Laboratories.
#ifndef GAUSS_INTEGRATION_HPP
#define GAUSS_INTEGRATION_HPP
#include "verdict.h"
namespace VERDICT_NAMESPACE
{
#define maxTotalNumberGaussPoints 27
#define maxNumberNodes 20
#define maxNumberGaussPoints 3
#define maxNumberGaussPointsTri 6
#define maxNumberGaussPointsTet 4
struct GaussIntegration
{
void get_signs_for_node_local_coord_hex(int node_id, double &sign_y1,
double &sign_y2, double &sign_y3);
//- to get the signs for coordinates of hex nodes in the local computational space
//constructors
void initialize(int n=2, int m=4, int dim=2, int tri=0);
//manipulators
void get_gauss_pts_and_weight();
//- get gauss point locations and weights
void get_tri_rule_pts_and_weight();
//- get integration points and weights for triangular rules
void calculate_shape_function_2d_tri();
//- calculate the shape functions and derivatives of shape functions
//- at integration points for 2D triangular elements
void calculate_shape_function_2d_quad();
//- calculate the shape functions and derivatives of shape functions
//- at gaussian points for 2D quad elements
void get_shape_func(double shape_function[], double dndy1_at_gauss_pts[], double dndy2_at_gauss_ptsp[], double gauss_weight[]);
//- get shape functions and the derivatives
void get_shape_func(double shape_function[], double dndy1_at_gauss_pts[],
double dndy2_at_gauss_pts[], double dndy3_at_gauss_pts[],
double gauss_weight[]);
//- get shape functions and the derivatives for 3D elements
void calculate_derivative_at_nodes(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes
void calculate_shape_function_3d_hex();
//- calculate shape functions and derivatives of shape functions
//- at gaussian points for 3D hex elements
void calculate_derivative_at_nodes_3d(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes],
double dndy3_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for hex elements
void calculate_derivative_at_nodes_2d_tri(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for triangular elements
void calculate_shape_function_3d_tet();
//- calculate shape functions and derivatives of shape functions
//- at integration points for 3D tet elements
void get_tet_rule_pts_and_weight();
//- get integration points and weights for tetrhedron rules
void calculate_derivative_at_nodes_3d_tet(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes],
double dndy3_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for tetrahedron elements
void get_node_local_coord_tet(int node_id, double &y1, double &y2,
double &y3, double &y4);
//- get nodal volume coordinates for tetrahedron element
int numberGaussPoints;
int numberNodes;
int numberDims;
double gaussPointY[maxNumberGaussPoints];
double gaussWeight[maxNumberGaussPoints];
double shapeFunction[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy1GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy2GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy3GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double totalGaussWeight[maxTotalNumberGaussPoints];
int totalNumberGaussPts;
double y1Area[maxNumberGaussPointsTri];
double y2Area[maxNumberGaussPointsTri];
double y1Volume[maxNumberGaussPointsTet];
double y2Volume[maxNumberGaussPointsTet];
double y3Volume[maxNumberGaussPointsTet];
double y4Volume[maxNumberGaussPointsTet];
void get_signs_for_node_local_coord_hex(
int node_id, double& sign_y1, double& sign_y2, double& sign_y3);
//- to get the signs for coordinates of hex nodes in the local computational space
// constructors
void initialize(int n = 2, int m = 4, int dim = 2, int tri = 0);
// manipulators
void get_gauss_pts_and_weight();
//- get gauss point locations and weights
void get_tri_rule_pts_and_weight();
//- get integration points and weights for triangular rules
void calculate_shape_function_2d_tri();
//- calculate the shape functions and derivatives of shape functions
//- at integration points for 2D triangular elements
void calculate_shape_function_2d_quad();
//- calculate the shape functions and derivatives of shape functions
//- at gaussian points for 2D quad elements
void get_shape_func(double shape_function[], double dndy1_at_gauss_pts[],
double dndy2_at_gauss_ptsp[], double gauss_weight[]);
//- get shape functions and the derivatives
void get_shape_func(double shape_function[], double dndy1_at_gauss_pts[],
double dndy2_at_gauss_pts[], double dndy3_at_gauss_pts[], double gauss_weight[]);
//- get shape functions and the derivatives for 3D elements
void calculate_derivative_at_nodes(
double dndy1_at_nodes[][maxNumberNodes], double dndy2_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes
void calculate_shape_function_3d_hex();
//- calculate shape functions and derivatives of shape functions
//- at gaussian points for 3D hex elements
void calculate_derivative_at_nodes_3d(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes], double dndy3_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for hex elements
void calculate_derivative_at_nodes_2d_tri(
double dndy1_at_nodes[][maxNumberNodes], double dndy2_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for triangular elements
void calculate_shape_function_3d_tet();
//- calculate shape functions and derivatives of shape functions
//- at integration points for 3D tet elements
void get_tet_rule_pts_and_weight();
//- get integration points and weights for tetrhedron rules
void calculate_derivative_at_nodes_3d_tet(double dndy1_at_nodes[][maxNumberNodes],
double dndy2_at_nodes[][maxNumberNodes], double dndy3_at_nodes[][maxNumberNodes]);
//- calculate shape function derivatives at nodes for tetrahedron elements
void get_node_local_coord_tet(int node_id, double& y1, double& y2, double& y3, double& y4);
//- get nodal volume coordinates for tetrahedron element
int numberGaussPoints;
int numberNodes;
int numberDims;
double gaussPointY[maxNumberGaussPoints];
double gaussWeight[maxNumberGaussPoints];
double shapeFunction[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy1GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy2GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double dndy3GaussPts[maxTotalNumberGaussPoints][maxNumberNodes];
double totalGaussWeight[maxTotalNumberGaussPoints];
int totalNumberGaussPts;
double y1Area[maxNumberGaussPointsTri];
double y2Area[maxNumberGaussPointsTri];
double y1Volume[maxNumberGaussPointsTet];
double y2Volume[maxNumberGaussPointsTet];
double y3Volume[maxNumberGaussPointsTet];
double y4Volume[maxNumberGaussPointsTet];
};
} // namespace verdict
#endif
#endif
This diff is collapsed.
/*=========================================================================
Module: V_HexMetric.hpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* V_HexMetric.hpp contains declarations of hex related shape quantities
......@@ -21,17 +21,19 @@
#ifndef VERDICT_HEX_METRIC_HPP
#define VERDICT_HEX_METRIC_HPP
#include <assert.h>
#include <math.h>
#include "verdict.h"
namespace verdict {
#include <cassert>
#include <cmath>
//
// Compute jacobian at each of the eight hex corner nodes
//
namespace VERDICT_NAMESPACE
{
//
// Compute jacobian at each of the eight hex corner nodes
//
template <typename T>
void hex_nodal_jacobians(const T coords[24], T Jdet8x[8]) {
void hex_nodal_jacobians(const T coords[24], T Jdet8x[8])
{
T x0 = coords[0];
T y0 = coords[1];
T z0 = coords[2];
......@@ -91,16 +93,23 @@ void hex_nodal_jacobians(const T coords[24], T Jdet8x[8]) {
T x6y7 = x6 * y7 - x7 * y6;
Jdet8x[0] = (-x1y3 + x1y4 - x3y4) * z0 + (x0y3 - x0y4 + x3y4) * z1 + (-x0y1 + x0y4 - x1y4) * z3 + (x0y1 - x0y3 + x1y3) * z4;
Jdet8x[1] = (-x1y2 + x1y5 - x2y5) * z0 + (x0y2 - x0y5 + x2y5) * z1 + (-x0y1 + x0y5 - x1y5) * z2 + (x0y1 - x0y2 + x1y2) * z5;
Jdet8x[2] = (-x2y3 + x2y6 - x3y6) * z1 + (x1y3 - x1y6 + x3y6) * z2 + (-x1y2 + x1y6 - x2y6) * z3 + (x1y2 - x1y3 + x2y3) * z6;
Jdet8x[3] = (-x2y3 + x2y7 - x3y7) * z0 + (x0y3 - x0y7 + x3y7) * z2 + (-x0y2 + x0y7 - x2y7) * z3 + (x0y2 - x0y3 + x2y3) * z7;
Jdet8x[4] = (-x4y5 + x4y7 - x5y7) * z0 + (x0y5 - x0y7 + x5y7) * z4 + (-x0y4 + x0y7 - x4y7) * z5 + (x0y4 - x0y5 + x4y5) * z7;
Jdet8x[5] = (-x4y5 + x4y6 - x5y6) * z1 + (x1y5 - x1y6 + x5y6) * z4 + (-x1y4 + x1y6 - x4y6) * z5 + (x1y4 - x1y5 + x4y5) * z6;
Jdet8x[6] = (-x5y6 + x5y7 - x6y7) * z2 + (x2y6 - x2y7 + x6y7) * z5 + (-x2y5 + x2y7 - x5y7) * z6 + (x2y5 - x2y6 + x5y6) * z7;
Jdet8x[7] = (-x4y6 + x4y7 - x6y7) * z3 + (x3y6 - x3y7 + x6y7) * z4 + (-x3y4 + x3y7 - x4y7) * z6 + (x3y4 - x3y6 + x4y6) * z7;
Jdet8x[0] = (-x1y3 + x1y4 - x3y4) * z0 + (x0y3 - x0y4 + x3y4) * z1 + (-x0y1 + x0y4 - x1y4) * z3 +
(x0y1 - x0y3 + x1y3) * z4;
Jdet8x[1] = (-x1y2 + x1y5 - x2y5) * z0 + (x0y2 - x0y5 + x2y5) * z1 + (-x0y1 + x0y5 - x1y5) * z2 +
(x0y1 - x0y2 + x1y2) * z5;
Jdet8x[2] = (-x2y3 + x2y6 - x3y6) * z1 + (x1y3 - x1y6 + x3y6) * z2 + (-x1y2 + x1y6 - x2y6) * z3 +
(x1y2 - x1y3 + x2y3) * z6;
Jdet8x[3] = (-x2y3 + x2y7 - x3y7) * z0 + (x0y3 - x0y7 + x3y7) * z2 + (-x0y2 + x0y7 - x2y7) * z3 +
(x0y2 - x0y3 + x2y3) * z7;
Jdet8x[4] = (-x4y5 + x4y7 - x5y7) * z0 + (x0y5 - x0y7 + x5y7) * z4 + (-x0y4 + x0y7 - x4y7) * z5 +
(x0y4 - x0y5 + x4y5) * z7;
Jdet8x[5] = (-x4y5 + x4y6 - x5y6) * z1 + (x1y5 - x1y6 + x5y6) * z4 + (-x1y4 + x1y6 - x4y6) * z5 +
(x1y4 - x1y5 + x4y5) * z6;
Jdet8x[6] = (-x5y6 + x5y7 - x6y7) * z2 + (x2y6 - x2y7 + x6y7) * z5 + (-x2y5 + x2y7 - x5y7) * z6 +
(x2y5 - x2y6 + x5y6) * z7;
Jdet8x[7] = (-x4y6 + x4y7 - x6y7) * z3 + (x3y6 - x3y7 + x6y7) * z4 + (-x3y4 + x3y7 - x4y7) * z6 +
(x3y4 - x3y6 + x4y6) * z7;
}
}
#endif
......@@ -3,13 +3,13 @@
Module: V_KnifeMetric.cpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* KnifeMetrics.cpp contains quality calculations for knives
......@@ -18,10 +18,8 @@
*
*/
#include "verdict.h"
#include "VerdictVector.hpp"
#include <memory.h>
#include "verdict.h"
namespace VERDICT_NAMESPACE
{
......@@ -37,116 +35,65 @@ namespace VERDICT_NAMESPACE
| 1\/ | |
| \ |
|_____\|_____|
4 5 6
4 5 6
(edge 3,5 is is a hidden line if you will)
if this is hard to visualize, consider a hex
with nodes 5 and 7 becoming the same node
*/
/*!
calculates the volume of a knife element
this is done by dividing the knife into 4 tets
and summing the volumes of each.
*/
double knife_volume( int num_nodes, double coordinates[][3] )
*/
double knife_volume(int num_nodes, const double coordinates[][3])
{
double volume = 0;
VerdictVector side1, side2, side3;
if (num_nodes == 7)
{
// divide the knife into 4 tets and calculate the volume
side1.set(
coordinates[1][0] - coordinates[0][0],
coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2]
);
side2.set(
coordinates[3][0] - coordinates[0][0],
coordinates[3][1] - coordinates[0][1],
coordinates[3][2] - coordinates[0][2]
);
side3.set(
coordinates[4][0] - coordinates[0][0],
coordinates[4][1] - coordinates[0][1],
coordinates[4][2] - coordinates[0][2]
);
side1.set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2]);
side2.set(coordinates[3][0] - coordinates[0][0], coordinates[3][1] - coordinates[0][1],
coordinates[3][2] - coordinates[0][2]);
side3.set(coordinates[4][0] - coordinates[0][0], coordinates[4][1] - coordinates[0][1],
coordinates[4][2] - coordinates[0][2]);
volume = side3 % (side1 * side2) / 6;
side1.set(
coordinates[5][0] - coordinates[1][0],
coordinates[5][1] - coordinates[1][1],
coordinates[5][2] - coordinates[1][2]
);
side2.set(
coordinates[3][0] - coordinates[1][0],
coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]
);
side3.set(
coordinates[4][0] - coordinates[1][0],
coordinates[4][1] - coordinates[1][1],
coordinates[4][2] - coordinates[1][2]
);
side1.set(coordinates[5][0] - coordinates[1][0], coordinates[5][1] - coordinates[1][1],
coordinates[5][2] - coordinates[1][2]);
side2.set(coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]);
side3.set(coordinates[4][0] - coordinates[1][0], coordinates[4][1] - coordinates[1][1],
coordinates[4][2] - coordinates[1][2]);
volume += side3 % (side1 * side2) / 6;
side1.set(
coordinates[2][0] - coordinates[1][0],
coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2]
);
side2.set(
coordinates[3][0] - coordinates[1][0],
coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]
);
side3.set(
coordinates[6][0] - coordinates[1][0],
coordinates[6][1] - coordinates[1][1],
coordinates[6][2] - coordinates[1][2]
);
side1.set(coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2]);
side2.set(coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]);
side3.set(coordinates[6][0] - coordinates[1][0], coordinates[6][1] - coordinates[1][1],
coordinates[6][2] - coordinates[1][2]);
volume += side3 % (side1 * side2) / 6;
side1.set(
coordinates[3][0] - coordinates[1][0],
coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]
);
side2.set(
coordinates[5][0] - coordinates[1][0],
coordinates[5][1] - coordinates[1][1],
coordinates[5][2] - coordinates[1][2]
);
side3.set(
coordinates[6][0] - coordinates[1][0],
coordinates[6][1] - coordinates[1][1],
coordinates[6][2] - coordinates[1][2]
);
side1.set(coordinates[3][0] - coordinates[1][0], coordinates[3][1] - coordinates[1][1],
coordinates[3][2] - coordinates[1][2]);
side2.set(coordinates[5][0] - coordinates[1][0], coordinates[5][1] - coordinates[1][1],
coordinates[5][2] - coordinates[1][2]);
side3.set(coordinates[6][0] - coordinates[1][0], coordinates[6][1] - coordinates[1][1],
coordinates[6][2] - coordinates[1][2]);
volume += side3 % (side1 * side2) / 6;
}
return (double)volume;
}
} // namespace verdict
......@@ -3,13 +3,13 @@
Module: V_PyramidMetric.cpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* PyramidMetrics.cpp contains quality calculations for Pyramids
......@@ -18,36 +18,35 @@
*
*/
#include "verdict.h"
#include "VerdictVector.hpp"
#include "verdict_defines.hpp"
#include <memory.h>
#include <vector>
#include <array>
#include <algorithm>
static const double SQRT2_HALVES = sqrt(2.0)/2.0;
#include "verdict.h"
extern double quad_equiangle_skew( int num_nodes, double coordinates[][3] );
extern double tri_equiangle_skew( int num_nodes, double coordinates[][3] );
#include <algorithm>
#include <array>
namespace VERDICT_NAMESPACE
{
extern double quad_equiangle_skew(int num_nodes, const double coordinates[][3]);
extern double tri_equiangle_skew(int num_nodes, const double coordinates[][3]);
static const double sqrt2_2 = std::sqrt(2.0) / 2.0;
// local methods
void make_pyramid_tets(double coordinates[][3], double tet1_coords[][3], double tet2_coords[][3],
double tet3_coords[][3], double tet4_coords[][3]);
void make_pyramid_faces(double coordinates[][3], double base[][3], double tri1[][3],
double tri2[][3], double tri3[][3], double tri4[][3]);
void make_pyramid_edges( VerdictVector edges[8], double coordinates[][3] );
double distance_point_to_pyramid_base( int num_nodes, double coordinates[][3], double &cos_angle);
double largest_pyramid_edge( double coordinates[][3] );
void make_pyramid_tets(const double coordinates[][3], double tet1_coords[][3],
double tet2_coords[][3], double tet3_coords[][3], double tet4_coords[][3]);
void make_pyramid_faces(const double coordinates[][3], double base[][3], double tri1[][3],
double tri2[][3], double tri3[][3], double tri4[][3]);
void make_pyramid_edges(VerdictVector edges[8], const double coordinates[][3]);
double distance_point_to_pyramid_base(
int num_nodes, const double coordinates[][3], double& cos_angle);
double largest_pyramid_edge(const double coordinates[][3]);
/*
the pyramid element
5
^
|\
|\
/| \\_
| \ \
| | \_ \_
......@@ -62,30 +61,28 @@ double largest_pyramid_edge( double coordinates[][3] );
2
a quadrilateral base and a pointy peak like a pyramid
*/
double pyramid_equiangle_skew( int num_nodes, double coordinates[][3] )
double pyramid_equiangle_skew(int /*num_nodes*/, const double coordinates[][3])
{
double base[4][3];
double tri1[3][3];
double tri2[3][3];
double tri3[3][3];
double tri4[3][3];
make_pyramid_faces(coordinates, base,tri1,tri2,tri3,tri4);
make_pyramid_faces(coordinates, base, tri1, tri2, tri3, tri4);
double quad_skew=quad_equiangle_skew( 4, base );
double tri1_skew=tri_equiangle_skew(3,tri1);
double tri2_skew=tri_equiangle_skew(3,tri2);
double tri3_skew=tri_equiangle_skew(3,tri3);
double tri4_skew=tri_equiangle_skew(3,tri4);
double quad_skew = quad_equiangle_skew(4, base);
double tri1_skew = tri_equiangle_skew(3, tri1);
double tri2_skew = tri_equiangle_skew(3, tri2);
double tri3_skew = tri_equiangle_skew(3, tri3);
double tri4_skew = tri_equiangle_skew(3, tri4);
double max_skew=quad_skew;
max_skew = max_skew > tri1_skew ? max_skew : tri1_skew;
max_skew = max_skew > tri2_skew ? max_skew : tri2_skew;
max_skew = max_skew > tri3_skew ? max_skew : tri3_skew;
max_skew = max_skew > tri4_skew ? max_skew : tri4_skew;
double max_skew = quad_skew;
max_skew = max_skew > tri1_skew ? max_skew : tri1_skew;
max_skew = max_skew > tri2_skew ? max_skew : tri2_skew;
max_skew = max_skew > tri3_skew ? max_skew : tri3_skew;
max_skew = max_skew > tri4_skew ? max_skew : tri4_skew;
return max_skew;
}
......@@ -95,19 +92,17 @@ double pyramid_equiangle_skew( int num_nodes, double coordinates[][3] )
the volume is calculated by dividing the pyramid into
2 tets and summing the volumes of the 2 tets.
*/
double pyramid_volume( int num_nodes, double coordinates[][3] )
*/
double pyramid_volume(int /*num_nodes*/, const double coordinates[][3])
{
double center_coords[3];
//calculate the center of the quads
center_coords[0] = (coordinates[0][0] + coordinates[1][0] +
coordinates[2][0] + coordinates[3][0]) / 4;
center_coords[1] = (coordinates[0][1] + coordinates[1][1] +
coordinates[2][1] + coordinates[3][1] ) / 4;
center_coords[2] = (coordinates[0][2] + coordinates[1][2] +
coordinates[2][2] + coordinates[3][2] ) / 4;
// calculate the center of the quads
center_coords[0] =
(coordinates[0][0] + coordinates[1][0] + coordinates[2][0] + coordinates[3][0]) / 4;
center_coords[1] =
(coordinates[0][1] + coordinates[1][1] + coordinates[2][1] + coordinates[3][1]) / 4;
center_coords[2] =
(coordinates[0][2] + coordinates[1][2] + coordinates[2][2] + coordinates[3][2]) / 4;
double tet_coords[4][4][3];
tet_coords[0][0][0] = coordinates[0][0];
......@@ -162,21 +157,15 @@ double pyramid_volume( int num_nodes, double coordinates[][3] )
tet_coords[3][3][1] = coordinates[4][1];
tet_coords[3][3][2] = coordinates[4][2];
double volume=0.0;
for(int t=0;t<4;t++)
double volume = 0.0;
for (int t = 0; t < 4; t++)
{
volume+=tet_volume(4,tet_coords[t]);
volume += tet_volume(4, tet_coords[t]);
}
return (double)volume;
}
double pyramid_jacobian( int num_nodes, double coordinates[][3] )
double pyramid_jacobian(int /*num_nodes*/, const double coordinates[][3])
{
// break the pyramid into four tets return the minimum jacobian of the two tets
double tet1_coords[4][3];
......@@ -197,7 +186,7 @@ double pyramid_jacobian( int num_nodes, double coordinates[][3] )
return p1 < p2 ? p1 : p2;
}
double pyramid_scaled_jacobian( int num_nodes, double coordinates[][3] )
double pyramid_scaled_jacobian(int /*num_nodes*/, const double coordinates[][3])
{
// break the pyramid into four tets return the minimum scaled jacobian of the tets
double tet1_coords[4][3];
......@@ -207,7 +196,7 @@ double pyramid_scaled_jacobian( int num_nodes, double coordinates[][3] )
make_pyramid_tets(coordinates, tet1_coords, tet2_coords, tet3_coords, tet4_coords);
std::array<double,4> scaled_jacob;
std::array<double, 4> scaled_jacob{};
scaled_jacob[0] = tet_scaled_jacobian(4, tet1_coords);
scaled_jacob[1] = tet_scaled_jacobian(4, tet2_coords);
scaled_jacob[2] = tet_scaled_jacobian(4, tet3_coords);
......@@ -218,13 +207,15 @@ double pyramid_scaled_jacobian( int num_nodes, double coordinates[][3] )
// scale the minimum scaled jacobian so that a perfect pyramid has
// a value of 1 and cap it to make sure it is not > 1.0 or < 0.0
if (*iter <= 0.0)
{
return 0.0;
}
double min_jac = (*iter)*2.0/sqrt(2.0);
double min_jac = (*iter) * 2.0 / std::sqrt(2.0);
return min_jac < 1.0 ? min_jac : 1.0 - (min_jac - 1.0);
}
double pyramid_shape( int num_nodes, double coordinates[][3] )
double pyramid_shape(int num_nodes, const double coordinates[][3])
{
// ideally there will be four equilateral triangles and one square.
// Test each face
......@@ -239,29 +230,37 @@ double pyramid_shape( int num_nodes, double coordinates[][3] )
double s1 = quad_shape(4, base);
if (s1 == 0.0)
{
return 0.0;
}
double cos_angle;
double dist_to_base = distance_point_to_pyramid_base(num_nodes, coordinates, cos_angle );
double dist_to_base = distance_point_to_pyramid_base(num_nodes, coordinates, cos_angle);
if (dist_to_base <= 0 || cos_angle <= 0.0)
{
return 0.0;
}
double longest_edge = largest_pyramid_edge(coordinates) * SQRT2_HALVES;
double longest_edge = largest_pyramid_edge(coordinates) * sqrt2_2;
if (dist_to_base < longest_edge)
dist_to_base = dist_to_base/longest_edge;
{
dist_to_base = dist_to_base / longest_edge;
}
else
dist_to_base = longest_edge/dist_to_base;
{
dist_to_base = longest_edge / dist_to_base;
}
double shape = s1 * cos_angle * dist_to_base;
return shape;
}
void make_pyramid_tets(double coordinates[][3], double tet1_coords[][3], double tet2_coords[][3],
double tet3_coords[][3], double tet4_coords[][3])
void make_pyramid_tets(const double coordinates[][3], double tet1_coords[][3],
double tet2_coords[][3], double tet3_coords[][3], double tet4_coords[][3])
{
// tet1
// tet1
tet1_coords[0][0] = coordinates[0][0];
tet1_coords[0][1] = coordinates[0][1];
tet1_coords[0][2] = coordinates[0][2];
......@@ -330,10 +329,10 @@ void make_pyramid_tets(double coordinates[][3], double tet1_coords[][3], double
tet4_coords[3][2] = coordinates[4][2];
}
void make_pyramid_faces(double coordinates[][3], double base[][3], double tri1[][3],
double tri2[][3], double tri3[][3], double tri4[][3])
void make_pyramid_faces(const double coordinates[][3], double base[][3], double tri1[][3],
double tri2[][3], double tri3[][3], double tri4[][3])
{
// base
// base
base[0][0] = coordinates[0][0];
base[0][1] = coordinates[0][1];
base[0][2] = coordinates[0][2];
......@@ -389,7 +388,7 @@ void make_pyramid_faces(double coordinates[][3], double base[][3], double tri1[]
tri3[2][1] = coordinates[4][1];
tri3[2][2] = coordinates[4][2];
// tri4
// tri4
tri4[0][0] = coordinates[3][0];
tri4[0][1] = coordinates[3][1];
tri4[0][2] = coordinates[3][2];
......@@ -403,52 +402,27 @@ void make_pyramid_faces(double coordinates[][3], double base[][3], double tri1[]
tri4[2][2] = coordinates[4][2];
}
void make_pyramid_edges( VerdictVector edges[8], double coordinates[][3] )
void make_pyramid_edges(VerdictVector edges[8], const double coordinates[][3])
{
edges[0].set(
coordinates[1][0] - coordinates[0][0],
coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2]
);
edges[1].set(
coordinates[2][0] - coordinates[1][0],
coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2]
);
edges[2].set(
coordinates[3][0] - coordinates[2][0],
coordinates[3][1] - coordinates[2][1],
coordinates[3][2] - coordinates[2][2]
);
edges[3].set(
coordinates[0][0] - coordinates[3][0],
coordinates[0][1] - coordinates[3][1],
coordinates[0][2] - coordinates[3][2]
);
edges[4].set(
coordinates[4][0] - coordinates[0][0],
coordinates[4][1] - coordinates[0][1],
coordinates[4][2] - coordinates[0][2]
);
edges[5].set(
coordinates[4][0] - coordinates[1][0],
coordinates[4][1] - coordinates[1][1],
coordinates[4][2] - coordinates[1][2]
);
edges[6].set(
coordinates[4][0] - coordinates[2][0],
coordinates[4][1] - coordinates[2][1],
coordinates[4][2] - coordinates[2][2]
);
edges[7].set(
coordinates[4][0] - coordinates[3][0],
coordinates[4][1] - coordinates[3][1],
coordinates[4][2] - coordinates[3][2]
);
edges[0].set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] - coordinates[0][1],
coordinates[1][2] - coordinates[0][2]);
edges[1].set(coordinates[2][0] - coordinates[1][0], coordinates[2][1] - coordinates[1][1],
coordinates[2][2] - coordinates[1][2]);
edges[2].set(coordinates[3][0] - coordinates[2][0], coordinates[3][1] - coordinates[2][1],
coordinates[3][2] - coordinates[2][2]);
edges[3].set(coordinates[0][0] - coordinates[3][0], coordinates[0][1] - coordinates[3][1],
coordinates[0][2] - coordinates[3][2]);
edges[4].set(coordinates[4][0] - coordinates[0][0], coordinates[4][1] - coordinates[0][1],
coordinates[4][2] - coordinates[0][2]);
edges[5].set(coordinates[4][0] - coordinates[1][0], coordinates[4][1] - coordinates[1][1],
coordinates[4][2] - coordinates[1][2]);
edges[6].set(coordinates[4][0] - coordinates[2][0], coordinates[4][1] - coordinates[2][1],
coordinates[4][2] - coordinates[2][2]);
edges[7].set(coordinates[4][0] - coordinates[3][0], coordinates[4][1] - coordinates[3][1],
coordinates[4][2] - coordinates[3][2]);
}
double largest_pyramid_edge( double coordinates[][3] )
double largest_pyramid_edge(const double coordinates[][3])
{
VerdictVector edges[8];
make_pyramid_edges(edges, coordinates);
......@@ -469,10 +443,11 @@ double largest_pyramid_edge( double coordinates[][3] )
max = std::max(max, l6);
max = std::max(max, l7);
return sqrt(max);
return std::sqrt(max);
}
double distance_point_to_pyramid_base( int num_nodes, double coordinates[][3], double &cos_angle )
double distance_point_to_pyramid_base(
int /*num_nodes*/, const double coordinates[][3], double& cos_angle)
{
VerdictVector a(coordinates[0]);
VerdictVector b(coordinates[1]);
......@@ -480,7 +455,7 @@ double distance_point_to_pyramid_base( int num_nodes, double coordinates[][3], d
VerdictVector d(coordinates[3]);
VerdictVector peak(coordinates[4]);
VerdictVector centroid = (a + b + c + d)/4.0;
VerdictVector centroid = (a + b + c + d) / 4.0;
VerdictVector t1 = b - a;
VerdictVector t2 = d - a;
......@@ -489,10 +464,9 @@ double distance_point_to_pyramid_base( int num_nodes, double coordinates[][3], d
VerdictVector pq = peak - centroid;
double distance = (pq % normal)/normal_length;
cos_angle = distance/pq.length();
double distance = (pq % normal) / normal_length;
cos_angle = distance / pq.length();
return distance;
}
} // namespace verdict
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
@PACKAGE_INIT@
# Our library dependencies (contains definitions for IMPORTED targets)
include("${CMAKE_CURRENT_LIST_DIR}/VerdictTargets.cmake")
check_required_components(Verdict)
......@@ -3,13 +3,13 @@
Module: VerdictVector.cpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* VerdictVector.cpp contains implementation of Vector operations
......@@ -18,64 +18,65 @@
*
*/
#include "verdict.h"
#include <math.h>
#include "VerdictVector.hpp"
#include <float.h>
#include "verdict.h"
#if defined(__BORLANDC__)
#pragma warn -8004 /* "assigned a value that is never used" */
#endif
#include <cmath>
namespace VERDICT_NAMESPACE
{
// scale the length of the vector to be the new_length
// unused
//VerdictVector &VerdictVector::length(const double new_length)
//{
// double len = this->length();
// xVal *= new_length / len;
// yVal *= new_length / len;
// zVal *= new_length / len;
// return *this;
//}
VerdictVector& VerdictVector::length(const double new_length)
{
double len = this->length();
xVal *= new_length / len;
yVal *= new_length / len;
zVal *= new_length / len;
return *this;
}
double VerdictVector::interior_angle(const VerdictVector &otherVector)
double VerdictVector::interior_angle(const VerdictVector& otherVector)
{
double cosAngle=0., angleRad=0., len1, len2=0.;
double cosAngle = 0., angleRad = 0., len1, len2 = 0.;
if (((len1 = this->length()) > 0) && ((len2 = otherVector.length()) > 0))
cosAngle = VerdictVector::Dot(*this , otherVector)/(len1 * len2);
{
cosAngle = VerdictVector::Dot(*this, otherVector) / (len1 * len2);
}
else
{
assert(len1 > 0);
assert(len2 > 0);
}
if ((cosAngle > 1.0) && (cosAngle < 1.0001))
{
cosAngle = 1.0;
angleRad = acos(cosAngle);
angleRad = std::acos(cosAngle);
}
else if (cosAngle < -1.0 && cosAngle > -1.0001)
{
cosAngle = -1.0;
angleRad = acos(cosAngle);
angleRad = std::acos(cosAngle);
}
else if (cosAngle >= -1.0 && cosAngle <= 1.0)
angleRad = acos(cosAngle);
{
angleRad = std::acos(cosAngle);
}
else
{
assert(cosAngle < 1.0001 && cosAngle > -1.0001);
}
return( (angleRad * 180.) / VERDICT_PI );
return ((angleRad * 180.) / VERDICT_PI);
}
VerdictVector::VerdictVector(const double xyz[3])
: xVal(xyz[0]), yVal(xyz[1]), zVal(xyz[2])
{}
VerdictVector::VerdictVector(const double xyz[3])
: xVal(xyz[0])
, yVal(xyz[1])
, zVal(xyz[2])
{
}
} // namespace verdict
This diff is collapsed.
/*=========================================================================
Module: V_EdgeMetric.cpp
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS,
the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
#include "gtest/gtest.h"
void InitializeTest( int argc, char **argv );
void CloseTest( void );
void InitializeTest(int argc, char** argv);
void CloseTest(void);
int main(int argc, char **argv)
int main(int argc, char** argv)
{
try
{
......@@ -11,15 +23,14 @@ int main(int argc, char **argv)
int return_status = RUN_ALL_TESTS();
return return_status;
}
catch(std::exception& e)
catch (std::exception& e)
{
printf("Exception Caught \"%s\": Test Failed\n", e.what());
return 1;
}
catch(...)
catch (...)
{
printf("Unexpected Exception Caught: Test Failed\n");
return 1;
}
}
This diff is collapsed.
/*=========================================================================
Module: v_vector.h
Copyright 2003,2006,2019 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
Under the terms of Contract DE-NA0003525 with NTESS, the U.S. Government retains certain rights in this software.
See LICENSE for details.
=========================================================================*/
/*
*
* v_vector.h contains simple vector operations
*
* This file is part of VERDICT
*
*/
// .SECTION Thanks
// Prior to its inclusion within VTK, this code was developed by the CUBIT
// project at Sandia National Laboratories.
#ifndef VERDICT_VECTOR
#define VERDICT_VECTOR
#include "verdict.h"
#include <math.h>
#include <assert.h>
namespace VERDICT_NAMESPACE
{
// computes the dot product of 3d vectors
//double dot_product( double vec1[], double vec2[] );
// computes the cross product
//double *cross_product( double vec1[], double vec2[], double answer[] = 0);
// computes the interior angle between 2 vectors in degrees
//double interior_angle ( double vec1[], double vec2[] );
// computes the length of a vector
//double length ( double vec[] );
//double length_squared (double vec[] );
inline double dot_product( double vec1[], double vec2[] )
{
double answer = vec1[0] * vec2[0] +
vec1[1] * vec2[1] +
vec1[2] * vec2[2];
return answer;
}
inline void normalize( double vec[] )
{
double x = sqrt( vec[0]*vec[0] +
vec[1]*vec[1] +
vec[2]*vec[2] );
vec[0] /= x;
vec[1] /= x;
vec[2] /= x;
}
inline double * cross_product( double vec1[], double vec2[], double answer[] )
{
answer[0] = vec1[1] * vec2[2] - vec1[2] * vec2[1];
answer[1] = vec1[2] * vec2[0] - vec1[0] * vec2[2];
answer[2] = vec1[0] * vec2[1] - vec1[1] * vec2[0];
return answer;
}
inline double length ( double vec[] )
{
return (sqrt ( vec[0] * vec[0] + vec[1] * vec[1] + vec[2] * vec[2] ));
}
inline double length_squared (double vec[] )
{
return (vec[0] * vec[0] + vec[1] * vec[1] + vec[2] * vec[2] );
}
inline double interior_angle( double vec1[], double vec2[] )
{
double len1, len2, cosAngle=0.0, angleRad=0.0;
if ( ((len1 = length(vec1)) > 0 ) && ((len2 = length(vec2)) > 0 ) )
{
cosAngle = dot_product(vec1, vec2) / (len1 * len2);
}
else
{
assert(len1 > 0);
assert(len2 > 0);
}
if ((cosAngle > 1.0) && (cosAngle < 1.0001))
{
cosAngle = 1.0;
angleRad = acos(cosAngle);
}
else if (cosAngle < -1.0 && cosAngle > -1.0001)
{
cosAngle = -1.0;
angleRad = acos(cosAngle);
}
else if (cosAngle >= -1.0 && cosAngle <= 1.0)
angleRad = acos(cosAngle);
else
{
assert(cosAngle < 1.0001 && cosAngle > -1.0001);
}
return( (angleRad * 180.) / VERDICT_PI );
}
} // namespace verdict
#endif
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment