Added a Tensor3-class to calculate principal values etc.

This commit is contained in:
Jacob Støren 2015-06-17 09:43:52 +02:00
parent 2cef1b8e3c
commit 30d41b9251
6 changed files with 800 additions and 0 deletions

View File

@ -0,0 +1,150 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions 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 <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "cvfBase.h"
#include "cafTensor3.h"
#include "cvfVector3.h"
#include "cvfMatrix3.h"
namespace caf {
//--------------------------------------------------------------------------------------------------
/// Compute the cofactors of the 3x3 matrix
/// Cofactor:
/// The determinant obtained by deleting the row and column of a given element and
/// preceded by a + or sign depending whether the element is in a + or position as follows:
///
/// + - +
/// - + -
/// + - +
//--------------------------------------------------------------------------------------------------
cvf::Mat3d cofactor3(const cvf::Mat3d& mx)
{
int detIdxi[2];
int detIdxj[2];
cvf::Mat3d cof;
double sign;
for (int i = 0; i < 3; i++)
{
detIdxi[0] = (i == 0) ? 1 : 0;
detIdxi[1] = (i == 2) ? 1 : 2;
for (int j = 0; j < 3; j++)
{
detIdxj[0] = (j == 0) ? 1 : 0;
detIdxj[1] = (j == 2) ? 1 : 2;
sign = (abs(j - i) == 1) ? -1 : 1;
cof(i, j) = sign * ( mx(detIdxi[0], detIdxj[0]) * mx(detIdxi[1], detIdxj[1])
- mx(detIdxi[0], detIdxj[1]) * mx(detIdxi[1], detIdxj[0]));
}
}
return cof;
}
//--------------------------------------------------------------------------------------------------
/// Compute the eigenvector of the matrix corresponding to the provided eigenvalue
/// The provided eigenvalue must be an actual eigenvalue of the matrix
//--------------------------------------------------------------------------------------------------
cvf::Vec3d eigenVector3(const cvf::Mat3d& mx, double eigenValue, bool* computedOk)
{
const double doubleThreshold = 1.0e-60;
if (computedOk) (*computedOk) = false;
cvf::Mat3d mxMinusEigv = mx;
for (int i = 0; i < 3; i++)
{
mxMinusEigv(i, i) -= eigenValue;
}
cvf::Mat3d cof = cofactor3(mxMinusEigv);
// Find largest absolute cofactor
int largestCof_i = -1;
int largestCof_j = -1;
double largestCof = 0.0;
for (int i = 0; i < 3; i++)
{
for (int j = 0; j < 3; j++)
{
double absCof = fabs(cof(i,j));
if (absCof > largestCof)
{
largestCof = absCof;
largestCof_i = i;
largestCof_j = j;
}
}
}
if (fabs(largestCof) < doubleThreshold) return cvf::Vec3d::ZERO;
// Find largest matrix element not in the max cofactor row/col
int largestMxElm_i = -1;
int largestMxElm_j = -1;
double largestMxElm = 0.0;
for (int i = 0; i < 3; i++)
{
if (i != largestCof_i)
{
for (int j = 0; j < 3; j++)
{
if (j != largestCof_j)
{
double absMxElm = fabs(mxMinusEigv(i,j));
if (absMxElm > largestMxElm)
{
largestMxElm = absMxElm;
largestMxElm_i = i;
largestMxElm_j = j;
}
}
}
}
}
// Check if largest coefficient is zero
if (fabs(largestMxElm) < doubleThreshold) return cvf::Vec3d::ZERO;
// Find last component index
int lastComp_j = 0;
for (int i = 0; i < 3; i++)
{
if ((i != largestCof_j) && (i != largestMxElm_j)) lastComp_j = i;
}
cvf::Vec3d eigenVector;
eigenVector[largestCof_j] = 1.0;
eigenVector[lastComp_j] = cof(largestCof_i, lastComp_j) / cof(largestCof_i, largestCof_j);
eigenVector[largestMxElm_j] = (-mxMinusEigv(largestMxElm_i, largestCof_j) - mxMinusEigv(largestMxElm_i, lastComp_j)*eigenVector[lastComp_j] )
/ mxMinusEigv(largestMxElm_i, largestMxElm_j);
if (computedOk) (*computedOk) = true;
return eigenVector;
}
}

