Script 'mail_helper' called by obssrc
Hello community,
here is the log from the commit of package verdict for openSUSE:Factory checked
in at 2023-11-10 12:33:23
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Comparing /work/SRC/openSUSE:Factory/verdict (Old)
and /work/SRC/openSUSE:Factory/.verdict.new.17445 (New)
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Package is "verdict"
Fri Nov 10 12:33:23 2023 rev:3 rq:1124784 version:1.4.2
Changes:
--------
--- /work/SRC/openSUSE:Factory/verdict/verdict.changes 2022-12-07
17:36:19.256984599 +0100
+++ /work/SRC/openSUSE:Factory/.verdict.new.17445/verdict.changes
2023-11-10 12:37:18.369333620 +0100
@@ -1,0 +2,11 @@
+Tue Oct 10 08:30:58 UTC 2023 - Dirk Müller <dmuel...@suse.com>
+
+- update to 1.4.2:
+ * Add interfaces to some verdict calculations that take double **
+ * Implementing length/area calculations and include higher order
+ nodes
+ * add verdict double ** interface for normalized in-radius
+ * More changes to tri normalized inradius to work for dimension 2
+ * fix tri normalized inradius to not return nan for collapsd tri
+
+-------------------------------------------------------------------
Old:
----
verdict-1.4.1.tar.gz
New:
----
verdict-1.4.2.tar.gz
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Other differences:
------------------
++++++ verdict.spec ++++++
--- /var/tmp/diff_new_pack.ocDN5h/_old 2023-11-10 12:37:19.069359528 +0100
+++ /var/tmp/diff_new_pack.ocDN5h/_new 2023-11-10 12:37:19.069359528 +0100
@@ -1,7 +1,7 @@
#
# spec file for package verdict
#
-# Copyright (c) 2022 SUSE LLC
+# Copyright (c) 2023 SUSE LLC
#
# All modifications and additions to the file contributed by third parties
# remain the property of their copyright owners, unless otherwise agreed
@@ -19,7 +19,7 @@
%define libname libverdict1_4
Name: verdict
-Version: 1.4.1
+Version: 1.4.2
Release: 0
Summary: Compute quality functions of 2 and 3-dimensional regions
License: BSD-3-Clause
++++++ verdict-1.4.1.tar.gz -> verdict-1.4.2.tar.gz ++++++
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/CMakeLists.txt
new/verdict-1.4.2/CMakeLists.txt
--- old/verdict-1.4.1/CMakeLists.txt 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/CMakeLists.txt 2023-07-19 18:21:57.000000000 +0200
@@ -7,7 +7,7 @@
set(CMAKE_CXX_STANDARD 11)
set(CMAKE_CXX_STANDARD_REQUIRED 1)
-project(verdict VERSION 1.4.1)
+project(verdict VERSION 1.4.2)
# Includes
include(GNUInstallDirs)
@@ -17,8 +17,7 @@
cmake_policy(SET CMP0063 NEW)
endif ()
-set(verdict_VERSION_FLAT
"${CMAKE_PROJECT_VERSION_MAJOR}${CMAKE_PROJECT_VERSION_MINOR}${CMAKE_PROJECT_VERSION_PATCH}")
-set(verdict_VERSION "${CMAKE_PROJECT_VERSION}")
+set(verdict_VERSION_FLAT
"${verdict_VERSION_MAJOR}${verdict_VERSION_MINOR}${verdict_VERSION_PATCH}")
option(VERDICT_BUILD_DOC "Build the 2007 Verdict User Manual" OFF)
option(VERDICT_MANGLE "Mangle verdict names for inclusion in a larger
library?" OFF)
@@ -53,7 +52,7 @@
${CMAKE_CURRENT_BINARY_DIR}/verdict_config.h
@ONLY)
-add_library(verdict ${verdict_SOURCES})
+add_library(verdict ${verdict_SOURCES} ${verdict_HEADERS}
${CMAKE_CURRENT_BINARY_DIR}/verdict_config.h)
target_include_directories(verdict PUBLIC
$<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}>
$<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}>)
if(UNIX)
@@ -68,7 +67,7 @@
if (NOT VERDICT_NO_LIBRARY_VERSION)
set_target_properties(verdict PROPERTIES
VERSION "${verdict_VERSION}"
- SOVERSION "${CMAKE_PROJECT_VERSION_MAJOR}.${CMAKE_PROJECT_VERSION_MINOR}")
+ SOVERSION "${verdict_VERSION_MAJOR}.${verdict_VERSION_MINOR}")
endif ()
if (NOT VERDICT_EXPORT_GROUP)
set(VERDICT_EXPORT_GROUP VerdictExport)
@@ -139,7 +138,7 @@
# Export Targets file
export(EXPORT ${VERDICT_EXPORT_GROUP}
- FILE "${CMAKE_CURRENT_BINARY_DIR}/cmake/VerdictTargets.cmake"
+ FILE "${CMAKE_CURRENT_BINARY_DIR}/VerdictTargets.cmake"
NAMESPACE Verdict::)
# Install Targets file
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/V_EdgeMetric.cpp
new/verdict-1.4.2/V_EdgeMetric.cpp
--- old/verdict-1.4.1/V_EdgeMetric.cpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/V_EdgeMetric.cpp 2023-07-19 18:21:57.000000000 +0200
@@ -27,12 +27,30 @@
/*!\brief Length of and edge.
* Length is calculated by taking the distance between the end nodes.
*/
-double edge_length(int /*num_nodes*/, const double coordinates[][3])
+double edge_length(int num_nodes, const double coordinates[][3])
{
+ double edge_length = 0.0;
- 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)(std::sqrt(x * x + y * y + z * z));
+ if (2 == num_nodes)
+ {
+ double x = coordinates[1][0] - coordinates[0][0];
+ double y = coordinates[1][1] - coordinates[0][1];
+ double z = coordinates[1][2] - coordinates[0][2];
+ edge_length = (double)(std::sqrt(x * x + y * y + z * z));
+ }
+ if (3 == num_nodes)
+ {
+ double x = coordinates[2][0] - coordinates[0][0];
+ double y = coordinates[2][1] - coordinates[0][1];
+ double z = coordinates[2][2] - coordinates[0][2];
+ edge_length += (double)(std::sqrt(x * x + y * y + z * z));
+
+ x = coordinates[2][0] - coordinates[1][0];
+ y = coordinates[2][1] - coordinates[1][1];
+ z = coordinates[2][2] - coordinates[1][2];
+ edge_length += (double)(std::sqrt(x * x + y * y + z * z));
+ }
+
+ return edge_length;
}
} // namespace verdict
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/V_HexMetric.cpp
new/verdict-1.4.2/V_HexMetric.cpp
--- old/verdict-1.4.1/V_HexMetric.cpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/V_HexMetric.cpp 2023-07-19 18:21:57.000000000 +0200
@@ -105,7 +105,7 @@
// Compute interior node
VerdictVector hex20_auxillary_node_coordinate(const double coordinates[][3])
{
- VerdictVector aux_node(0, 0, 0);
+ VerdictVector aux_node(0.0, 0.0, 0.0);
for (int i = 0; i < 8; i++)
{
VerdictVector tmp_vec(coordinates[i][0], coordinates[i][1],
coordinates[i][2]);
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/V_QuadMetric.cpp
new/verdict-1.4.2/V_QuadMetric.cpp
--- old/verdict-1.4.1/V_QuadMetric.cpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/V_QuadMetric.cpp 2023-07-19 18:21:57.000000000 +0200
@@ -24,6 +24,7 @@
#include "verdict_defines.hpp"
#include <algorithm>
+#include <vector>
namespace VERDICT_NAMESPACE
{
@@ -738,18 +739,103 @@
jacobian at quad center
*/
-double quad_area(int /*num_nodes*/, const double coordinates[][3])
+double quad_area(int num_nodes, const double coordinates[][3])
{
- double corner_areas[4];
- signed_corner_areas(corner_areas, coordinates);
+ if (4 == num_nodes)
+ {
+ double corner_areas[4];
+ signed_corner_areas(corner_areas, coordinates);
- double area = 0.25 * (corner_areas[0] + corner_areas[1] + corner_areas[2] +
corner_areas[3]);
+ double area = 0.25 * (corner_areas[0] + corner_areas[1] + corner_areas[2]
+ corner_areas[3]);
- if (area > 0)
+ if (area > 0)
+ {
+ return (double)std::min(area, VERDICT_DBL_MAX);
+ }
+ return (double)std::max(area, -VERDICT_DBL_MAX);
+ }
+ else
{
- return (double)std::min(area, VERDICT_DBL_MAX);
+ double area = 0;
+ double tmp_coords[4][3];
+
+ if (5 == num_nodes)
+ {
+ std::vector< std::vector<int> > tri_conn = { {0,1}, {1,2}, {2,3}, {3,0}
};
+
+ //center node 4
+ tmp_coords[2][0] = coordinates[4][0];
+ tmp_coords[2][1] = coordinates[4][1];
+ tmp_coords[2][2] = coordinates[4][2];
+
+ for (std::vector<int> v : tri_conn)
+ {
+ tmp_coords[0][0] = coordinates[v[0]][0];
+ tmp_coords[0][1] = coordinates[v[0]][1];
+ tmp_coords[0][2] = coordinates[v[0]][2];
+ tmp_coords[1][0] = coordinates[v[1]][0];
+ tmp_coords[1][1] = coordinates[v[1]][1];
+ tmp_coords[1][2] = coordinates[v[1]][2];
+ area += tri_area(3, tmp_coords);
+ }
+ }
+ else if (8 == num_nodes)
+ {
+ std::vector< std::vector<int> > tri_conn =
+ { {0,4,7}, {4,1,5}, {5,2,6}, {6,3,7} };
+
+ for (std::vector<int> v : tri_conn)
+ {
+ tmp_coords[0][0] = coordinates[v[0]][0];
+ tmp_coords[0][1] = coordinates[v[0]][1];
+ tmp_coords[0][2] = coordinates[v[0]][2];
+ tmp_coords[1][0] = coordinates[v[1]][0];
+ tmp_coords[1][1] = coordinates[v[1]][1];
+ tmp_coords[1][2] = coordinates[v[1]][2];
+ tmp_coords[2][0] = coordinates[v[2]][0];
+ tmp_coords[2][1] = coordinates[v[2]][1];
+ tmp_coords[2][2] = coordinates[v[2]][2];
+ area += tri_area(3, tmp_coords);
+ }
+
+ //interior quad 4567
+ tmp_coords[0][0] = coordinates[4][0];
+ tmp_coords[0][1] = coordinates[4][1];
+ tmp_coords[0][2] = coordinates[4][2];
+ tmp_coords[1][0] = coordinates[5][0];
+ tmp_coords[1][1] = coordinates[5][1];
+ tmp_coords[1][2] = coordinates[5][2];
+ tmp_coords[2][0] = coordinates[6][0];
+ tmp_coords[2][1] = coordinates[6][1];
+ tmp_coords[2][2] = coordinates[6][2];
+ tmp_coords[3][0] = coordinates[7][0];
+ tmp_coords[3][1] = coordinates[7][1];
+ tmp_coords[3][2] = coordinates[7][2];
+ area += quad_area(4, tmp_coords);
+ }
+ else if (9 == num_nodes)
+ {
+ std::vector< std::vector<int> > tri_conn =
+ { {0,4}, {4,1}, {1,5}, {5,2}, {2,6}, {6,3}, {3,7}, {7,0} };
+
+ //quad center node
+ tmp_coords[2][0] = coordinates[8][0];
+ tmp_coords[2][1] = coordinates[8][1];
+ tmp_coords[2][2] = coordinates[8][2];
+
+ for (std::vector<int> v : tri_conn)
+ {
+ tmp_coords[0][0] = coordinates[v[0]][0];
+ tmp_coords[0][1] = coordinates[v[0]][1];
+ tmp_coords[0][2] = coordinates[v[0]][2];
+ tmp_coords[1][0] = coordinates[v[1]][0];
+ tmp_coords[1][1] = coordinates[v[1]][1];
+ tmp_coords[1][2] = coordinates[v[1]][2];
+ area += tri_area(3, tmp_coords);
+ }
+ }
+ return area;
}
- return (double)std::max(area, -VERDICT_DBL_MAX);
}
/*!
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/V_TetMetric.cpp
new/verdict-1.4.2/V_TetMetric.cpp
--- old/verdict-1.4.1/V_TetMetric.cpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/V_TetMetric.cpp 2023-07-19 18:21:57.000000000 +0200
@@ -278,20 +278,16 @@
minimum of the jacobian divided by the lengths of 3 edge vectors
*/
-double tet_scaled_jacobian(int /*num_nodes*/, const double coordinates[][3])
+
+template <typename CoordsContainerType>
+double tet_scaled_jacobian_impl(int /*num_nodes*/, const CoordsContainerType
coordinates)
{
- const VerdictVector side0(coordinates[1][0] - coordinates[0][0],
- coordinates[1][1] - coordinates[0][1], coordinates[1][2] -
coordinates[0][2]);
- const VerdictVector side1(coordinates[2][0] - coordinates[1][0],
- coordinates[2][1] - coordinates[1][1], coordinates[2][2] -
coordinates[1][2]);
- const VerdictVector side2(coordinates[0][0] - coordinates[2][0],
- coordinates[0][1] - coordinates[2][1], coordinates[0][2] -
coordinates[2][2]);
- const VerdictVector side3(coordinates[3][0] - coordinates[0][0],
- coordinates[3][1] - coordinates[0][1], coordinates[3][2] -
coordinates[0][2]);
- const VerdictVector side4(coordinates[3][0] - coordinates[1][0],
- coordinates[3][1] - coordinates[1][1], coordinates[3][2] -
coordinates[1][2]);
- const VerdictVector side5(coordinates[3][0] - coordinates[2][0],
- coordinates[3][1] - coordinates[2][1], coordinates[3][2] -
coordinates[2][2]);
+ const VerdictVector side0{coordinates[0], coordinates[1]};
+ const VerdictVector side1{coordinates[1], coordinates[2]};
+ const VerdictVector side2{coordinates[2], coordinates[0]};
+ const VerdictVector side3{coordinates[0], coordinates[3]};
+ const VerdictVector side4{coordinates[1], coordinates[3]};
+ const VerdictVector side5{coordinates[2], coordinates[3]};
const double jacobi = side3 % (side2 * side0);
@@ -303,8 +299,8 @@
const double side4_length_squared = side4.length_squared();
const double side5_length_squared = side5.length_squared();
- const double length_squared[4] = { side0_length_squared *
side2_length_squared *
- side3_length_squared,
+ const double length_squared[4] = {
+ side0_length_squared * side2_length_squared * side3_length_squared,
side0_length_squared * side1_length_squared * side4_length_squared,
side1_length_squared * side2_length_squared * side5_length_squared,
side3_length_squared * side4_length_squared * side5_length_squared };
@@ -336,6 +332,15 @@
return (double)(sqrt2 * jacobi / length_product);
}
+double tet_scaled_jacobian(int num_nodes, const double coordinates[][3])
+{
+ return tet_scaled_jacobian_impl(num_nodes, coordinates);
+}
+
+double tet_scaled_jacobian_from_loc_ptrs(int num_nodes, const double * const *
coordinates)
+{
+ return tet_scaled_jacobian_impl(num_nodes, coordinates);
+}
/*!
The radius ratio of a tet
@@ -401,19 +406,13 @@
conditioned. Previously, this would only happen when the volume was small
and positive, but now ill-conditioned inverted tetrahedra are also
included.
*/
-double tet_aspect_ratio(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet_aspect_ratio_impl(int /*num_nodes*/, const CoordsContainerType
coordinates)
{
// Determine side vectors
- VerdictVector ab, bc, ac, ad, bd, cd;
-
- ab.set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- ac.set(coordinates[2][0] - coordinates[0][0], coordinates[2][1] -
coordinates[0][1],
- coordinates[2][2] - coordinates[0][2]);
-
- ad.set(coordinates[3][0] - coordinates[0][0], coordinates[3][1] -
coordinates[0][1],
- coordinates[3][2] - coordinates[0][2]);
+ const VerdictVector ab{coordinates[0], coordinates[1]};
+ const VerdictVector ac{coordinates[0], coordinates[2]};
+ const VerdictVector ad{coordinates[0], coordinates[3]};
double detTet = ab % (ac * ad);
@@ -422,27 +421,22 @@
return (double)VERDICT_DBL_MAX;
}
- bc.set(coordinates[2][0] - coordinates[1][0], coordinates[2][1] -
coordinates[1][1],
- coordinates[2][2] - coordinates[1][2]);
-
- bd.set(coordinates[3][0] - coordinates[1][0], coordinates[3][1] -
coordinates[1][1],
- coordinates[3][2] - coordinates[1][2]);
-
- cd.set(coordinates[3][0] - coordinates[2][0], coordinates[3][1] -
coordinates[2][1],
- coordinates[3][2] - coordinates[2][2]);
-
- double ab2 = ab.length_squared();
- double bc2 = bc.length_squared();
- double ac2 = ac.length_squared();
- double ad2 = ad.length_squared();
- double bd2 = bd.length_squared();
- double cd2 = cd.length_squared();
+ VerdictVector bc{coordinates[1], coordinates[2]};
+ VerdictVector bd{coordinates[1], coordinates[3]};
+ VerdictVector cd{coordinates[2], coordinates[3]};
+
+ const double ab2 = ab.length_squared();
+ const double bc2 = bc.length_squared();
+ const double ac2 = ac.length_squared();
+ const double ad2 = ad.length_squared();
+ const double bd2 = bd.length_squared();
+ const double cd2 = cd.length_squared();
double A = ab2 > bc2 ? ab2 : bc2;
double B = ac2 > ad2 ? ac2 : ad2;
double C = bd2 > cd2 ? bd2 : cd2;
double D = A > B ? A : B;
- double hm = D > C ? std::sqrt(D) : std::sqrt(C);
+ const double hm = D > C ? std::sqrt(D) : std::sqrt(C);
bd = ab * bc;
A = bd.length();
@@ -458,6 +452,15 @@
return fix_range(aspect_ratio);
}
+double tet_aspect_ratio(int num_nodes, const double coordinates[][3])
+{
+ return tet_aspect_ratio_impl(num_nodes, coordinates);
+}
+double tet_aspect_ratio_from_loc_ptrs(int num_nodes, const double * const
*coordinates)
+{
+ return tet_aspect_ratio_impl(num_nodes, coordinates);
+}
+
/*!
the aspect gamma of a tet
@@ -954,25 +957,21 @@
1/6 * jacobian at a corner node
*/
-double tet_volume(int num_nodes, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet_volume_impl(int num_nodes, const CoordsContainerType coordinates)
{
// Determine side vectors
- VerdictVector side0, side2, side3;
if (4 == num_nodes)
{
- side2.set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- side0.set(coordinates[2][0] - coordinates[0][0], coordinates[2][1] -
coordinates[0][1],
- coordinates[2][2] - coordinates[0][2]);
-
- side3.set(coordinates[3][0] - coordinates[0][0], coordinates[3][1] -
coordinates[0][1],
- coordinates[3][2] - coordinates[0][2]);
+ const VerdictVector side2{coordinates[0], coordinates[1]};
+ const VerdictVector side0{coordinates[0], coordinates[2]};
+ const VerdictVector side3{coordinates[0], coordinates[3]};
return calculate_tet_volume_using_sides(side0, side2, side3);
}
else
{
VerdictVector tet_pts[15];
+ VerdictVector side0, side2, side3;
// create a vector for each point
for (int k = 0; k < num_nodes; k++)
@@ -981,7 +980,7 @@
}
// determine center point of the higher-order nodes
- VerdictVector centroid(0, 0, 0);
+ VerdictVector centroid(0.0, 0.0, 0.0);
for (int k = 4; k < num_nodes; k++)
{
centroid += VerdictVector(coordinates[k][0], coordinates[k][1],
coordinates[k][2]);
@@ -1133,6 +1132,16 @@
return 0;
}
+double tet_volume(int num_nodes, const double coordinates[][3])
+{
+ return tet_volume_impl(num_nodes, coordinates);
+}
+
+double tet_volume_from_loc_ptrs(int num_nodes, const double * const
*coordinates)
+{
+ return tet_volume_impl(num_nodes, coordinates);
+}
+
/*!
the condition of a tet
@@ -1143,30 +1152,20 @@
conditioned. Previously, this would only happen when the volume was small
and positive, but now ill-conditioned inverted tetrahedra are also
included.
*/
-double tet_condition(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet_condition_impl(int /*num_nodes*/, const CoordsContainerType
coordinates)
{
- double condition, term1, term2, det;
-
- VerdictVector side0, side2, side3;
-
- side0.set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- side2.set(coordinates[0][0] - coordinates[2][0], coordinates[0][1] -
coordinates[2][1],
- coordinates[0][2] - coordinates[2][2]);
-
- side3.set(coordinates[3][0] - coordinates[0][0], coordinates[3][1] -
coordinates[0][1],
- coordinates[3][2] - coordinates[0][2]);
-
- VerdictVector c_1, c_2, c_3;
-
- c_1 = side0;
- c_2 = (-2 * side2 - side0) / sqrt3;
- c_3 = (3 * side3 + side2 - side0) / sqrt6;
-
- term1 = c_1 % c_1 + c_2 % c_2 + c_3 % c_3;
- term2 = (c_1 * c_2) % (c_1 * c_2) + (c_2 * c_3) % (c_2 * c_3) + (c_1 * c_3)
% (c_1 * c_3);
- det = c_1 % (c_2 * c_3);
+ const VerdictVector side0{coordinates[0], coordinates[1]};
+ const VerdictVector side2{coordinates[2], coordinates[0]};
+ const VerdictVector side3{coordinates[0], coordinates[3]};
+
+ const VerdictVector c_1 = side0;
+ const VerdictVector c_2 = (-2 * side2 - side0) / sqrt3;
+ const VerdictVector c_3 = (3 * side3 + side2 - side0) / sqrt6;
+
+ const double term1 = c_1 % c_1 + c_2 % c_2 + c_3 % c_3;
+ const double term2 = (c_1 * c_2) % (c_1 * c_2) + (c_2 * c_3) % (c_2 * c_3) +
(c_1 * c_3) % (c_1 * c_3);
+ const double det = c_1 % (c_2 * c_3);
if (std::abs(det) <= VERDICT_DBL_MIN)
{
@@ -1174,12 +1173,19 @@
}
else
{
- condition = std::sqrt(term1 * term2) / (3.0 * det);
+ return std::sqrt(term1 * term2) / (3.0 * det);
}
+}
- return (double)condition;
+double tet_condition(int num_nodes, const double coordinates[][3])
+{
+ return tet_condition_impl(num_nodes, coordinates);
}
+double tet_condition_from_loc_ptrs(int num_nodes, const double * const
*coordinates)
+{
+ return tet_condition_impl(num_nodes, coordinates);
+}
/*!
the jacobian of a tet
@@ -1517,9 +1523,10 @@
return char_length / denominator;
}
-VerdictVector tet10_auxillary_node_coordinate(const double coordinates[][3])
+template <typename CoordsContainerType>
+VerdictVector tet10_auxillary_node_coordinate(const CoordsContainerType
coordinates)
{
- VerdictVector aux_node(0, 0, 0);
+ VerdictVector aux_node(0.0, 0.0, 0.0);
for (int i = 4; i < 10; i++)
{
VerdictVector tmp_vec(coordinates[i][0], coordinates[i][1],
coordinates[i][2]);
@@ -1530,7 +1537,8 @@
return aux_node;
}
-double tet10_min_inradius(const double coordinates[][3], int begin_index, int
end_index)
+template <typename CoordsContainerType>
+double tet10_min_inradius(const CoordsContainerType coordinates, int
begin_index, int end_index)
{
double min_tetinradius = VERDICT_DBL_MAX;
@@ -1582,7 +1590,9 @@
return min_tetinradius;
}
-double calculate_tet4_outer_radius(const double coordinates[][3])
+
+template <typename CoordsContainerType>
+double calculate_tet4_outer_radius(const CoordsContainerType coordinates)
{
verdict::VerdictVector nE[4];
for (int i{ 0 }; i < 4; i++)
@@ -1604,7 +1614,8 @@
return outer_radius;
}
-double tet10_normalized_inradius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet10_normalized_inradius(const CoordsContainerType coordinates)
{
double min_inradius_for_subtet_with_parent_node =
tet10_min_inradius(coordinates, 0, 3);
double min_inradius_for_subtet_with_no_parent_node =
tet10_min_inradius(coordinates, 4, 11);
@@ -1621,7 +1632,8 @@
return fix_range(norm_inrad);
}
-double tet4_normalized_inradius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet4_normalized_inradius(const CoordsContainerType coordinates)
{
double tet10_coords[10][3];
for (int i = 0; i < 4; i++)
@@ -1643,7 +1655,8 @@
return tet10_normalized_inradius(tet10_coords);
}
-double tet_normalized_inradius(int num_nodes, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet_normalized_inradius_impl(int num_nodes, const CoordsContainerType
coordinates)
{
if (num_nodes == 4)
{
@@ -1656,16 +1669,22 @@
return 0.0;
}
-double tet4_mean_ratio(const double coordinates[][3])
+double tet_normalized_inradius(int num_nodes, const double coordinates[][3])
{
- const VerdictVector side0(coordinates[1][0] - coordinates[0][0],
- coordinates[1][1] - coordinates[0][1], coordinates[1][2] -
coordinates[0][2]);
+ return tet_normalized_inradius_impl(num_nodes, coordinates);
+}
- const VerdictVector side2(coordinates[0][0] - coordinates[2][0],
- coordinates[0][1] - coordinates[2][1], coordinates[0][2] -
coordinates[2][2]);
+double tet_normalized_inradius_from_loc_ptrs(int num_nodes, const double *
const *coordinates)
+{
+ return tet_normalized_inradius_impl(num_nodes, coordinates);
+}
- const VerdictVector side3(coordinates[3][0] - coordinates[0][0],
- coordinates[3][1] - coordinates[0][1], coordinates[3][2] -
coordinates[0][2]);
+template <typename CoordsContainerType>
+double tet4_mean_ratio(const CoordsContainerType coordinates)
+{
+ const VerdictVector side0{coordinates[0], coordinates[1]};
+ const VerdictVector side2{coordinates[2], coordinates[0]};
+ const VerdictVector side3{coordinates[0], coordinates[3]};
const double tetVolume = calculate_tet_volume_using_sides(side0, side2,
side3);
if (std::abs( tetVolume ) < VERDICT_DBL_MIN)
@@ -1673,14 +1692,9 @@
return 0.0;
}
- const VerdictVector side1(coordinates[2][0] - coordinates[1][0],
- coordinates[2][1] - coordinates[1][1], coordinates[2][2] -
coordinates[1][2]);
-
- const VerdictVector side4(coordinates[3][0] - coordinates[1][0],
- coordinates[3][1] - coordinates[1][1], coordinates[3][2] -
coordinates[1][2]);
-
- const VerdictVector side5(coordinates[3][0] - coordinates[2][0],
- coordinates[3][1] - coordinates[2][1], coordinates[3][2] -
coordinates[2][2]);
+ const VerdictVector side1{coordinates[1], coordinates[2]};
+ const VerdictVector side4{coordinates[1], coordinates[3]};
+ const VerdictVector side5{coordinates[2], coordinates[3]};
const double side0_length_squared = side0.length_squared();
const double side1_length_squared = side1.length_squared();
@@ -1696,7 +1710,8 @@
return 6 * std::pow(2, 0.5) * tetVolume / std::pow(sum, 3. / 2.);
}
-double tet10_mean_ratio(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet10_mean_ratio(const CoordsContainerType coordinates)
{
double min_tet_mean_ratio = VERDICT_DBL_MAX;
@@ -1747,7 +1762,8 @@
return min_tet_mean_ratio;
}
-double tet_mean_ratio(int num_nodes, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tet_mean_ratio_impl(int num_nodes, const CoordsContainerType
coordinates)
{
if (num_nodes == 4)
{
@@ -1759,4 +1775,14 @@
}
return 0.0;
}
+
+double tet_mean_ratio(int num_nodes, const double coordinates[][3])
+{
+ return tet_mean_ratio_impl(num_nodes, coordinates);
+}
+
+double tet_mean_ratio_from_loc_ptrs(int num_nodes, const double * const *
coordinates)
+{
+ return tet_mean_ratio_impl(num_nodes, coordinates);
+}
} // namespace verdict
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/V_TriMetric.cpp
new/verdict-1.4.2/V_TriMetric.cpp
--- old/verdict-1.4.1/V_TriMetric.cpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/V_TriMetric.cpp 2023-07-19 18:21:57.000000000 +0200
@@ -24,6 +24,7 @@
#include "verdict_defines.hpp"
#include <algorithm>
+#include <array>
namespace VERDICT_NAMESPACE
{
@@ -149,27 +150,23 @@
what is now called "v_tri_aspect_frobenius"
*/
-double tri_aspect_ratio(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_aspect_ratio_impl(int /*num_nodes*/, const CoordsContainerType
coordinates, const int dimension)
{
// three vectors for each side
- VerdictVector a(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- VerdictVector b(coordinates[2][0] - coordinates[1][0], coordinates[2][1] -
coordinates[1][1],
- coordinates[2][2] - coordinates[1][2]);
-
- VerdictVector c(coordinates[0][0] - coordinates[2][0], coordinates[0][1] -
coordinates[2][1],
- coordinates[0][2] - coordinates[2][2]);
-
- double a1 = a.length();
- double b1 = b.length();
- double c1 = c.length();
+ const VerdictVector a{coordinates[0], coordinates[1], dimension};
+ const VerdictVector b{coordinates[1], coordinates[2], dimension};
+ const VerdictVector c{coordinates[2], coordinates[0], dimension};
+
+ const double a1 = a.length();
+ const double b1 = b.length();
+ const double c1 = c.length();
double hm = a1 > b1 ? a1 : b1;
hm = hm > c1 ? hm : c1;
- VerdictVector ab = a * b;
- double denominator = ab.length();
+ const VerdictVector ab = a * b;
+ const double denominator = ab.length();
if (denominator < VERDICT_DBL_MIN)
{
@@ -177,9 +174,7 @@
}
else
{
- double aspect_ratio;
- aspect_ratio = aspect_ratio_normal_coeff * hm * (a1 + b1 + c1) /
denominator;
-
+ const double aspect_ratio = aspect_ratio_normal_coeff * hm * (a1 + b1 +
c1) / denominator;
if (aspect_ratio > 0)
{
return (double)std::min(aspect_ratio, VERDICT_DBL_MAX);
@@ -188,6 +183,14 @@
}
}
+double tri_aspect_ratio(int num_nodes, const double coordinates[][3])
+{
+ return tri_aspect_ratio_impl(num_nodes, coordinates, 3);
+}
+double tri_aspect_ratio_from_loc_ptrs(int num_nodes, const double * const
*coordinates, const int dimension)
+{
+ return tri_aspect_ratio_impl(num_nodes, coordinates, dimension);
+}
/*!
the radius ratio of a triangle
@@ -279,26 +282,124 @@
0.5 * jacobian at a node
*/
-double tri_area(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_area_impl(int num_nodes, const CoordsContainerType coordinates,
const int dimension)
{
- // two vectors for two sides
- VerdictVector side1(coordinates[1][0] - coordinates[0][0], coordinates[1][1]
- coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
+ if (3 == num_nodes)
+ {
+ // two vectors for two sides
+ const VerdictVector side1{ coordinates[0], coordinates[1], dimension };
+ const VerdictVector side3{ coordinates[0], coordinates[2], dimension };
+
+ // the cross product of the two vectors representing two sides of the
+ // triangle
+ const VerdictVector tmp = side1 * side3;
+
+ // return the magnitude of the vector divided by two
+ const double area = 0.5 * tmp.length();
+ if (area > 0)
+ {
+ return (double)std::min(area, VERDICT_DBL_MAX);
+ }
+ return (double)std::max(area, -VERDICT_DBL_MAX);
+ }
+ else
+ {
+ const double *tmp_coords[3];
+ double tri_area = 0.0;
- VerdictVector side3(coordinates[2][0] - coordinates[0][0], coordinates[2][1]
- coordinates[0][1],
- coordinates[2][2] - coordinates[0][2]);
+ if (6 == num_nodes)
+ {
+ // 035
+ tmp_coords[0] = coordinates[0];
+ tmp_coords[1] = coordinates[3];
+ tmp_coords[2] = coordinates[5];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 314
+ tmp_coords[0] = coordinates[3];
+ tmp_coords[1] = coordinates[1];
+ tmp_coords[2] = coordinates[4];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 425
+ tmp_coords[0] = coordinates[4];
+ tmp_coords[1] = coordinates[2];
+ tmp_coords[2] = coordinates[5];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 345
+ tmp_coords[0] = coordinates[3];
+ tmp_coords[1] = coordinates[4];
+ tmp_coords[2] = coordinates[5];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+ }
+ else if (7 == num_nodes)
+ {
+ //center node 6
+ tmp_coords[2] = coordinates[6];
+
+ // 036
+ tmp_coords[0] = coordinates[0];
+ tmp_coords[1] = coordinates[3];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 316
+ tmp_coords[0] = coordinates[3];
+ tmp_coords[1] = coordinates[1];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 146
+ tmp_coords[0] = coordinates[1];
+ tmp_coords[1] = coordinates[4];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 426
+ tmp_coords[0] = coordinates[4];
+ tmp_coords[1] = coordinates[2];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 256
+ tmp_coords[0] = coordinates[2];
+ tmp_coords[1] = coordinates[5];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 506
+ tmp_coords[0] = coordinates[5];
+ tmp_coords[1] = coordinates[0];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+ }
+ else if( 4 == num_nodes )
+ {
+ //center node 3
+ tmp_coords[2] = coordinates[3];
- // the cross product of the two vectors representing two sides of the
- // triangle
- VerdictVector tmp = side1 * side3;
+ // 013
+ tmp_coords[0] = coordinates[0];
+ tmp_coords[1] = coordinates[1];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
- // return the magnitude of the vector divided by two
- double area = 0.5 * tmp.length();
- if (area > 0)
- {
- return (double)std::min(area, VERDICT_DBL_MAX);
+ // 123
+ tmp_coords[0] = coordinates[1];
+ tmp_coords[1] = coordinates[2];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+
+ // 203
+ tmp_coords[0] = coordinates[2];
+ tmp_coords[1] = coordinates[0];
+ tri_area += tri_area_impl(3, tmp_coords, dimension);
+ }
+ return tri_area;
}
- return (double)std::max(area, -VERDICT_DBL_MAX);
+}
+double tri_area(int num_nodes, const double coordinates[][3])
+{
+ return tri_area_impl(num_nodes, coordinates, 3);
+}
+
+double tri_area_from_loc_ptrs(int num_nodes, const double * const
*coordinates, const int dimension)
+{
+ return tri_area_impl(num_nodes, coordinates, dimension);
}
/*!
@@ -463,54 +564,56 @@
Condition number of the jacobian matrix at any corner
*/
-double tri_condition(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_condition_impl(int /*num_nodes*/, const CoordsContainerType
coordinates, const int dimension)
{
- VerdictVector v1(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- VerdictVector v2(coordinates[2][0] - coordinates[0][0], coordinates[2][1] -
coordinates[0][1],
- coordinates[2][2] - coordinates[0][2]);
-
- VerdictVector tri_normal = v1 * v2;
- double areax2 = tri_normal.length();
+ const VerdictVector v1{coordinates[0], coordinates[1], dimension};
+ const VerdictVector v2{coordinates[0], coordinates[2], dimension};
+ const VerdictVector tri_normal = v1 * v2;
+ const double areax2 = tri_normal.length();
if (areax2 == 0.0)
{
return (double)VERDICT_DBL_MAX;
}
- double condition = (double)(((v1 % v1) + (v2 % v2) - (v1 % v2)) / (areax2 *
sqrt3));
+ const double condition = (double)(((v1 % v1) + (v2 % v2) - (v1 % v2)) /
(areax2 * sqrt3));
return (double)std::min(condition, VERDICT_DBL_MAX);
}
+double tri_condition(int num_nodes, const double coordinates[][3])
+{
+ return tri_condition_impl(num_nodes, coordinates, 3);
+}
+
+double tri_condition_from_loc_ptrs(int num_nodes, const double * const *
coordinates, const int dimension)
+{
+ return tri_condition_impl(num_nodes, coordinates, dimension);
+}
+
/*!
The scaled jacobian of a tri
minimum of the jacobian divided by the lengths of 2 edge vectors
*/
-double tri_scaled_jacobian(int /*num_nodes*/, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_scaled_jacobian_impl(int /*num_nodes*/, const CoordsContainerType
coordinates, const int dimension)
{
- VerdictVector first, second;
- double jacobian;
-
- VerdictVector edge[3];
- edge[0].set(coordinates[1][0] - coordinates[0][0], coordinates[1][1] -
coordinates[0][1],
- coordinates[1][2] - coordinates[0][2]);
-
- edge[1].set(coordinates[2][0] - coordinates[0][0], coordinates[2][1] -
coordinates[0][1],
- coordinates[2][2] - coordinates[0][2]);
+ const VerdictVector edge[3] =
+ {
+ {coordinates[0], coordinates[1], dimension},
+ {coordinates[0], coordinates[2], dimension},
+ {coordinates[1], coordinates[2], dimension},
+ };
- edge[2].set(coordinates[2][0] - coordinates[1][0], coordinates[2][1] -
coordinates[1][1],
- coordinates[2][2] - coordinates[1][2]);
- first = edge[1] - edge[0];
- second = edge[2] - edge[0];
+ const VerdictVector first = edge[1] - edge[0];
+ const VerdictVector second = edge[2] - edge[0];
- VerdictVector cross = first * second;
- jacobian = cross.length();
+ const VerdictVector cross = first * second;
+ double jacobian = cross.length();
- double max_edge_length_product;
- max_edge_length_product = std::max(edge[0].length() * edge[1].length(),
+ const double max_edge_length_product = std::max(edge[0].length() *
edge[1].length(),
std::max(edge[1].length() * edge[2].length(), edge[0].length() *
edge[2].length()));
if (max_edge_length_product < VERDICT_DBL_MIN)
@@ -528,6 +631,16 @@
return (double)std::max(jacobian, -VERDICT_DBL_MAX);
}
+double tri_scaled_jacobian(int num_nodes, const double coordinates[][3])
+{
+ return tri_scaled_jacobian_impl(num_nodes, coordinates, 3);
+}
+
+double tri_scaled_jacobian_from_loc_ptrs(int num_nodes, const double * const *
coordinates, const int dimension)
+{
+ return tri_scaled_jacobian_impl(num_nodes, coordinates, dimension);
+}
+
/*!
The shape of a tri
@@ -649,10 +762,11 @@
return (double)distortion;
}
- else if (num_nodes == 6)
+ else if (num_nodes >= 6)
{
total_number_of_gauss_points = 6;
number_of_gauss_points = 6;
+ num_nodes = 6; //seven nodes not handled
}
distortion = VERDICT_DBL_MAX;
@@ -829,7 +943,8 @@
return (double)std::max(distortion, -VERDICT_DBL_MAX);
}
-double tri_inradius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_inradius(const CoordsContainerType coordinates)
{
double sp = 0.0;
VerdictVector sides[3];
@@ -849,7 +964,8 @@
return ir;
}
-double tri6_min_inradius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri6_min_inradius(const CoordsContainerType coordinates, const int
dimension)
{
static int SUBTRI_NODES[4][3] = { { 0, 3, 5 }, { 3, 1, 4 }, { 5, 4, 2 }, {
3, 4, 5 } };
double min_inrad = VERDICT_DBL_MAX;
@@ -861,7 +977,7 @@
int idx = SUBTRI_NODES[i][j];
subtri_coords[j][0] = coordinates[idx][0];
subtri_coords[j][1] = coordinates[idx][1];
- subtri_coords[j][2] = coordinates[idx][2];
+ subtri_coords[j][2] = dimension == 2 ? 0.0 : coordinates[idx][2];
}
double subtri_inrad = tri_inradius(subtri_coords);
if (subtri_inrad < min_inrad)
@@ -872,7 +988,8 @@
return min_inrad;
}
-double calculate_tri3_outer_radius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double calculate_tri3_outer_radius(const CoordsContainerType coordinates,
const int dimension)
{
double sp = 0.0;
VerdictVector sides[3];
@@ -880,8 +997,8 @@
for (int i = 0; i < 3; i++)
{
int j = (i + 1) % 3;
- sides[i].set(coordinates[j][0] - coordinates[i][0], coordinates[j][1] -
coordinates[i][1],
- coordinates[j][2] - coordinates[i][2]);
+ const VerdictVector thisSide{coordinates[i], coordinates[j], dimension};
+ sides[i] = thisSide;
slen[i] = sides[i].length();
sp += slen[i];
}
@@ -895,50 +1012,80 @@
return cr;
}
-double tri6_normalized_inradius(const double coordinates[][3])
+static double fix_range(double v)
{
- double min_inradius_for_subtri = tri6_min_inradius(coordinates);
- double outer_radius = calculate_tri3_outer_radius(coordinates);
+ if (std::isnan(v))
+ {
+ return VERDICT_DBL_MAX;
+ }
+ if (v >= VERDICT_DBL_MAX)
+ {
+ return VERDICT_DBL_MAX;
+ }
+ if (v <= -VERDICT_DBL_MAX)
+ {
+ return -VERDICT_DBL_MAX;
+ }
+ return v;
+}
+
+template <typename CoordsContainerType>
+double tri6_normalized_inradius(const CoordsContainerType coordinates, const
int dimension)
+{
+ double min_inradius_for_subtri = tri6_min_inradius(coordinates, dimension);
+ double outer_radius = calculate_tri3_outer_radius(coordinates, dimension);
double normalized_inradius = 4.0 * min_inradius_for_subtri / outer_radius;
- return normalized_inradius;
+ return fix_range(normalized_inradius);
}
-double tri3_normalized_inradius(const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri3_normalized_inradius(const CoordsContainerType coordinates, const
int dimension)
{
- double tri6_coords[6][3];
+ std::array<const double*, 6> tri6_coords;
for (int i = 0; i < 3; i++)
- {
- tri6_coords[i][0] = coordinates[i][0];
- tri6_coords[i][1] = coordinates[i][1];
- tri6_coords[i][2] = coordinates[i][2];
- }
+ tri6_coords[i] = coordinates[i];
static int eidx[3][2] = { { 0, 1 }, { 1, 2 }, { 2, 0 } };
+ double tri6_midnode_coords[3][3];
for (int i = 3; i < 6; i++)
{
int i0 = eidx[i - 3][0];
int i1 = eidx[i - 3][1];
- tri6_coords[i][0] = (coordinates[i0][0] + coordinates[i1][0]) * 0.5;
- tri6_coords[i][1] = (coordinates[i0][1] + coordinates[i1][1]) * 0.5;
- tri6_coords[i][2] = (coordinates[i0][2] + coordinates[i1][2]) * 0.5;
+
+ tri6_midnode_coords[i-3][0] = (coordinates[i0][0] + coordinates[i1][0]) *
0.5;
+ tri6_midnode_coords[i-3][1] = (coordinates[i0][1] + coordinates[i1][1]) *
0.5;
+ tri6_midnode_coords[i-3][2] = (coordinates[i0][2] + coordinates[i1][2]) *
0.5;
+ tri6_coords[i] = tri6_midnode_coords[i-3];
}
- return tri6_normalized_inradius(tri6_coords);
+ return tri6_normalized_inradius(tri6_coords.data(), dimension);
}
//! Calculates the minimum normalized inner radius of a tri6
/** Ratio of the minimum subtri inner radius to tri outer radius*/
/* Currently, supports tri 6 and 3.*/
-double tri_normalized_inradius(int num_nodes, const double coordinates[][3])
+template <typename CoordsContainerType>
+double tri_normalized_inradius_impl(int num_nodes, const CoordsContainerType
coordinates, const int dimension)
{
if (num_nodes == 3)
{
- return tri3_normalized_inradius(coordinates);
+ return tri3_normalized_inradius(coordinates, dimension);
}
else if (num_nodes == 6)
{
- return tri6_normalized_inradius(coordinates);
+ return tri6_normalized_inradius(coordinates, dimension);
}
return 0.0;
}
+
+double tri_normalized_inradius(int num_nodes, const double coordinates[][3])
+{
+ return tri_normalized_inradius_impl(num_nodes, coordinates, 3);
+}
+
+double tri_normalized_inradius_from_loc_ptrs(int num_nodes, const double *
const *coordinates, const int dimension)
+{
+ return tri_normalized_inradius_impl(num_nodes, coordinates, dimension);
+}
+
} // namespace verdict
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/VerdictVector.hpp
new/verdict-1.4.2/VerdictVector.hpp
--- old/verdict-1.4.1/VerdictVector.hpp 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/VerdictVector.hpp 2023-07-19 18:21:57.000000000 +0200
@@ -45,9 +45,14 @@
//- Constructor: create vector from tuple
VerdictVector(const VerdictVector& tail, const VerdictVector& head);
+ VerdictVector(const double *tail, const double *head, int dimension);
+ VerdictVector(const double *tail, const double *head);
//- Constructor for a VerdictVector starting at tail and pointing
//- to head.
+ template <typename ARG1, typename ARG2, typename ARG3> VerdictVector(ARG1,
ARG2, ARG3) = delete;
+ //- define this template to avoid ambiguity between the (double, double,
double) and (double *, double *, int) constructors
+
VerdictVector(const VerdictVector& copy_from); //- Copy Constructor
//- Heading: Set and Inquire Functions
@@ -267,6 +272,20 @@
{
}
+inline VerdictVector::VerdictVector(const double *tail, const double *head,
int dimension)
+ : xVal{head[0] - tail[0]}
+ , yVal{head[1] - tail[1]}
+ , zVal{dimension == 2 ? 0.0 : head[2] - tail[2]}
+{
+}
+
+inline VerdictVector::VerdictVector(const double *tail, const double *head)
+ : xVal{head[0] - tail[0]}
+ , yVal{head[1] - tail[1]}
+ , zVal{head[2] - tail[2]}
+{
+}
+
inline VerdictVector::VerdictVector(const VerdictVector& tail, const
VerdictVector& head)
: xVal(head.xVal - tail.xVal)
, yVal(head.yVal - tail.yVal)
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/unittests/verdict.test.cpp
new/verdict-1.4.2/unittests/verdict.test.cpp
--- old/verdict-1.4.1/unittests/verdict.test.cpp 2022-10-04
18:52:36.000000000 +0200
+++ new/verdict-1.4.2/unittests/verdict.test.cpp 2023-07-19
18:21:57.000000000 +0200
@@ -265,6 +265,23 @@
runtest(testcase);
}
+TEST(verdict, edge_calc_4)
+{
+ test_case testcase = { "edge_calc_4", { { verdict::edge_length,
1.414213562373 } }, 3,
+ { { 0, 0, 0 }, { 1, 0, 0 }, { 0.5, 0.5, 0 } } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, edge_calc_5)
+{
+ test_case testcase = { "edge_calc_5", { { verdict::edge_length, 1.0 } }, 3,
+ { { 0, 0, 0 }, { 1, 0, 0 }, { 0.5, 0, 0 } } };
+
+ runtest(testcase);
+}
+
+
TEST(verdict, wedge_simple)
{
test_case testcase = { "wedge_simple",
@@ -324,11 +341,42 @@
runtest(testcase);
}
+TEST(verdict, tri4_area)
+{
+ test_case testcase = { "tri4_area",
+ { { verdict::tri_area, 0.8660254037844 }, },
+ 4,
+ {
+ { -0.5, 0.5, 0.5 },
+ { 0.5, -0.5, 0.5 },
+ { -0.5, -0.5, -0.5 },
+ { -0.166667, -0.166667, 0.166667 }
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, tri4_area_b)
+{
+ test_case testcase = { "tri4_area_b",
+ { { verdict::tri_area, 0.9006928267085 }, },
+ 4,
+ {
+ { -0.5, 0.5, 0.5 },
+ { 0.5, -0.5, 0.5 },
+ { -0.5, -0.5, -0.5 },
+ { -0.0666666666667, -0.0666666666667, 0.166667 }
+ } };
+
+ runtest(testcase);
+}
+
+
// tri_distortion is not well covered, we need to test it with a six noded
triangle */
TEST(verdict, tri_six_nodes)
{
test_case testcase = { "tri_six_nodes",
- { { verdict::tri_area, 0.54772255751 }, { verdict::tri_aspect_ratio,
1.3268079265 },
+ { { verdict::tri_area, 0.5554471166 }, { verdict::tri_aspect_ratio,
1.3268079265 },
{ verdict::tri_condition, 1.1173381066 }, { verdict::tri_distortion,
-3.3198288345 },
{ verdict::tri_minimum_angle, 45.406778838 }, {
verdict::tri_maximum_angle, 85.823122909 },
{ verdict::tri_shape, 0.89498424344 },
@@ -343,6 +391,42 @@
runtest(testcase);
}
+TEST(verdict, tri_seven_nodes)
+{
+ test_case testcase = { "tri_seven_nodes",
+ { { verdict::tri_area, 0.5477225575 }, { verdict::tri_aspect_ratio,
1.3268079265 },
+ { verdict::tri_condition, 1.117338106 }, { verdict::tri_distortion,
1.0000001069 },
+ { verdict::tri_minimum_angle, 45.406778838 }, {
verdict::tri_maximum_angle, 85.823122909 },
+ { verdict::tri_shape, 0.89498424344 },
+ /* 8 */ { verdict::tri_edge_ratio, 1.4005493428 },
+ /* 9 */ { verdict::tri_aspect_frobenius, 1.1173381066 },
+ /* 10 */ { verdict::tri_equiangle_skew, 0.24322035270 },
+ /* 11 */ { verdict::tri_normalized_inradius, 0.0 } },
+ 7,
+ { { 0, 0, 0 }, { 1, 0, 0.2 }, { 0, 1, 0.4 }, { 0.5, 0.0, 0.1 }, { 0.5, 0.5,
0.3 },
+ { 0.0, 0.5, 0.2 }, { 0.3333333333, 0.3333333333, 0.2 } } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, tri_seven_nodes_ho_nodes_moved)
+{
+ test_case testcase = { "tri_seven_nodes_ho_nodes_moved",
+ { { verdict::tri_area, 0.6215416837 }, { verdict::tri_aspect_ratio,
1.3268079265 },
+ { verdict::tri_condition, 1.117338106 }, { verdict::tri_distortion,
0.6514797121 },
+ { verdict::tri_minimum_angle, 45.406778838 }, {
verdict::tri_maximum_angle, 85.823122909 },
+ { verdict::tri_shape, 0.89498424344 },
+ /* 8 */ { verdict::tri_edge_ratio, 1.4005493428 },
+ /* 9 */ { verdict::tri_aspect_frobenius, 1.1173381066 },
+ /* 10 */ { verdict::tri_equiangle_skew, 0.24322035270 },
+ /* 11 */ { verdict::tri_normalized_inradius, 0.0 } },
+ 7,
+ { { 0, 0, 0 }, { 1, 0, 0.2 }, { 0, 1, 0.4 }, { 0.6, -0.1, 0.0 }, { 0.5, 0.6,
0.3 },
+ { 0.1, 0.5, 0.2 }, { 0.3333333333, 0.3333333333, 0.1 } } };
+
+ runtest(testcase);
+}
+
TEST(verdict, quad_simple1)
{
test_case testcase = { "quad_simple1",
@@ -492,6 +576,116 @@
runtest(testcase);
}
+TEST(verdict, quad5_area)
+{
+ test_case testcase = { "quad5_area",
+ { { verdict::quad_area, 21.5843025975 }, },
+ 5,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {3.125, 0.7, -1.0375 }
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, quad5_area_b)
+{
+ test_case testcase = { "quad5_area_b",
+ { { verdict::quad_area, 21.7405993317 }, },
+ 5,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {3.125, 0.7, -0.7875 }
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, quad8_area)
+{
+ test_case testcase = { "quad8_area",
+ { { verdict::quad_area, 21.5843485480 }, },
+ 8,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {1, -1.75, -1.075},
+ {4.5, -0.25, -1.125},
+ {5.25, 3.15, -1},
+ {1.75, 1.65, -0.95}
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, quad8_area_b)
+{
+ test_case testcase = { "quad8_area_b",
+ { { verdict::quad_area, 21.3830774682 }, },
+ 8,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {1, -1.75, -0.875},
+ {4.7, -0.25, -1.125},
+ {5.45, 3.35, -1 },
+ {2.15, 1.65, -0.95}
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, quad9_area)
+{
+ test_case testcase = { "quad9_area",
+ { { verdict::quad_area, 21.5843025975 }, },
+ 9,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {1, -1.75, -1.075},
+ {4.5, -0.25, -1.125},
+ {5.25, 3.15, -1 },
+ {1.75, 1.65, -0.95},
+ {3.125, 0.7, -1.0375 }
+ } };
+
+ runtest(testcase);
+}
+
+TEST(verdict, quad9_area_b)
+{
+ test_case testcase = { "quad9_area_b",
+ { { verdict::quad_area, 21.6059362945 }, },
+ 9,
+ {
+ { -1, -1, -1 },
+ { 3, -2.5, -1.15 },
+ { 6, 2, -1.1 },
+ { 4.5, 4.3, -0.9 },
+ {1, -1.75, -0.875},
+ {4.7, -0.25, -1.125},
+ {5.45, 3.35, -1 },
+ {2.15, 1.65, -0.95},
+ {3.125, 0.7, -1.2875 }
+ } };
+
+ runtest(testcase);
+}
+
TEST(verdict, tet_test1)
{
test_case testcase = { "tet_test1",
@@ -935,6 +1129,23 @@
runtest(testcase);
}
+TEST(verdict, tri_normalized_inradius_collapsed_tri6)
+{
+ test_case testcase = { "tri_normalized_inradius_good_tri6",
+ { { verdict::tri_normalized_inradius, verdict::VERDICT_DBL_MAX } }, 6,
+ {
+ {0.0, 0.0, 0.0},
+ {0.0, 0.0, 0.0},
+ {0.0, 0.0, 0.0},
+ {0.0, 0.0, 0.0},
+ {0.0, 0.0, 0.0},
+ {0.0, 0.0, 0.0},
+ } };
+
+ runtest(testcase);
+}
+
+
TEST(verdict, tet8_volume)
{
test_case testcase = { "Volume of tet8 element",
diff -urN '--exclude=CVS' '--exclude=.cvsignore' '--exclude=.svn'
'--exclude=.svnignore' old/verdict-1.4.1/verdict.h new/verdict-1.4.2/verdict.h
--- old/verdict-1.4.1/verdict.h 2022-10-04 18:52:36.000000000 +0200
+++ new/verdict-1.4.2/verdict.h 2023-07-19 18:21:57.000000000 +0200
@@ -283,6 +283,7 @@
length and the inradius of the tetrahedron
Reference --- P. Frey and P.-L. George, Meshing, Hermes (2000). */
VERDICT_EXPORT double tet_aspect_ratio(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tet_aspect_ratio_from_loc_ptrs(int num_nodes, const
double * const *coordinates);
//! Calculates tet aspect gamma metric.
/** Srms**3 / (8.479670*V) where Srms = sqrt(Sum(Si**2)/6), Si = edge length.
@@ -310,6 +311,7 @@
Reference --- V. N. Parthasarathy et al, A comparison of tetrahedron
quality measures, Finite Elem. Anal. Des., Vol 15(1993), 255-261. */
VERDICT_EXPORT double tet_volume(int num_nodes, const double coordinates[][3]);
+VERDICT_EXPORT double tet_volume_from_loc_ptrs(int num_nodes, const double *
const *coordinates);
//! Calculates tet condition metric.
/** Condition number of the Jacobian matrix at any corner.
@@ -317,6 +319,7 @@
Optimization of the Jacobian Matrix Norm and Associated Quantities,
Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */
VERDICT_EXPORT double tet_condition(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tet_condition_from_loc_ptrs(int num_nodes, const double
* const *coordinates);
//! Calculates tet jacobian.
/** Minimum pointwise volume at any corner.
@@ -331,6 +334,7 @@
Optimization of the Jacobian Matrix Norm and Associated Quantities,
Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */
VERDICT_EXPORT double tet_scaled_jacobian(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tet_scaled_jacobian_from_loc_ptrs(int num_nodes, const
double * const * coordinates);
//! Calculates tet mean ratio.
/** Ratio of tet volume to volume of an equilateral tet with the same RMS edge
length
@@ -338,11 +342,13 @@
meshes, IJNME 2010 82:843-867 Reference 2 --- Danial Ibanez - PhD Thesis,
Conformal Mesh
Adaptation on Heterogeneous Supercomputers */
VERDICT_EXPORT double tet_mean_ratio(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tet_mean_ratio_from_loc_ptrs(int num_nodes, const double
* const *coordinates);
//! Calculates the minimum normalized inner radius of a tet
/** Ratio of the minimum subtet inner radius to tet outer radius*/
/* Currently supports tetra 10 and 4.*/
VERDICT_EXPORT double tet_normalized_inradius(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tet_normalized_inradius_from_loc_ptrs(int num_nodes,
const double * const *coordinates);
//! Calculates tet shape metric.
/** 3/Mean Ratio of weighted Jacobian matrix.
@@ -615,6 +621,7 @@
Reference --- P. P. Pebay & T. J. Baker, Analysis of Triangle Quality
Measures, AMS Math. Comp., 2003, 72(244):1817-1839 */
VERDICT_EXPORT double tri_aspect_ratio(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tri_aspect_ratio_from_loc_ptrs(int num_nodes, const
double * const * coordinates, const int dimension = 3);
//! Calculates triangle metric.
/** radius ratio
@@ -629,6 +636,7 @@
//! Calculates triangle metric.
/** Maximum included angle in triangle */
VERDICT_EXPORT double tri_area(int num_nodes, const double coordinates[][3]);
+VERDICT_EXPORT double tri_area_from_loc_ptrs(int num_nodes, const double *
const *coordinates, const int dimension = 3);
//! Calculates triangle metric.
/** Minimum included angle in triangle */
@@ -644,6 +652,7 @@
Optimization of the Jacobian Matrix Norm and Associated Quantities,
Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */
VERDICT_EXPORT double tri_condition(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tri_condition_from_loc_ptrs(int num_nodes, const double
* const *coordinates, const int dimension = 3);
//! Calculates triangle metric.
/** Minimum Jacobian divided by the lengths of 2 edge vectors.
@@ -651,6 +660,7 @@
Optimization of the Jacobian Matrix Norm and Associated Quantities,
Intl. J. Numer. Meth. Engng. 2000, 48:1165-1185. */
VERDICT_EXPORT double tri_scaled_jacobian(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tri_scaled_jacobian_from_loc_ptrs(int num_nodes, const
double * const *coordinates, const int dimension=3);
//! Calculates triangle metric.
/** Min( J, 1/J ), where J is determinant of weighted Jacobian matrix.
@@ -684,6 +694,7 @@
/** Ratio of the minimum subtet inner radius to tet outer radius*/
/* Currently supports tri 6 and 3.*/
VERDICT_EXPORT double tri_normalized_inradius(int num_nodes, const double
coordinates[][3]);
+VERDICT_EXPORT double tri_normalized_inradius_from_loc_ptrs(int num_nodes,
const double * const *coordinates, const int dimension=3);
} // namespace verdict
#endif /* __verdict_h */
