From 5556fb39b141d67690c70a27c2c7cb0e8d2c30e9 Mon Sep 17 00:00:00 2001 From: Kristian Bendiksen Date: Wed, 22 Jan 2020 12:26:41 +0100 Subject: [PATCH] #4669 Add unit test for RigHexGradientTools. --- .../RigFemPartResultsCollection.cpp | 2 +- .../GeoMechDataModel/RigHexGradientTools.cpp | 14 +- .../UnitTests/CMakeLists_files.cmake | 1 + .../UnitTests/RigHexGradientTools-Test.cpp | 424 ++++++++++++++++++ 4 files changed, 433 insertions(+), 8 deletions(-) create mode 100644 ApplicationCode/UnitTests/RigHexGradientTools-Test.cpp diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp index dbca303d90..d33c6c0837 100644 --- a/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp +++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigFemPartResultsCollection.cpp @@ -1007,7 +1007,7 @@ RigFemScalarResultFrames* RigFemPartResultsCollection::calculateStressGradient( RigFemScalarResultFrames* inputResultFrames = this->findOrLoadScalarResult( partIndex, - RigFemResultAddress( resVarAddr.resultPosType, "SE", componentName.toStdString() ) ); + RigFemResultAddress( resVarAddr.resultPosType, "ST", componentName.toStdString() ) ); RigFemScalarResultFrames* dataFramesX = m_femPartResults[partIndex]->createScalarResult( RigFemResultAddress( resVarAddr.resultPosType, resVarAddr.fieldName, componentName.toStdString() + "-x" ) ); diff --git a/ApplicationCode/GeoMech/GeoMechDataModel/RigHexGradientTools.cpp b/ApplicationCode/GeoMech/GeoMechDataModel/RigHexGradientTools.cpp index 7471ed4e7b..00f76232ba 100644 --- a/ApplicationCode/GeoMech/GeoMechDataModel/RigHexGradientTools.cpp +++ b/ApplicationCode/GeoMech/GeoMechDataModel/RigHexGradientTools.cpp @@ -29,13 +29,13 @@ std::array RigHexGradientTools::gradients( const std::array gradientsUVW; gradientsUVW[0] = cvf::Vec3d( forwardFD( values, 0, 1 ), forwardFD( values, 0, 3 ), forwardFD( values, 0, 4 ) ); - gradientsUVW[1] = cvf::Vec3d( backwardFD( values, 1, 0 ), forwardFD( values, 1, 2 ), forwardFD( values, 1, 5 ) ); - gradientsUVW[2] = cvf::Vec3d( backwardFD( values, 2, 3 ), backwardFD( values, 2, 1 ), forwardFD( values, 2, 6 ) ); - gradientsUVW[3] = cvf::Vec3d( forwardFD( values, 3, 2 ), backwardFD( values, 3, 0 ), forwardFD( values, 3, 7 ) ); - gradientsUVW[4] = cvf::Vec3d( forwardFD( values, 4, 5 ), forwardFD( values, 4, 7 ), backwardFD( values, 4, 0 ) ); - gradientsUVW[5] = cvf::Vec3d( backwardFD( values, 5, 4 ), forwardFD( values, 5, 6 ), backwardFD( values, 5, 1 ) ); - gradientsUVW[6] = cvf::Vec3d( backwardFD( values, 6, 7 ), backwardFD( values, 6, 5 ), backwardFD( values, 6, 2 ) ); - gradientsUVW[7] = cvf::Vec3d( forwardFD( values, 7, 6 ), backwardFD( values, 7, 4 ), backwardFD( values, 7, 3 ) ); + gradientsUVW[1] = cvf::Vec3d( backwardFD( values, 0, 1 ), forwardFD( values, 1, 2 ), forwardFD( values, 1, 5 ) ); + gradientsUVW[2] = cvf::Vec3d( backwardFD( values, 3, 2 ), backwardFD( values, 1, 2 ), forwardFD( values, 2, 6 ) ); + gradientsUVW[3] = cvf::Vec3d( forwardFD( values, 3, 2 ), backwardFD( values, 0, 3 ), forwardFD( values, 3, 7 ) ); + gradientsUVW[4] = cvf::Vec3d( forwardFD( values, 4, 5 ), forwardFD( values, 4, 7 ), backwardFD( values, 0, 4 ) ); + gradientsUVW[5] = cvf::Vec3d( backwardFD( values, 4, 5 ), forwardFD( values, 5, 6 ), backwardFD( values, 1, 5 ) ); + gradientsUVW[6] = cvf::Vec3d( backwardFD( values, 7, 6 ), backwardFD( values, 5, 6 ), backwardFD( values, 2, 6 ) ); + gradientsUVW[7] = cvf::Vec3d( forwardFD( values, 7, 6 ), backwardFD( values, 4, 7 ), backwardFD( values, 3, 7 ) ); std::array NC; NC[0] = {-1, -1, -1}; diff --git a/ApplicationCode/UnitTests/CMakeLists_files.cmake b/ApplicationCode/UnitTests/CMakeLists_files.cmake index bbca23a011..56e8c2e527 100644 --- a/ApplicationCode/UnitTests/CMakeLists_files.cmake +++ b/ApplicationCode/UnitTests/CMakeLists_files.cmake @@ -63,6 +63,7 @@ ${CMAKE_CURRENT_LIST_DIR}/RiaSummaryCurveAnalyzer-Test.cpp ${CMAKE_CURRENT_LIST_DIR}/RiaStdStringTools-Test.cpp ${CMAKE_CURRENT_LIST_DIR}/RifWellMeasurementReader-Test.cpp ${CMAKE_CURRENT_LIST_DIR}/RiaDateStringParser-Test.cpp +${CMAKE_CURRENT_LIST_DIR}/RigHexGradientTools-Test.cpp ) if (RESINSIGHT_ENABLE_GRPC) diff --git a/ApplicationCode/UnitTests/RigHexGradientTools-Test.cpp b/ApplicationCode/UnitTests/RigHexGradientTools-Test.cpp new file mode 100644 index 0000000000..34e7f22347 --- /dev/null +++ b/ApplicationCode/UnitTests/RigHexGradientTools-Test.cpp @@ -0,0 +1,424 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2011-2012 Statoil ASA, Ceetron AS +// +// ResInsight is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. +// +// See the GNU General Public License at +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// + +#include "gtest/gtest.h" + +#include "RigHexGradientTools.h" + +#include "../cafHexInterpolator/cafHexInterpolator.h" + +#include "cvfMatrix4.h" + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForIdentityElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + std::array cornerValues = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + for ( int i = 0; i < 8; i++ ) + { + EXPECT_EQ( 0.0, gradients[i].x() ); + EXPECT_EQ( 0.0, gradients[i].y() ); + EXPECT_EQ( 0.0, gradients[i].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForScaledIdentityElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + for ( auto& v : hexCorners ) + v *= 2.5; + + std::array cornerValues = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + for ( int i = 0; i < 8; i++ ) + { + EXPECT_EQ( 0.0, gradients[i].x() ); + EXPECT_EQ( 0.0, gradients[i].y() ); + EXPECT_EQ( 0.0, gradients[i].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForTranslatedIdentityElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + for ( auto& v : hexCorners ) + v += {3.2, 9.5, -20.3}; + + std::array cornerValues = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + for ( int i = 0; i < 8; i++ ) + { + EXPECT_EQ( 0.0, gradients[i].x() ); + EXPECT_EQ( 0.0, gradients[i].y() ); + EXPECT_EQ( 0.0, gradients[i].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForRotatedIdentityElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + cvf::Mat4d rot = cvf::Mat4d::fromRotation( {1, 1, 0}, 3.0 ); + for ( auto& v : hexCorners ) + v.transformPoint( rot ); + + std::array cornerValues = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + for ( int i = 0; i < 8; i++ ) + { + EXPECT_EQ( 0.0, gradients[i].x() ); + EXPECT_EQ( 0.0, gradients[i].y() ); + EXPECT_EQ( 0.0, gradients[i].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + // Set a higher value in the first corner + std::array cornerValues = {2.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + // The adjacent corners should a non-zero gradient in the direction + // towards the first corner + EXPECT_EQ( -0.5, gradients[0].x() ); + EXPECT_EQ( -0.5, gradients[0].y() ); + EXPECT_EQ( -0.5, gradients[0].z() ); + + EXPECT_EQ( 0.5, gradients[1].x() ); + EXPECT_EQ( 0.0, gradients[1].y() ); + EXPECT_EQ( 0.0, gradients[1].z() ); + + EXPECT_EQ( 0.0, gradients[3].x() ); + EXPECT_EQ( 0.5, gradients[3].y() ); + EXPECT_EQ( 0.0, gradients[3].z() ); + + EXPECT_EQ( 0.0, gradients[4].x() ); + EXPECT_EQ( 0.0, gradients[4].y() ); + EXPECT_EQ( 0.5, gradients[4].z() ); + + // Non-adjacent should be unaffected + std::array independentCorners = {2, 5, 7, 6}; + for ( auto c : independentCorners ) + { + EXPECT_EQ( 0.0, gradients[c].x() ); + EXPECT_EQ( 0.0, gradients[c].y() ); + EXPECT_EQ( 0.0, gradients[c].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GradientsForLongElement ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 2.0, -1.0, 1.0 ), + cvf::Vec3d( 2.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + // Set a higher value in the first corner + std::array cornerValues = {2.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; + + std::array gradients = RigHexGradientTools::gradients( hexCorners, cornerValues ); + + // The adjacent corners should a non-zero gradient in the direction + // towards the first corner + EXPECT_EQ( -1.0 / 3.0, gradients[0].x() ); + EXPECT_EQ( -0.5, gradients[0].y() ); + EXPECT_EQ( -0.5, gradients[0].z() ); + + EXPECT_EQ( 1.0 / 3.0, gradients[1].x() ); + EXPECT_EQ( 0.0, gradients[1].y() ); + EXPECT_EQ( 0.0, gradients[1].z() ); + + EXPECT_EQ( 0.0, gradients[3].x() ); + EXPECT_EQ( 0.5, gradients[3].y() ); + EXPECT_EQ( 0.0, gradients[3].z() ); + + EXPECT_EQ( 0.0, gradients[4].x() ); + EXPECT_EQ( 0.0, gradients[4].y() ); + EXPECT_EQ( 0.5, gradients[4].z() ); + + // Non-adjacent should be unaffected + std::array independentCorners = {2, 5, 7, 6}; + for ( auto c : independentCorners ) + { + EXPECT_EQ( 0.0, gradients[c].x() ); + EXPECT_EQ( 0.0, gradients[c].y() ); + EXPECT_EQ( 0.0, gradients[c].z() ); + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GenerateJacobianForIdentity ) +{ + // Identity element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + for ( int i = 0; i < 8; i++ ) + { + cvf::Mat3d actualJacobian = caf::HexInterpolator::jacobi( hexCorners, hexCorners[i] ); + + // Calculate expected jacobian + double dx = 2.0; + double dy = 2.0; + double dz = 2.0; + + double di = 2.0; + double dj = 2.0; + double dk = 2.0; + + cvf::Mat3d expectedJacobian( dx / di, 0.0, 0.0, 0.0, dy / dj, 0.0, 0.0, 0.0, dz / dk ); + + cvf::Mat3d identityMat = cvf::Mat3d::IDENTITY; + + for ( int col = 0; col < 3; col++ ) + { + for ( int row = 0; row < 3; row++ ) + { + EXPECT_EQ( expectedJacobian.rowCol( row, col ), actualJacobian.rowCol( row, col ) ); + EXPECT_EQ( identityMat.rowCol( row, col ), actualJacobian.rowCol( row, col ) ); + } + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GenerateJacobianForIdentityScaled ) +{ + // Identity element + std::array normalizedCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, -1.0, -1.0 ), + cvf::Vec3d( 1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, -1.0, 1.0 ), + cvf::Vec3d( 1.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + double scale = 2.5; + + std::array hexCorners; + for ( int i = 0; i < 8; i++ ) + { + hexCorners[i] = normalizedCorners[i] * scale; + } + + for ( int i = 0; i < 8; i++ ) + { + cvf::Mat3d actualJacobian = caf::HexInterpolator::jacobi( hexCorners, normalizedCorners[i] ); + + // Calculate expected jacobian + double dx = 5.0; + double dy = 5.0; + double dz = 5.0; + + double di = 2.0; + double dj = 2.0; + double dk = 2.0; + + cvf::Mat3d expectedJacobian( dx / di, 0.0, 0.0, 0.0, dy / dj, 0.0, 0.0, 0.0, dz / dk ); + + cvf::Mat3d identityMat = cvf::Mat3d::IDENTITY; + + for ( int col = 0; col < 3; col++ ) + { + for ( int row = 0; row < 3; row++ ) + { + EXPECT_EQ( expectedJacobian.rowCol( row, col ), actualJacobian.rowCol( row, col ) ); + EXPECT_EQ( identityMat.rowCol( row, col ) * scale, actualJacobian.rowCol( row, col ) ); + } + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GenerateJacobianForLongElement ) +{ + // Try a more complex element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, 1.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 1.0 ), + cvf::Vec3d( 2.0, -1.0, 1.0 ), + cvf::Vec3d( 2.0, 1.0, 1.0 ), + cvf::Vec3d( -1.0, 1.0, 1.0 )}; + + cvf::Vec3d corner0( -1.0, -1.0, -1.0 ); + + cvf::Mat3d actualJacobian = caf::HexInterpolator::jacobi( hexCorners, corner0 ); + + // Calculate expected jacobian + double dx = 3.0; + double dy = 2.0; + double dz = 2.0; + + double di = 2.0; + double dj = 2.0; + double dk = 2.0; + + cvf::Mat3d expectedJacobian( dx / di, 0.0, 0.0, 0.0, dy / dj, 0.0, 0.0, 0.0, dz / dk ); + + for ( int col = 0; col < 3; col++ ) + { + for ( int row = 0; row < 3; row++ ) + { + EXPECT_EQ( expectedJacobian.rowCol( row, col ), actualJacobian.rowCol( row, col ) ); + } + } +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigHexGradientTools, GenerateJacobianForRotatedElement ) +{ + // Try a more complex element + std::array hexCorners = {cvf::Vec3d( -1.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, -1.0, -1.0 ), + cvf::Vec3d( 2.0, 3.0, -1.0 ), + cvf::Vec3d( -1.0, 3.0, -1.0 ), + cvf::Vec3d( -1.0, -1.0, 4.0 ), + cvf::Vec3d( 2.0, -1.0, 4.0 ), + cvf::Vec3d( 2.0, 3.0, 4.0 ), + cvf::Vec3d( -1.0, 3.0, 4.0 )}; + + cvf::Vec3d corner0( -1.0, -1.0, -1.0 ); + + cvf::Mat4d rot = cvf::Mat4d::fromRotation( {1, 0, 0}, cvf::PI_D / 2.0 ); + for ( auto& v : hexCorners ) + v.transformPoint( rot ); + + cvf::Mat3d actualJacobian = caf::HexInterpolator::jacobi( hexCorners, corner0 ); + + // Calculate expected jacobian + cvf::Vec3d dx_rot = ( hexCorners[1] - hexCorners[0] ); + cvf::Vec3d dy_rot = ( hexCorners[3] - hexCorners[0] ); + cvf::Vec3d dz_rot = ( hexCorners[4] - hexCorners[0] ); + + double di = 2.0; + double dj = 2.0; + double dk = 2.0; + + cvf::Mat3d expectedJacobian( dx_rot.x() / di, + dx_rot.y() / di, + dx_rot.z() / di, + dy_rot.x() / dj, + dy_rot.y() / dj, + dy_rot.z() / dj, + dz_rot.x() / dk, + dz_rot.y() / dk, + dz_rot.z() / dk ); + expectedJacobian.transpose(); + + for ( int row = 0; row < 3; row++ ) + { + for ( int col = 0; col < 3; col++ ) + { + EXPECT_EQ( expectedJacobian.rowCol( row, col ), actualJacobian.rowCol( row, col ) ); + } + } +}