View File

@ -0,0 +1,58 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions 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 <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once
#include "cvfBase.h"
#include "cvfVector3.h"
namespace caf
{
template< typename S>
class Tensor3
{
S m_tensor[6]; // SXX, SYY, SZZ, SXY, SYZ, SZX
public:
Tensor3() {}
Tensor3(S sxx, S syy, S szz, S sxy, S syz, S szx);
Tensor3(const Tensor3& other);
inline Tensor3& operator=(const Tensor3& rhs);
bool equals(const Tensor3& mat) const;
bool operator==(const Tensor3& rhs) const;
bool operator!=(const Tensor3& rhs) const;
enum TensorComponentEnum { SXX, SYY, SZZ, SXY, SYZ, SZX };
inline S& operator[](TensorComponentEnum comp);
inline S operator[](TensorComponentEnum comp) const;
void setFromInternalLayout(S* tensorData);
void setFromAbaqusLayout(S* tensorData);
cvf::Vec3f calculatePrincipals(cvf::Vec3f principalDirections[3]);
float calculateVonMises();
};
typedef Tensor3<float> Ten3f;
}
#include "cafTensor3.inl"

View File

@ -0,0 +1,324 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions 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 <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "cvfAssert.h"
#include "cvfMath.h"
#include "cvfSystem.h"
#include <algorithm>
#include "cvfMatrix3.h"
namespace caf {
//----------------------------------------------------------------------------------------------------
/// Copy constructor
//----------------------------------------------------------------------------------------------------
template <typename S>
inline Tensor3<S>::Tensor3(const Tensor3& other)
{
cvf::System::memcpy(m_tensor, sizeof(m_tensor), other.m_tensor, sizeof(other.m_tensor));
}
//----------------------------------------------------------------------------------------------------
/// Constructor with explicit initialization of all tensor elements.
///
//----------------------------------------------------------------------------------------------------
template <typename S>
Tensor3<S>::Tensor3(S sxx, S syy, S szz, S sxy, S syz, S szx)
{
m_tensor[0] = sxx;
m_tensor[1] = syy;
m_tensor[2] = szz;
m_tensor[3] = sxy;
m_tensor[4] = syz;
m_tensor[5] = szx;
}
//----------------------------------------------------------------------------------------------------
/// Assignment operator
//----------------------------------------------------------------------------------------------------
template <typename S>
inline Tensor3<S>& Tensor3<S>::operator=(const Tensor3& obj)
{
cvf::System::memcpy(m_tensor, sizeof(m_tensor), obj.m_tensor, sizeof(obj.m_tensor));
return *this;
}
//----------------------------------------------------------------------------------------------------
/// Check if matrices are equal using exact comparisons.
//----------------------------------------------------------------------------------------------------
template<typename S>
bool Tensor3<S>::equals(const Tensor3& ten) const
{
for (int i = 0; i < 6; i++)
{
if (m_tensor[i] != ten.m_tensor[i]) return false;
}
return true;
}
//----------------------------------------------------------------------------------------------------
/// Comparison operator. Checks for equality using exact comparisons.
//----------------------------------------------------------------------------------------------------
template <typename S>
bool Tensor3<S>::operator==(const Tensor3& rhs) const
{
return this->equals(rhs);
}
//----------------------------------------------------------------------------------------------------
/// Comparison operator. Checks for not equal using exact comparisons.
//----------------------------------------------------------------------------------------------------
template <typename S>
bool Tensor3<S>::operator!=(const Tensor3& rhs) const
{
int i;
for (i = 0; i < 6; i++)
{
if (m_tensor[i] != rhs.m_tensor[i]) return true;
}
return false;
}
//--------------------------------------------------------------------------------------------------
/// Get modifiable component 0,1,2. E.g. x = v[0];
//--------------------------------------------------------------------------------------------------
template<typename S>
inline S Tensor3<S>::operator[](TensorComponentEnum index) const
{
CVF_TIGHT_ASSERT(index >= 0);
CVF_TIGHT_ASSERT(index < 6);
return m_tensor[index];
}
//--------------------------------------------------------------------------------------------------
/// Get const component 0,1,2. E.g. x = v[0];
//--------------------------------------------------------------------------------------------------
template<typename S>
inline S& Tensor3<S>::operator[](TensorComponentEnum index)
{
CVF_TIGHT_ASSERT(index >= 0);
CVF_TIGHT_ASSERT(index < 6);
return m_tensor[index];
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
template< typename S>
void Tensor3<S>::setFromInternalLayout(S* tensorData)
{
m_tensor[0] = tensorData[0];
m_tensor[1] = tensorData[1];
m_tensor[2] = tensorData[2];
m_tensor[3] = tensorData[3];
m_tensor[4] = tensorData[4];
m_tensor[5] = tensorData[5];
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
template< typename S>
void Tensor3<S>::setFromAbaqusLayout(S* tensorData)
{
m_tensor[0] = tensorData[0];
m_tensor[1] = tensorData[1];
m_tensor[2] = tensorData[2];
m_tensor[3] = tensorData[3];
m_tensor[4] = tensorData[5];
m_tensor[5] = tensorData[4];
}
cvf::Mat3d cofactor3(const cvf::Mat3d& mx);
cvf::Vec3d eigenVector3(const cvf::Mat3d& mx, double eigenValue, bool* computedOk);
//--------------------------------------------------------------------------------------------------
/// Compute principal values and optionally the principal directions
/// The tensor must be laid out as follows: SXX, SYY, SZZ, SXY, SYZ, SZX
//--------------------------------------------------------------------------------------------------
template<typename S>
cvf::Vec3f Tensor3<S>::calculatePrincipals( cvf::Vec3f principalDirections[3])
{
CVF_TIGHT_ASSERT(m_tensor);
const float floatThreshold = 1.0e-30f;
const double doubleThreshold = 1.0e-60;
cvf::Vec3f principalValues;
// Init return arrays to invalid
principalValues[0] = std::numeric_limits<float>::infinity();
principalValues[1] = std::numeric_limits<float>::infinity();
principalValues[2] = std::numeric_limits<float>::infinity();
if (principalDirections)
{
principalDirections[0] = cvf::Vec3f::ZERO;
principalDirections[1] = cvf::Vec3f::ZERO;
principalDirections[2] = cvf::Vec3f::ZERO;
}
// Return if we have an undefined component
int i;
for (i = 0; i < 6; i++)
{
if (m_tensor[i] == std::numeric_limits<S>::infinity())
{
return principalValues;
}
}
// Return 0, 0, 0 if all components are zero
bool isAllTensCompsZero = true;
for (i = 0; i < 6; i++)
{
if (!(abs(m_tensor[i]) < floatThreshold))
{
isAllTensCompsZero = false;
break;
}
}
if (isAllTensCompsZero)
{
return cvf::Vec3f::ZERO;
}
double SXX = m_tensor[0], SYY = m_tensor[1], SZZ = m_tensor[2];
double SXY = m_tensor[3], SYZ = m_tensor[4], SZX = m_tensor[5];
double pressure = -(SXX + SYY + SZZ)/3.0;
// Normally we would solve the eigenvalues by solving the 3'rd degree equation:
// -sigma^3 + A*sigma^2 - B*sigma + C = 0
// in which A, B, and C are the invariants of the stress tensor.
// http://www.engapplets.vt.edu/Mohr/java/nsfapplets/MohrCircles2-3D/Theory/theory.htm
// But the roots(eigenvalues) are calculated by transforming the above equation into
// s**3 + aa*s + b = 0 and using the trignometric solution.
// See crc standard mathematical tables 19th edition pp. 103-104.
SXX += pressure;
SYY += pressure;
SZZ += pressure;
double S1, S2, S3;
double AA, BB, CC, DD, angleP;
AA = SXY*SXY + SYZ*SYZ + SZX*SZX - SXX*SYY - SYY*SZZ - SXX*SZZ;
BB = SXX * SYZ * SYZ
+ SYY * SZX * SZX
+ SZZ * SXY * SXY
- SXX * SYY * SZZ
- 2.0 * SXY * SYZ * SZX;
if (fabs(AA) < doubleThreshold)
{
S1 = 0.0;
S2 = 0.0;
S3 = 0.0;
}
else
{
CC = -sqrt(27.0/AA) * BB * 0.5 / AA;
if (CC > 1.0) CC = 1.0;
else if (CC < -1.0) CC = -1.0;
angleP = acos(CC)/3.0;
DD = 2.0*sqrt(AA/3.0);
S1 = DD*cos(angleP);
S2 = DD*cos(angleP + 4.0*cvf::PI_D/3.0);
S3 = DD*cos(angleP + 2.0*cvf::PI_D/3.0);
}
int idxPMin = 2;
int idxPMid = 1;
int idxPMax = 0;
double principalsd[3];
principalsd[idxPMax] = (S1 - pressure);
principalsd[idxPMid] = (S2 - pressure);
principalsd[idxPMin] = (S3 - pressure);
// Sort the principals if we have no Z component in the tensor at all
if ((m_tensor[2] == 0.0f) && (m_tensor[4] == 0.0f) && (m_tensor[5] == 0.0f))
{
if (fabs(principalsd[idxPMin]) > fabs(principalsd[idxPMid])) std::swap(idxPMin, idxPMid);
if (fabs(principalsd[idxPMin]) > fabs(principalsd[idxPMax])) std::swap(idxPMin, idxPMax);
if (principalsd[idxPMax] < principalsd[idxPMid]) std::swap(idxPMax, idxPMid);
principalsd[idxPMin] = 0;
}
// Calculate the principal directions if needed
if (principalDirections)
{
cvf::Mat3d T;
T(0,0) = m_tensor[0]; T(0,1) = m_tensor[3]; T(0,2) = m_tensor[5];
T(1,0) = m_tensor[3]; T(1,1) = m_tensor[1]; T(1,2) = m_tensor[4];
T(2,0) = m_tensor[5]; T(2,1) = m_tensor[4]; T(2,2) = m_tensor[2];
principalDirections[0] = cvf::Vec3f(eigenVector3(T, principalsd[idxPMax], NULL));
principalDirections[0].normalize();
principalDirections[1] = cvf::Vec3f(eigenVector3(T, principalsd[idxPMid], NULL));
principalDirections[1].normalize();
principalDirections[2] = cvf::Vec3f(eigenVector3(T, principalsd[idxPMin], NULL));
principalDirections[2].normalize();
}
principalValues[0] = (float)principalsd[idxPMax];
principalValues[1] = (float)principalsd[idxPMid];
principalValues[2] = (float)principalsd[idxPMin];
return principalValues;
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
template< typename S>
float caf::Tensor3<S>::calculateVonMises()
{
return (float) sqrt( ( (m_tensor[0]*m_tensor[0] + m_tensor[1]*m_tensor[1] + m_tensor[2]*m_tensor[2]) ) +
( -(m_tensor[0]*m_tensor[1] + m_tensor[1]*m_tensor[2] + m_tensor[0]*m_tensor[2]) ) +
( 3*(m_tensor[3]*m_tensor[3] + m_tensor[4]*m_tensor[4] + m_tensor[5]*m_tensor[5]) ) );
}
}

View File

@ -0,0 +1,32 @@
cmake_minimum_required (VERSION 2.8)
project ( cafTensor_UnitTests )
set(RI_VIZ_FWK_ROOT ../../../Fwk/VizFwk CACHE PATH "Path to VizFwk")
set(RI_GTEST_ROOT .. CACHE PATH "Path to folder containing gtest folder")
set(RI_SRC_ROOT ../../cafTensor CACHE PATH "Path to the code to test")
set(RI_TEST_FILE "" CACHE FILEPATH "Path to test file")
include(${RI_VIZ_FWK_ROOT}/CMake/Utils/ceeDetermineCompilerFlags.cmake)
add_subdirectory(${RI_VIZ_FWK_ROOT}/LibCore buildVizFwk)
add_definitions( -DTEST_FILE="${RI_TEST_FILE}")
include_directories(${RI_VIZ_FWK_ROOT}/LibCore)
include_directories(${RI_GTEST_ROOT})
include_directories(${RI_SRC_ROOT})
set( UNIT_TEST_CPP_SOURCES
main.cpp
cafTensor_UnitTests.cpp
${RI_SRC_ROOT}/cafTensor3.cpp
${RI_SRC_ROOT}/cafTensor3.h
${RI_SRC_ROOT}/cafTensor3.inl
${RI_GTEST_ROOT}/gtest/gtest-all.cpp
)
add_executable( ${PROJECT_NAME} ${UNIT_TEST_CPP_SOURCES} )
target_link_libraries( ${PROJECT_NAME} LibCore)

View File

@ -0,0 +1,194 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015- Statoil ASA
// Copyright (C) 2015- Ceetron Solutions 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 <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "cafTensor3.h"
#include "gtest/gtest.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, BasicTests)
{
caf::Ten3f T1;
caf::Ten3f T2(1, 2, 3, 4, 5, 6);
caf::Ten3f T3(T2);
EXPECT_EQ(1, T2[caf::Ten3f::SXX]);
EXPECT_EQ(2, T2[caf::Ten3f::SYY]);
EXPECT_EQ(3, T2[caf::Ten3f::SZZ]);
EXPECT_EQ(4, T2[caf::Ten3f::SXY]);
EXPECT_EQ(5, T2[caf::Ten3f::SYZ]);
EXPECT_EQ(6, T2[caf::Ten3f::SZX]);
T1 = T2;
EXPECT_EQ(1, T1[caf::Ten3f::SXX]);
EXPECT_EQ(2, T1[caf::Ten3f::SYY]);
EXPECT_EQ(3, T1[caf::Ten3f::SZZ]);
EXPECT_EQ(4, T1[caf::Ten3f::SXY]);
EXPECT_EQ(5, T1[caf::Ten3f::SYZ]);
EXPECT_EQ(6, T1[caf::Ten3f::SZX]);
EXPECT_TRUE(T2 == T3);
EXPECT_TRUE(T1 == T2);
EXPECT_TRUE(T1.equals(T2));
EXPECT_FALSE(T1 != T2);
T1[caf::Ten3f::SXX] = 7;
EXPECT_TRUE(T1 != T2);
EXPECT_FALSE(T1 == T2);
EXPECT_FALSE(T1.equals(T2));
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, setFromNativeArray)
{
float tensData[6] = {11,12,13,14,15,16};
caf::Ten3f T1;
T1.setFromAbaqusLayout(tensData);
EXPECT_EQ(11, T1[caf::Ten3f::SXX]);
EXPECT_EQ(12, T1[caf::Ten3f::SYY]);
EXPECT_EQ(13, T1[caf::Ten3f::SZZ]);
EXPECT_EQ(14, T1[caf::Ten3f::SXY]);
EXPECT_EQ(15, T1[caf::Ten3f::SZX]);
EXPECT_EQ(16, T1[caf::Ten3f::SYZ]);
caf::Ten3f T2;
T2.setFromInternalLayout(tensData);
EXPECT_EQ(11, T2[caf::Ten3f::SXX]);
EXPECT_EQ(12, T2[caf::Ten3f::SYY]);
EXPECT_EQ(13, T2[caf::Ten3f::SZZ]);
EXPECT_EQ(14, T2[caf::Ten3f::SXY]);
EXPECT_EQ(15, T2[caf::Ten3f::SYZ]);
EXPECT_EQ(16, T2[caf::Ten3f::SZX]);
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, zero)
{
caf::Ten3f T0(0,0,0,0,0,0);
cvf::Vec3f pDirs[3];
cvf::Vec3f p0 = T0.calculatePrincipals(pDirs);
EXPECT_TRUE(p0 == cvf::Vec3f::ZERO);
EXPECT_TRUE(pDirs[0] == cvf::Vec3f::ZERO);
EXPECT_TRUE(pDirs[1] == cvf::Vec3f::ZERO);
EXPECT_TRUE(pDirs[2] == cvf::Vec3f::ZERO);
float vm = T0.calculateVonMises();
EXPECT_EQ(0.0f, vm );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, undef)
{
float inf = std::numeric_limits<float>::infinity();
caf::Ten3f T0(0,0,0,0,0,inf);
cvf::Vec3f pDirs[3];
cvf::Vec3f p0 = T0.calculatePrincipals(pDirs);
EXPECT_TRUE(p0 == cvf::Vec3f(inf, inf, inf));
EXPECT_TRUE(pDirs[0] == cvf::Vec3f::ZERO);
EXPECT_TRUE(pDirs[1] == cvf::Vec3f::ZERO);
EXPECT_TRUE(pDirs[2] == cvf::Vec3f::ZERO);
float vm = T0.calculateVonMises();
EXPECT_EQ(inf, vm );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, realTensors1)
{
caf::Ten3f T0(80,50,20,40,45,50);
cvf::Vec3f pDirs[3];
cvf::Vec3f p0 = T0.calculatePrincipals(pDirs);
EXPECT_NEAR( 143.8f, p0[0], 0.1 );
EXPECT_NEAR( 23.6f, p0[1], 0.1 );
EXPECT_NEAR( -17.4f, p0[2], 0.1 );
float vm = T0.calculateVonMises();
EXPECT_NEAR(145.2f, vm, 0.1 );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, realTensors2)
{
caf::Ten3f T0(20,50,80,50,45,40);
cvf::Vec3f pDirs[3];
cvf::Vec3f p0 = T0.calculatePrincipals(pDirs);
EXPECT_NEAR( 143.8f, p0[0], 0.1 );
EXPECT_NEAR( 23.9f, p0[1], 0.1 );
EXPECT_NEAR( -17.6f, p0[2], 0.1 );
float vm = T0.calculateVonMises();
EXPECT_NEAR(145.2f, vm, 0.1 );
}
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
TEST(cafTensor3Test, realTensors3)
{
caf::Ten3f T0(10,20,30,0,0,0);
cvf::Vec3f pDirs[3];
cvf::Vec3f p0 = T0.calculatePrincipals(pDirs);
EXPECT_NEAR( 30.0f, p0[0], 0.1 );
EXPECT_NEAR( 20.0f, p0[1], 0.1 );
EXPECT_NEAR( 10.0f, p0[2], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[0][0], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[0][1], 0.1 );
EXPECT_NEAR( 1.0f, pDirs[0][2], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[1][0], 0.1);
EXPECT_NEAR( 1.0f, pDirs[1][1], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[1][2], 0.1 );
EXPECT_NEAR( 1.0f, pDirs[2][0], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[2][1], 0.1 );
EXPECT_NEAR( 0.0f, pDirs[2][2], 0.1 );
float vm = T0.calculateVonMises();
EXPECT_NEAR(17.3f, vm, 0.1 );
}

View File

@ -0,0 +1,42 @@
/////////////////////////////////////////////////////////////////////////////////
//
// Copyright (C) 2015 - Statoil ASA
// Copyright (C) 2015 - Ceetron Solutions 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 <http://www.gnu.org/licenses/gpl.html>
// for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#include "cvfBase.h"
#include "gtest/gtest.h"
#include <stdio.h>
#include "cvfTrace.h"
//--------------------------------------------------------------------------------------------------
///
//--------------------------------------------------------------------------------------------------
int main(int argc, char **argv)
{
cvf::Assert::setReportMode(cvf::Assert::CONSOLE);
testing::InitGoogleTest(&argc, argv);
int result = RUN_ALL_TESTS();
std::cout << "Please press <Enter> to close the window.";
std::cin.get();
return result;
}