diff --git a/applications/test/DiagTensor/Make/files b/applications/test/DiagTensor/Make/files new file mode 100644 index 0000000000..a3e36220b3 --- /dev/null +++ b/applications/test/DiagTensor/Make/files @@ -0,0 +1,3 @@ +Test-DiagTensor.C + +EXE = $(FOAM_USER_APPBIN)/Test-DiagTensor diff --git a/applications/test/tensor/Make/options b/applications/test/DiagTensor/Make/options similarity index 100% rename from applications/test/tensor/Make/options rename to applications/test/DiagTensor/Make/options diff --git a/applications/test/DiagTensor/Test-DiagTensor.C b/applications/test/DiagTensor/Test-DiagTensor.C new file mode 100644 index 0000000000..48b68be74c --- /dev/null +++ b/applications/test/DiagTensor/Test-DiagTensor.C @@ -0,0 +1,439 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-DiagTensor + +Description + Tests for \c DiagTensor constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "Tensor.H" +#include "SymmTensor.H" +#include "SphericalTensor.H" +#include "DiagTensor.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of DiagTensor, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const DiagTensor dT(Zero); + Info<< dT << endl; + } + { + Info<< "# Construct given VectorSpace of the same rank:" << nl; + const VectorSpace, Type, 3> V(Zero); + const DiagTensor dT(V); + Info<< dT << endl; + } + { + Info<< "# Construct given the three components:" << nl; + const DiagTensor dT + ( + Type(1), + Type(5), + Type(-9) + ); + Info<< dT << endl; + } + { + Info<< "# Copy construct:" << nl; + const DiagTensor dT(Zero); + const DiagTensor copydT(dT); + Info<< dT << tab << copydT << endl; + } +} + + +// Execute each member function of DiagTensor, and print output +template +void test_member_funcs(Type) +{ + DiagTensor dT(Type(1), Type(5), Type(-9)); + const DiagTensor cdT(Type(-9), Type(5), Type(1)); + + Info<< "# Operand: " << nl + << " DiagTensor = " << dT << endl; + + + { + Info<< "# Component access:" << nl; + + DiagTensor cpdT(dT.xx(), dT.yy(), dT.zz()); + cmp(" 'DiagTensor' access:", dT, cpdT); + + const DiagTensor cpcdT(cdT.xx(), cdT.yy(), cdT.zz()); + cmp(" 'const DiagTensor' access:", cdT, cpcdT); + } +} + + +// Execute each global function of DiagTensor, and print output +template +void test_global_funcs(Type) +{ + const Tensor T + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(-1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const DiagTensor dT(Type(1), Type(5), Type(-9)); + + Info<< "# Operands: " << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << nl + << " DiagTensor = " << dT << endl; + + + cmp(" Trace = ", tr(dT), Type(-3)); + cmp(" Spherical part = ", sph(dT), SphericalTensor(tr(dT)/Type(3))); + cmp(" Determinant = ", det(dT), Type(-44.99999999999999)); + cmp + ( + " Inverse = ", + inv(dT), + DiagTensor(Type(1), Type(0.2), Type(-0.11111111)) + ); + cmp + ( + " Diagonal of Tensor = ", + diag(T), + DiagTensor(Type(-1), Type(5), Type(-9)) + ); + cmp + ( + " Diagonal of SymmTensor = ", + diag(sT), + DiagTensor(Type(-1), Type(5), Type(-9)) + ); +} + + +// Execute each global operator of DiagTensor, and print output +template +void test_global_opers(Type) +{ + const Tensor T + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(-1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const DiagTensor dT(Type(1), Type(5), Type(-9)); + const SphericalTensor spT(Type(1)); + const Vector v(Type(3), Type(2), Type(1)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << nl + << " DiagTensor = " << dT << nl + << " SphericalTensor = " << spT << nl + << " Vector = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Sum of DiagTensor-Tensor = ", + (dT + T), + Tensor + ( + Type(0), Type(2), Type(-3), + Type(4), Type(10), Type(-6), + Type(7), Type(8), Type(-18) + ) + ); + cmp + ( + " Sum of Tensor-DiagTensor = ", + (T + dT), + Tensor + ( + Type(0), Type(2), Type(-3), + Type(4), Type(10), Type(-6), + Type(7), Type(8), Type(-18) + ) + ); + cmp + ( + " Subtract Tensor from DiagTensor = ", + (dT - T), + Tensor + ( + Type(2), Type(-2), Type(3), + Type(-4), Type(0), Type(6), + Type(-7), Type(-8), Type(0) + ) + ); + cmp + ( + " Subtract DiagTensor from Tensor = ", + (T - dT), + Tensor + ( + Type(-2), Type(2), Type(-3), + Type(4), Type(0), Type(-6), + Type(7), Type(8), Type(0) + ) + ); + cmp + ( + " Division of Type by DiagTensor = ", + (x/dT), + DiagTensor(Type(4), Type(0.8), Type(-0.44444444)) + ); + cmp + ( + " Division of DiagTensor by Type = ", + (dT/x), + DiagTensor(Type(0.25), Type(1.25), Type(-2.25)) + ); + cmp + ( + " Division of Vector by DiagTensor = ", + (v/dT), + Vector(Type(3), Type(0.4), Type(-0.11111111)) + ); + cmp + ( + " Inner-product of DiagTensor-DiagTensor = ", + (dT & dT), + DiagTensor(Type(1), Type(25), Type(81)) + ); + cmp + ( + " Inner-product of DiagTensor-Tensor = ", + (dT & T), + Tensor + ( + Type(-1), Type(2), Type(-3), + Type(20), Type(25), Type(-30), + Type(-63), Type(-72), Type(81) + ) + ); + cmp + ( + " Inner-product of Tensor-DiagTensor = ", + (T & dT), + Tensor + ( + Type(-1), Type(10), Type(27), + Type(4), Type(25), Type(54), + Type(7), Type(40), Type(81) + ) + ); + cmp + ( + " Inner-product of DiagTensor-Vector = ", + (dT & v), + Vector(Type(3), Type(10), Type(-9)) + ); + cmp + ( + " Inner-product of Vector-DiagTensor = ", + (v & dT), + Vector(Type(3), Type(10), Type(-9)) + ); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main() +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "DiagTensor", + "DiagTensor", + "DiagTensor" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/test/SphericalTensor/Make/files b/applications/test/SphericalTensor/Make/files new file mode 100644 index 0000000000..0a2c967665 --- /dev/null +++ b/applications/test/SphericalTensor/Make/files @@ -0,0 +1,3 @@ +Test-SphericalTensor.C + +EXE = $(FOAM_USER_APPBIN)/Test-SphericalTensor diff --git a/applications/test/SphericalTensor/Make/options b/applications/test/SphericalTensor/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/SphericalTensor/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/SphericalTensor/Test-SphericalTensor.C b/applications/test/SphericalTensor/Test-SphericalTensor.C new file mode 100644 index 0000000000..79f2b9d738 --- /dev/null +++ b/applications/test/SphericalTensor/Test-SphericalTensor.C @@ -0,0 +1,348 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-SphericalTensor + +Description + Tests for \c SphericalTensor constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "Tensor.H" +#include "SymmTensor.H" +#include "SphericalTensor.H" +#include "DiagTensor.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of SphericalTensor, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const SphericalTensor spT(Zero); + Info<< spT << endl; + } + { + Info<< "# Construct given VectorSpace of the same rank:" << nl; + const VectorSpace, Type, 1> V(Zero); + const SphericalTensor spT(V); + Info<< spT << endl; + } + { + Info<< "# Construct given the component:" << nl; + const SphericalTensor spT(Type(1)); + Info<< spT << endl; + } + { + Info<< "# Copy construct:" << nl; + const SphericalTensor spT(Zero); + const SphericalTensor copyspT(spT); + Info<< spT << tab << copyspT << endl; + } +} + + +// Execute each member function of SphericalTensor, and print output +template +void test_member_funcs(Type) +{ + SphericalTensor spT(Type(1)); + const SphericalTensor cspT(Type(-9)); + + Info<< "# Operand: " << nl + << " SphericalTensor = " << spT << endl; + + + { + Info<< "# Component access:" << nl; + + SphericalTensor cpspT(spT.ii()); + cmp(" 'SphericalTensor' access:", spT, cpspT); + + const SphericalTensor cpcspT(cspT.ii()); + cmp(" 'const SphericalTensor' access:", cspT, cpcspT); + } + { + Info<< "# SphericalTensor operations:" << nl; + + Info<< " Transpose:" << nl; + cmp(" 'SphericalTensor'.T():", spT.T(), spT); + } +} + + +// Execute each global function of SphericalTensor, and print output +template +void test_global_funcs(Type) +{ + const SphericalTensor spT(Type(5)); + + Info<< "# Operand: " << nl + << " SphericalTensor = " << spT << endl; + + + cmp(" Trace = ", tr(spT), Type(15)); + cmp(" Spherical part = ", sph(spT), spT); + cmp(" Determinant = ", det(spT), Type(124.99999999999994)); + cmp + ( + " Inverse = ", + inv(spT), + SphericalTensor(Type(0.2)) + ); + cmp(" Square of Frobenius norm = ", magSqr(spT), Type(75)); + cmp(" Max component = ", cmptMax(spT), Type(5)); + cmp(" Min component = ", cmptMax(spT), Type(5)); + cmp(" Sum of components = ", cmptSum(spT), Type(15)); + cmp(" Arithmetic average of components = ", cmptAv(spT), Type(5)); +} + + +// Execute each global operator of SphericalTensor, and print output +template +void test_global_opers(Type) +{ + const Tensor T + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(-1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const DiagTensor dT(Type(1), Type(5), Type(-9)); + const SphericalTensor spT(Type(-2)); + const Vector v(Type(3), Type(2), Type(1)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << nl + << " DiagTensor = " << dT << nl + << " SphericalTensor = " << spT << nl + << " Vector = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Division of Type by SpTensor = ", + (x/spT), + SphericalTensor(Type(-2)) + ); + cmp + ( + " Division of SpTensor by Type = ", + (spT/x), + SphericalTensor(Type(-0.5)) + ); + cmp + ( + " Inner-product of SpTensor-SpTensor = ", + (spT & spT), + SphericalTensor(Type(4)) + ); + cmp + ( + " Inner-product of SpTensor-Vector = ", + (spT & v), + Vector(Type(-6), Type(-4), Type(-2)) // Column-vector + ); + cmp + ( + " Inner-product of Vector-SpTensor = ", + (v & spT), + Vector(Type(-6), Type(-4), Type(-2)) // Row-vector + ); + cmp(" D-inner-product of SpTensor-SpTensor = ", (spT && spT), Type(12)); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main() +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "SphericalTensor", + "SphericalTensor", + "SphericalTensor" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/test/SphericalTensor2D/Make/files b/applications/test/SphericalTensor2D/Make/files new file mode 100644 index 0000000000..2b492a2f8c --- /dev/null +++ b/applications/test/SphericalTensor2D/Make/files @@ -0,0 +1,3 @@ +Test-SphericalTensor2D.C + +EXE = $(FOAM_USER_APPBIN)/Test-SphericalTensor2D diff --git a/applications/test/SphericalTensor2D/Make/options b/applications/test/SphericalTensor2D/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/SphericalTensor2D/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/SphericalTensor2D/Test-SphericalTensor2D.C b/applications/test/SphericalTensor2D/Test-SphericalTensor2D.C new file mode 100644 index 0000000000..00d17610dc --- /dev/null +++ b/applications/test/SphericalTensor2D/Test-SphericalTensor2D.C @@ -0,0 +1,331 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-SphericalTensor2D + +Description + Tests for \c SphericalTensor2D constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "Tensor2D.H" +#include "SymmTensor2D.H" +#include "SphericalTensor2D.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of SphericalTensor2D, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const SphericalTensor2D spT(Zero); + Info<< spT << endl; + } + { + Info<< "# Construct given VectorSpace of the same rank:" << nl; + const VectorSpace, Type, 1> V(Zero); + const SphericalTensor2D spT(V); + Info<< spT << endl; + } + { + Info<< "# Construct given the component:" << nl; + const SphericalTensor2D spT(Type(1)); + Info<< spT << endl; + } + { + Info<< "# Copy construct:" << nl; + const SphericalTensor2D spT(Zero); + const SphericalTensor2D copyspT(spT); + Info<< spT << tab << copyspT << endl; + } +} + + +// Execute each member function of SphericalTensor2D, and print output +template +void test_member_funcs(Type) +{ + SphericalTensor2D spT(Type(1)); + const SphericalTensor2D cspT(Type(-9)); + + Info<< "# Operand: " << nl + << " SphericalTensor2D = " << spT << endl; + + + { + Info<< "# Component access:" << nl; + + SphericalTensor2D cpspT(spT.ii()); + cmp(" 'SphericalTensor2D' access:", spT, cpspT); + + const SphericalTensor2D cpcspT(cspT.ii()); + cmp(" 'const SphericalTensor2D' access:", cspT, cpcspT); + } +} + + +// Execute each global function of SphericalTensor2D, and print output +template +void test_global_funcs(Type) +{ + const SphericalTensor2D spT(Type(5)); + + Info<< "# Operand: " << nl + << " SphericalTensor2D = " << spT << endl; + + + cmp(" Trace = ", tr(spT), Type(10)); + cmp(" Spherical part = ", sph(spT), spT); + cmp(" Determinant = ", det(spT), Type(24.99999999999994)); + cmp + ( + " Inverse = ", + inv(spT), + SphericalTensor2D(Type(0.2)) + ); +} + + +// Execute each global operator of SphericalTensor2D, and print output +template +void test_global_opers(Type) +{ + const Tensor2D T + ( + Type(-1), Type(2), + Type(4), Type(5) + ); + const SymmTensor2D sT + ( + Type(-1), Type(2), + Type(5) + ); + const SphericalTensor2D spT(Type(-2)); + const Vector2D v(Type(3), Type(2)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor2D = " << T << nl + << " SymmTensor2D = " << sT << nl + << " SphericalTensor2D = " << spT << nl + << " Vector2D = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Division of Type by SpTensor2D = ", + (x/spT), + SphericalTensor2D(Type(-2)) + ); + cmp + ( + " Division of SpTensor2D by Type = ", + (spT/x), + SphericalTensor2D(Type(-0.5)) + ); + cmp + ( + " Inner-product of SpTensor2D-SpTensor2D = ", + (spT & spT), + SphericalTensor2D(Type(4)) + ); + cmp + ( + " Inner-product of SpTensor2D-Vector2D = ", + (spT & v), + Vector2D(Type(-6), Type(-4)) // Column-vector + ); + cmp + ( + " Inner-product of Vector2D-SpTensor2D = ", + (v & spT), + Vector2D(Type(-6), Type(-4)) // Row-vector + ); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main() +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "SphericalTensor2D", + "SphericalTensor2D", + "SphericalTensor2D" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/test/SymmTensor/Make/files b/applications/test/SymmTensor/Make/files new file mode 100644 index 0000000000..186a997815 --- /dev/null +++ b/applications/test/SymmTensor/Make/files @@ -0,0 +1,3 @@ +Test-SymmTensor.C + +EXE = $(FOAM_USER_APPBIN)/Test-SymmTensor diff --git a/applications/test/SymmTensor/Make/options b/applications/test/SymmTensor/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/SymmTensor/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/SymmTensor/Test-SymmTensor.C b/applications/test/SymmTensor/Test-SymmTensor.C new file mode 100644 index 0000000000..5e48ef7c34 --- /dev/null +++ b/applications/test/SymmTensor/Test-SymmTensor.C @@ -0,0 +1,576 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-SymmTensor + +Description + Tests for \c SymmTensor constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Eigen decomposition tests for \c symmTensor, i.e. SymmTensor. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "symmTensor.H" +#include "transform.H" +#include "Random.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Create a random symmTensor +symmTensor makeRandomContainer(Random& rnd) +{ + symmTensor A(Zero); + std::generate(A.begin(), A.end(), [&]{ return rnd.GaussNormal(); }); + return A; +} + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of SymmTensor, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const SymmTensor sT(Zero); + Info<< sT << endl; + } + { + Info<< "# Construct given VectorSpace of the same rank:" << nl; + const VectorSpace, Type, 6> M(Zero); + const SymmTensor sT(M); + Info<< sT << endl; + } + { + Info<< "# Construct given SphericalTensor:" << nl; + const SphericalTensor Sp(Type(5)); + const SymmTensor sT(Sp); + Info<< sT << endl; + } + { + Info<< "# Construct given the six components:" << nl; + const SymmTensor sT + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + Info<< sT << endl; + } + { + Info<< "# Copy construct:" << nl; + const SymmTensor sT(Zero); + const SymmTensor copysT(sT); + Info<< sT << tab << copysT << endl; + } +} + + +// Execute each member function of SymmTensor, and print output +template +void test_member_funcs(Type) +{ + SymmTensor sT + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const SymmTensor csT + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + + Info<< "# Operand: " << nl + << " SymmTensor = " << sT << endl; + + + { + Info<< "# Component access:" << nl; + + SymmTensor cpsT + ( + sT.xx(), sT.xy(), sT.xz(), + sT.yy(), sT.yz(), + sT.zz() + ); + cmp(" 'SymmTensor' access:", sT, cpsT); + cmp(" xy()=yx():", sT.xy(), sT.yx()); + cmp(" xz()=zx():", sT.xz(), sT.zx()); + cmp(" yz()=zy():", sT.yz(), sT.zy()); + + const SymmTensor cpcsT + ( + csT.xx(), csT.xy(), csT.xz(), + csT.yy(), csT.yz(), + csT.zz() + ); + cmp(" 'const SymmTensor' access:", csT, cpcsT); + cmp(" xy()=yx():", sT.xy(), sT.yx()); + cmp(" xz()=zx():", sT.xz(), sT.zx()); + cmp(" yz()=zy():", sT.yz(), sT.zy()); + } + { + Info<< "# Diagonal access:" << nl; + + cmp + ( + " 'SymmTensor'.diag():", + sT.diag(), + Vector(Type(1), Type(5), Type(-9)) + ); + cmp + ( + " 'const SymmTensor'.diag():", + csT.diag(), + Vector(Type(1), Type(5), Type(-9)) + ); + + + Info<< "# Diagonal manipulation:" << nl; + + sT.diag(Vector(Type(-10), Type(-15), Type(-20))); + cmp + ( + " 'SymmTensor'.diag('Vector'):", + sT.diag(), + Vector(Type(-10), Type(-15), Type(-20)) + ); + } + { + Info<< "# Tensor operations:" << nl; + + Info<< " Transpose:" << nl; + cmp(" 'SymmTensor'.T():", sT.T(), sT); + } + { + Info<< "# Member operators:" << nl; + + sT = SphericalTensor(Type(5)); + cmp + ( + " Assign to a SphericalTensor:", + sT, + SymmTensor + ( + Type(5), Zero, Zero, + Type(5), Zero, + Type(5) + ) + ); + } +} + + +// Execute each global function of SymmTensor, and print output +template +void test_global_funcs(Type) +{ + const SymmTensor sT + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + + Info<< "# Operand: " << nl + << " SymmTensor = " << sT << nl << endl; + + + cmp(" Trace = ", tr(sT), Type(-3)); + cmp(" Spherical part = ", sph(sT), SphericalTensor(tr(sT)/Type(3))); + cmp(" Symmetric part = ", symm(sT), sT); + cmp(" Twice the symmetric part = ", twoSymm(sT), 2*sT); + cmp + ( + " Deviatoric part = ", + dev(sT), + SymmTensor + ( + Type(2), Type(2), Type(-3), + Type(6), Type(-6), + Type(-8) + ) + ); + cmp(" Two-third deviatoric part = ", dev2(sT), sT - 2*sph(sT)); + cmp(" Determinant = ", det(sT), Type(-17.999999999999996)); + cmp + ( + " Cofactor tensor = ", + cof(sT), + SymmTensor + ( + Type(-81), Type(36), Type(3), + Type(-18), Type(0), + Type(1) + ) + ); + cmp + ( + " Inverse = ", + inv(sT, det(sT)), + SymmTensor + ( + Type(4.5), Type(-2), Type(-0.16666667), + Type(1), Type(0), + Type(-0.05555556) + ), + 1e-8, + 1e-8 + ); + cmp + ( + " Inverse (another) = ", + inv(sT), + SymmTensor + ( + Type(4.5), Type(-2), Type(-0.16666667), + Type(1), Type(0), + Type(-0.05555556) + ), + 1e-8, + 1e-8 + ); + cmp(" First invariant = ", invariantI(sT), Type(-3)); + cmp(" Second invariant = ", invariantII(sT), Type(-98)); + cmp(" Third invariant = ", invariantIII(sT), Type(-17.999999999999996)); + cmp + ( + " Inner-product with self = ", + innerSqr(sT), + SymmTensor + ( + Type(14), Type(30), Type(12), + Type(65), Type(18), + Type(126) + ) + ); + cmp(" Square of Frobenius norm = ", magSqr(sT), Type(205)); +} + + +// Execute each global operator of SymmTensor, and print output +template +void test_global_opers(Type) +{ + const Tensor T + ( + Type(1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const SphericalTensor spT(Type(1)); + const Vector v(Type(3), Type(2), Type(1)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << nl + << " SphericalTensor = " << spT << nl + << " Vector = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Sum of SpTensor-SymmTensor = ", + (spT + sT), + SymmTensor + ( + Type(2), Type(2), Type(-3), + Type(6), Type(-6), + Type(-8) + ) + ); + cmp + ( + " Sum of SymmTensor-SpTensor = ", + (sT + spT), + SymmTensor + ( + Type(2), Type(2), Type(-3), + Type(6), Type(-6), + Type(-8) + ) + ); + cmp + ( + " Subtract SymmTensor from SpTensor = ", + (spT - sT), + SymmTensor + ( + Type(0), Type(-2), Type(3), + Type(-4), Type(6), + Type(10) + ) + ); + cmp + ( + " Subtract SpTensor from SymmTensor = ", + (sT - spT), + SymmTensor + ( + Type(0), Type(2), Type(-3), + Type(4), Type(-6), + Type(-10) + ) + ); + cmp + ( + " Hodge dual of a SymmTensor", + *sT, + Vector(Type(-6), Type(3), Type(2)) + ); + cmp + ( + " Division of a SymmTensor by a Type", + sT/x, + SymmTensor + ( + Type(0.25), Type(0.5), Type(-0.75), + Type(1.25), Type(-1.5), + Type(-2.25) + ) + ); + cmp + ( + " Inner-product of SymmTensor-SymmTensor = ", + (sT & sT), + Tensor + ( + Type(14), Type(30), Type(12), + Type(30), Type(65), Type(18), + Type(12), Type(18), Type(126) + ) + ); + cmp + ( + " Inner-product of SpTensor-SymmTensor = ", + (spT & sT), + SymmTensor + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ) + ); + cmp + ( + " Inner-product of SymmTensor-SpTensor = ", + (sT & spT), + SymmTensor + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ) + ); + cmp + ( + " Inner-product of SymmTensor-Vector = ", + (sT & v), + Vector(Type(4), Type(10), Type(-30)) // Column-vector + ); + cmp + ( + " Inner-product of Vector-SymmTensor = ", + (v & sT), + Vector(Type(4), Type(10), Type(-30)) // Row-vector + ); + cmp(" D-inner-product of SymmTensor-SymmTensor = ", (sT && sT), Type(205)); + cmp(" D-inner-product of SymmTensor-SpTensor = ", (sT && spT), Type(-3)); + cmp(" D-inner-product of SpTensor-SymmTensor = ", (spT && sT), Type(-3)); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main() +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "SymmTensor", + "SymmTensor", + "SymmTensor" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; + + return 0; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/test/SymmTensor2D/Make/files b/applications/test/SymmTensor2D/Make/files new file mode 100644 index 0000000000..ce2bf56e3f --- /dev/null +++ b/applications/test/SymmTensor2D/Make/files @@ -0,0 +1,3 @@ +Test-SymmTensor2D.C + +EXE = $(FOAM_USER_APPBIN)/Test-SymmTensor2D diff --git a/applications/test/SymmTensor2D/Make/options b/applications/test/SymmTensor2D/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/SymmTensor2D/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/SymmTensor2D/Test-SymmTensor2D.C b/applications/test/SymmTensor2D/Test-SymmTensor2D.C new file mode 100644 index 0000000000..236516fb04 --- /dev/null +++ b/applications/test/SymmTensor2D/Test-SymmTensor2D.C @@ -0,0 +1,538 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2019-2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-SymmTensor2D + +Description + Tests for \c SymmTensor2D constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Eigen decomposition tests for \c symmTensor2D, i.e. SymmTensor2D. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "symmTensor2D.H" +#include "transform.H" +#include "Random.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Create a random symmTensor2D +symmTensor2D makeRandomContainer(Random& rnd) +{ + symmTensor2D A(Zero); + std::generate(A.begin(), A.end(), [&]{ return rnd.GaussNormal(); }); + return A; +} + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of SymmTensor2D, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const SymmTensor2D sT(Zero); + Info<< sT << endl; + } + { + Info<< "# Construct given VectorSpace:" << nl; + const VectorSpace, Type, 3> V(Zero); + const SymmTensor2D sT(V); + Info<< sT << endl; + } + { + Info<< "# Construct given SphericalTensor2D:" << nl; + const SphericalTensor2D Sp(Type(5)); + const SymmTensor2D sT(Sp); + Info<< sT << endl; + } + { + Info<< "# Construct given the three components:" << nl; + const SymmTensor2D sT + ( + Type(1), Type(2), + Type(4) + ); + Info<< sT << endl; + } + { + Info<< "# Copy construct:" << nl; + const SymmTensor2D S(Type(1), Type(2), Type(3)); + const SymmTensor2D sT(S); + Info<< sT << endl; + } +} + + +// Execute each member function of SymmTensor2D, and print output +template +void test_member_funcs(Type) +{ + SymmTensor2D sT(Type(1), Type(2), Type(-3)); + const SymmTensor2D csT(Type(-3), Type(2), Type(1)); + + + { + Info<< "# Component access:" << nl; + + SymmTensor2D cpsT + ( + sT.xx(), sT.xy(), + sT.yy() + ); + cmp(" 'SymmTensor2D' access:", sT, cpsT); + cmp(" xy()=yx():", sT.xy(), sT.yx()); + + const SymmTensor2D cpcsT + ( + csT.xx(), csT.xy(), + csT.yy() + ); + cmp(" 'const SymmTensor2D' access:", csT, cpcsT); + cmp(" xy()=yx():", sT.xy(), sT.yx()); + } + { + Info<< "# Diagonal access:" << nl; + + cmp + ( + " 'SymmTensor2D'.diag():", + sT.diag(), + Vector2D(Type(1), Type(-3)) + ); + cmp + ( + " 'const SymmTensor2D'.diag():", + csT.diag(), + Vector2D(Type(-3), Type(1)) + ); + + + Info<< "# Diagonal manipulation:" << nl; + + sT.diag(Vector2D(Type(-10), Type(-15))); + cmp + ( + " 'SymmTensor2D'.diag('Vector2D'):", + sT.diag(), + Vector2D(Type(-10), Type(-15)) + ); + } + { + Info<< "# Tensor operations:" << nl; + + Info<< " Transpose:" << nl; + cmp(" 'SymmTensor2D'.T():", sT.T(), sT); + } + { + Info<< "# Member operators:" << nl; + + sT = SphericalTensor2D(Type(5)); + cmp + ( + " Assign to a SphericalTensor2D:", + sT, + SymmTensor2D + ( + Type(5), Zero, + Type(5) + ) + ); + } +} + + +// Execute each global function of SymmTensor2D, and print output +template +void test_global_funcs(Type) +{ + const SymmTensor2D sT(Type(1), Type(2), Type(-3)); + const Vector2D v(Type(-3), Type(2)); + + Info<< "# Operands: " << nl + << " SymmTensor2D = " << sT << nl + << " Vector2D = " << v << endl; + + + cmp(" Trace = ", tr(sT), Type(-2)); + cmp(" Spherical part = ", sph(sT), SphericalTensor2D(tr(sT)/Type(2))); + cmp(" Symmetric part = ", symm(sT), sT); + cmp(" Twice the symmetric part = ", twoSymm(sT), 2*sT); + cmp + ( + " Deviatoric part = ", + dev(sT), + SymmTensor2D + ( + Type(2), Type(2), + Type(-2) + ) + ); + cmp(" Two-third deviatoric part = ", dev2(sT), sT - 2*sph(sT)); + cmp(" Determinant = ", det(sT), Type(-7.000000000000001)); + cmp + ( + " Cofactor tensor = ", + cof(sT), + SymmTensor2D + ( + Type(-3), Type(-2), + Type(1) + ) + ); + cmp + ( + " Inverse = ", + inv(sT, det(sT)), + SymmTensor2D + ( + Type(0.42857143), Type(0.28571429), + Type(-0.14285714) + ), + 1e-8, + 1e-6 + ); + cmp + ( + " Inverse (another) = ", + inv(sT), + SymmTensor2D + ( + Type(0.42857143), Type(0.28571429), + Type(-0.14285714) + ), + 1e-8, + 1e-6 + ); + cmp(" First invariant = ", invariantI(sT), Type(-2)); + cmp(" Second invariant = ", invariantII(sT), Type(-7)); + cmp + ( + " Inner-product with self = ", + innerSqr(sT), + SymmTensor2D + ( + Type(5), Type(-4), + Type(13) + ) + ); + cmp(" Square of Frobenius norm = ", magSqr(sT), Type(17.999999999999996)); + cmp + ( + " Outer-product of a Vector2D with itself = ", + sqr(v), + SymmTensor2D + ( + Type(9), Type(-6), + Type(4) + ) + ); +} + + +// Execute each global operator of SymmTensor2D, and print output +template +void test_global_opers(Type) +{ + const Tensor2D T + ( + Type(1), Type(-2), + Type(3), Type(-4) + ); + const SymmTensor2D sT + ( + Type(1), Type(2), + Type(-4) + ); + const SphericalTensor2D spT(Type(-1)); + const Vector2D v(Type(3), Type(-2)); + const Type x(-4); + + Info<< "# Operands:" << nl + << " Tensor2D = " << T << nl + << " SymmTensor2D = " << sT << nl + << " SphericalTensor2D = " << spT << nl + << " Vector2D = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Sum of SpTensor2D-SymmTensor2D = ", + (spT + sT), + SymmTensor2D + ( + Type(0), Type(2), + Type(-5) + ) + ); + cmp + ( + " Sum of SymmTensor2D-SpTensor2D = ", + (sT + spT), + SymmTensor2D + ( + Type(0), Type(2), + Type(-5) + ) + ); + cmp + ( + " Subtract SymmTensor2D from SpTensor2D = ", + (spT - sT), + SymmTensor2D + ( + Type(-2), Type(-2), + Type(3) + ) + ); + cmp + ( + " Subtract SpTensor2D from SymmTensor2D = ", + (sT - spT), + SymmTensor2D + ( + Type(2), Type(2), + Type(-3) + ) + ); + cmp + ( + " Division of a SymmTensor2D by a Type", + sT/x, + SymmTensor2D + ( + Type(-0.25), Type(-0.5), + Type(1) + ) + ); + cmp + ( + " Inner-product of SymmTensor2D-SymmTensor2D = ", + (sT & sT), + Tensor2D + ( + Type(5), Type(-6), + Type(-6), Type(20) + ) + ); + cmp + ( + " Inner-product of SpTensor2D-SymmTensor2D = ", + (spT & sT), + SymmTensor2D + ( + Type(-1), Type(-2), + Type(4) + ) + ); + cmp + ( + " Inner-product of SymmTensor2D-SpTensor2D = ", + (sT & spT), + SymmTensor2D + ( + Type(-1), Type(-2), + Type(4) + ) + ); + cmp + ( + " Inner-product of SymmTensor2D-Vector2D = ", + (sT & v), + Vector2D(Type(-1), Type(14)) // Column-vector + ); + cmp + ( + " Inner-product of Vector2D-SymmTensor2D = ", + (v & sT), + Vector2D(Type(-1), Type(14)) // Row-vector + ); + cmp(" D-inner-prod of SymmTensor2D-SymmTensor2D = ", (sT && sT), Type(25)); + cmp(" D-inner-prod of SymmTensor2D-SpTensor2D = ", (sT && spT), Type(3)); + cmp(" D-inner-prod of SpTensor2D-SymmTensor2D = ", (spT && sT), Type(3)); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "SymmTensor2D", + "SymmTensor2D", + "SymmTensor2D" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// ************************************************************************* // diff --git a/applications/test/Tensor/Make/files b/applications/test/Tensor/Make/files new file mode 100644 index 0000000000..c0464c2461 --- /dev/null +++ b/applications/test/Tensor/Make/files @@ -0,0 +1,3 @@ +Test-Tensor.C + +EXE = $(FOAM_USER_APPBIN)/Test-Tensor diff --git a/applications/test/Tensor/Make/options b/applications/test/Tensor/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/Tensor/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/Tensor/Test-Tensor.C b/applications/test/Tensor/Test-Tensor.C new file mode 100644 index 0000000000..9bc1066726 --- /dev/null +++ b/applications/test/Tensor/Test-Tensor.C @@ -0,0 +1,964 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2018 OpenFOAM Foundation + Copyright (C) 2019-2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-Tensor + +Description + Tests for \c Tensor constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Eigen decomposition tests for \c tensor, i.e. Tensor. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "tensor.H" +#include "transform.H" +#include "Random.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Create a random tensor +tensor makeRandomContainer(Random& rnd) +{ + tensor A(Zero); + std::generate(A.begin(), A.end(), [&]{ return rnd.GaussNormal(); }); + return A; +} + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of Tensor, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const Tensor T(Zero); + Info<< T << endl; + } + { + Info<< "# Construct given MatrixSpace of the same rank:" << nl; + const MatrixSpace, Type, 3, 3> M(Zero); + const Tensor T(M); + Info<< T << endl; + } + { + Info<< "# Construct given VectorSpace of the same rank:" << nl; + const VectorSpace, Type, 9> V(Zero); + const Tensor T(V); + Info<< T << endl; + } + { + Info<< "# Construct given SphericalTensor:" << nl; + const SphericalTensor Sp(Type(5)); + const Tensor T(Sp); + Info<< T << endl; + } + { + Info<< "# Construct given SymmTensor:" << nl; + const SymmTensor S + ( + Type(1), Type(2), Type(3), + Type(5), Type(6), + Type(9) + ); + const Tensor T(S); + Info<< T << endl; + } + { + Info<< "# Construct given triad of row vectors," + << " optionally treated as transposed (ie, column vectors)" << nl; + const Vector> vecs + ( + Vector(Type(1), Type(2), Type(3)), + Vector(Type(4), Type(5), Type(6)), + Vector(Type(7), Type(8), Type(9)) + ); + const Tensor T(vecs); + Info<< T << nl; + + const Tensor transposedT(vecs, true); + Info<< transposedT << endl; + } + { + Info<< "# Construct given the three row vectors," + << " optionally treated as transposed (ie, column vectors)" << nl; + const Vector a(Type(1), Type(2), Type(3)); + const Vector b(Type(4), Type(5), Type(6)); + const Vector c(Type(7), Type(8), Type(9)); + const Tensor T(a, b, c); + Info<< T << nl; + + const Tensor transposedT(a, b, c, true); + Info<< transposedT << endl; + } + { + Info<< "# Construct given the nine components:" << nl; + const Tensor T + ( + Type(1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + Info<< T << endl; + } + { + Info<< "# Copy construct:" << nl; + const Tensor T(Zero); + const Tensor copyT(T); + Info<< T << tab << copyT << endl; + } +} + + +// Execute each member function of Tensor, and print output +template +void test_member_funcs(Type) +{ + Tensor T + ( + Type(1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + Tensor Tbak = T; + const Tensor cT + ( + Type(-9), Type(8), Type(7), + Type(-6), Type(5), Type(4), + Type(-3), Type(2), Type(1) + ); + + Info<< "# Operand: " << nl + << " Tensor = " << T << endl; + + + { + Info<< "# Component access:" << nl; + + Tensor cpT + ( + T.xx(), T.xy(), T.xz(), + T.yx(), T.yy(), T.yz(), + T.zx(), T.zy(), T.zz() + ); + cmp(" 'Tensor' access:", T, cpT); + + const Tensor cpcT + ( + cT.xx(), cT.xy(), cT.xz(), + cT.yx(), cT.yy(), cT.yz(), + cT.zx(), cT.zy(), cT.zz() + ); + cmp(" 'const Tensor' access:", cT, cpcT); + } + { + Info<< "# Column-vector access:" << nl; + cmp(" cx():", T.cx(), Vector(Type(1), Type(4), Type(7))); + cmp(" cy():", T.cy(), Vector(Type(2), Type(5), Type(8))); + cmp(" cz():", T.cz(), Vector(Type(-3), Type(-6), Type(-9))); + cmp(" col(0):", T.col(0), Vector(Type(1), Type(4), Type(7))); + cmp(" col(1):", T.col(1), Vector(Type(2), Type(5), Type(8))); + cmp(" col(2):", T.col(2), Vector(Type(-3), Type(-6), Type(-9))); + cmp + ( + " col<0>:", + T.template col<0>(), + Vector(Type(1), Type(4), Type(7)) + ); + cmp + ( + " col<1>:", + T.template col<1>(), + Vector(Type(2), Type(5), Type(8)) + ); + cmp + ( + " col<2>:", + T.template col<2>(), + Vector(Type(-3), Type(-6), Type(-9)) + ); + // Compilation error: Info << " col<3> = " << T.col<3>() << nl; + + + Info<< "# Column-vector manipulation:" << nl; + T.col(1, Vector(Type(0), Type(1), Type(99))); + cmp + ( + " col(1, Vector):", + T.col(1), + Vector(Type(0), Type(1), Type(99)) + ); + + T.cols + ( + Vector(Type(1), Type(1), Type(1)), + Vector(Type(-1), Type(1), Type(2)), + Vector(Type(1), Type(1), Type(3)) + ); + cmp + ( + " cols(Vectors):", + T, + Tensor + ( + Type(1), Type(-1), Type(1), + Type(1), Type(1), Type(1), + Type(1), Type(2), Type(3) + ) + ); + } + { + Info<< "# Row-vector access:" << nl; + T = Tbak; + cmp(" x():", T.x(), Vector(Type(1), Type(2), Type(-3))); + cmp(" y():", T.y(), Vector(Type(4), Type(5), Type(-6))); + cmp(" z():", T.z(), Vector(Type(7), Type(8), Type(-9))); + cmp(" row(0):", T.row(0), Vector(Type(1), Type(2), Type(-3))); + cmp(" row(1):", T.row(1), Vector(Type(4), Type(5), Type(-6))); + cmp(" row(2):", T.row(2), Vector(Type(7), Type(8), Type(-9))); + cmp + ( + " row<0>:", + T.template row<0>(), + Vector(Type(1), Type(2), Type(-3)) + ); + cmp + ( + " row<1>:", + T.template row<1>(), + Vector(Type(4), Type(5), Type(-6)) + ); + cmp + ( + " row<2>:", + T.template row<2>(), + Vector(Type(7), Type(8), Type(-9)) + ); + // Compilation error: Info << " row<3> = " << T.row<3>() << nl; + + + Info<< "# Row-vector manipulation:" << nl; + T.row(1, Vector(Type(0), Type(1), Type(99))); + cmp + ( + " row(1, Vector):", + T.row(1), + Vector(Type(0), Type(1), Type(99)) + ); + + T.rows + ( + Vector(Type(1), Type(1), Type(1)), + Vector(Type(-1), Type(1), Type(2)), + Vector(Type(1), Type(1), Type(3)) + ); + cmp + ( + " rows(Vectors):", + T, + Tensor + ( + Type(1), Type(1), Type(1), + Type(-1), Type(1), Type(2), + Type(1), Type(1), Type(3) + ) + ); + } + { + Info<< "# Diagonal access:" << nl; + + T = Tbak; + cmp + ( + " 'Tensor'.diag():", + T.diag(), + Vector(Type(1), Type(5), Type(-9)) + ); + cmp + ( + " 'const Tensor'.diag():", + cT.diag(), + Vector(Type(-9), Type(5), Type(1)) + ); + + + Info<< "# Diagonal manipulation:" << nl; + + T.diag(Vector(Type(-10), Type(-15), Type(-20))); + cmp + ( + " 'Tensor'.diag('Vector'):", + T.diag(), + Vector(Type(-10), Type(-15), Type(-20)) + ); + } + { + Info<< "# Tensor operations:" << nl; + + T = Tbak; + cmp(" Transpose:", T, (T.T()).T()); + cmp // Singular matrix + ( + " Inverse:", + T.inv(), + Tensor + ( + Type(-4.50359963e+15), Type(9.00719925e+15), Type(-4.50359963e+15), + Type(9.00719925e+15), Type(-1.80143985e+16), Type(9.00719925e+15), + Type(4.50359963e+15), Type(-9.00719925e+15), Type(4.50359963e+15) + ) + ); + cmp + ( + " Inner-product:", + T.inner(T), + Tensor + ( + Type(-12), Type(-12), Type(12), + Type(-18), Type(-15), Type(12), + Type(-24), Type(-18), Type(12) + ) + ); + cmp + ( + " Schur-product:", + T.schur(T), + Tensor + ( + Type(1), Type(4), Type(9), + Type(16), Type(25), Type(36), + Type(49), Type(64), Type(81) + ) + ); + } + { + Info<< "# Member operators:" << nl; + + T = Tbak; + T &= T; + cmp + ( + " Assign inner-product of this with another Tensor:", + T, + (Tbak & Tbak) + ); + + T = VectorSpace, Type, 9>(Zero); + cmp + ( + " Assign to an equivalent vector space:", + T, + Tensor(Zero) + ); + + T = SphericalTensor(Type(5)); + cmp + ( + " Assign to a SphericalTensor:", + T, + Tensor + ( + Type(5), Zero, Zero, + Zero, Type(5), Zero, + Zero, Zero, Type(5) + ) + ); + + T = SymmTensor + ( + Type(1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + cmp + ( + " Assign to a SymmTensor:", + T, + Tensor + ( + Type(1), Type(2), Type(-3), + Type(2), Type(5), Type(-6), + Type(-3), Type(-6), Type(-9) + ) + ); + + T = Vector> + ( + Vector(Type(-1), Type(2), Type(3)), + Vector(Type(-4), Type(5), Type(6)), + Vector(Type(4), Type(-5), Type(6)) + ); + cmp + ( + " Assign to a triad of row vectors:", + T, + Tensor + ( + Type(-1), Type(2), Type(3), + Type(-4), Type(5), Type(6), + Type(4), Type(-5), Type(6) + ) + ); + } +} + + +// Execute each global function of Tensor, and print output +template +void test_global_funcs(Type) +{ + const Tensor T + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(-1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + + Info<< "# Operands: " << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << endl; + + + cmp(" Trace = ", tr(T), Type(-5)); + cmp(" Spherical part = ", sph(T), SphericalTensor(tr(T)/Type(3))); + cmp + ( + " Symmetric part = ", + symm(T), + SymmTensor + ( + Type(-1), Type(3), Type(2), + Type(5), Type(1), + Type(-9) + ) + ); + cmp + ( + " Twice the symmetric part = ", + twoSymm(T), + SymmTensor + ( + Type(-2), Type(6), Type(4), + Type(10), Type(2), + Type(-18) + ) + ); + cmp + ( + " Skew-symmetric part = ", + skew(T), + Tensor + ( + Type(0), Type(-1), Type(-5), + Type(1), Type(0), Type(-7), + Type(5), Type(7), Type(0) + ) + ); + /* + // Complex-type is not supported for this function. + cmp + ( + " Skew-symmetric part of a SymmTensor = ", + skew(sT), + Tensor(Zero) + ); + */ + cmp + ( + " Deviatoric part = ", + dev(T), + Tensor + ( + Type(0.66666667), Type(2), Type(-3), + Type(4), Type(6.66666667), Type(-6), + Type(7), Type(8), Type(-7.33333333) + ), + 1e-6, + 1e-6 + ); + cmp(" Two-third deviatoric part = ", dev2(T), T - 2*sph(T)); + cmp(" Determinant = ", det(T), Type(-6.000000000000005)); + cmp + ( + " Cofactor tensor = ", + cof(T), + Tensor + ( + Type(3), Type(-6), Type(-3), + Type(-6), Type(30), Type(22), + Type(3), Type(-18), Type(-13) + ) + ); + cmp + ( + " Inverse = ", + inv(T, det(T)), + Tensor + ( + Type(-0.5), Type(1), Type(-0.5), + Type(1), Type(-5), Type(3), + Type(0.5), Type(-3.66666667), Type(2.16666667) + ), + 1e-8, + 1e-8 + ); + cmp + ( + " Inverse (another) = ", + inv(T), + Tensor + ( + Type(-0.5), Type(1), Type(-0.5), + Type(1), Type(-5), Type(3), + Type(0.5), Type(-3.66666667), Type(2.16666667) + ), + 1e-8, + 1e-8 + ); + cmp + ( + " Inverse (another) = ", + T.inv(), + Tensor + ( + Type(-0.5), Type(1), Type(-0.5), + Type(1), Type(-5), Type(3), + Type(0.5), Type(-3.66666667), Type(2.16666667) + ), + 1e-8, + 1e-8 + ); + cmp(" First invariant = ", invariantI(T), Type(-5)); + cmp(" Second invariant = ", invariantII(T), Type(20)); + cmp(" Third invariant = ", invariantIII(T), Type(-6.000000000000005)); +} + + +// Execute each global operator of Tensor, and print output +template +void test_global_opers(Type) +{ + const Tensor T + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ); + const SymmTensor sT + ( + Type(-1), Type(2), Type(-3), + Type(5), Type(-6), + Type(-9) + ); + const SphericalTensor spT(Type(1)); + const Vector v(Type(3), Type(2), Type(1)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor = " << T << nl + << " SymmTensor = " << sT << nl + << " SphericalTensor = " << spT << nl + << " Vector = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Sum of SpTensor-Tensor = ", + (spT + T), + Tensor + ( + Type(0), Type(2), Type(-3), + Type(4), Type(6), Type(-6), + Type(7), Type(8), Type(-8) + ) + ); + cmp + ( + " Sum of Tensor-SpTensor = ", + (T + spT), + Tensor + ( + Type(0), Type(2), Type(-3), + Type(4), Type(6), Type(-6), + Type(7), Type(8), Type(-8) + ) + ); + cmp + ( + " Sum of SymmTensor-Tensor = ", + (sT + T), + Tensor + ( + Type(-2), Type(4), Type(-6), + Type(6), Type(10), Type(-12), + Type(4), Type(2), Type(-18) + ) + ); + cmp + ( + " Sum of Tensor-SymmTensor = ", + (T + sT), + Tensor + ( + Type(-2), Type(4), Type(-6), + Type(6), Type(10), Type(-12), + Type(4), Type(2), Type(-18) + ) + ); + cmp + ( + " Subtract Tensor from SpTensor = ", + (spT - T), + Tensor + ( + Type(2), Type(-2), Type(3), + Type(-4), Type(-4), Type(6), + Type(-7), Type(-8), Type(10) + ) + ); + cmp + ( + " Subtract SpTensor from Tensor = ", + (T - spT), + Tensor + ( + Type(-2), Type(2), Type(-3), + Type(4), Type(4), Type(-6), + Type(7), Type(8), Type(-10) + ) + ); + cmp + ( + " Subtract Tensor from SymmTensor = ", + (sT - T), + Tensor + ( + Type(0), Type(0), Type(0), + Type(-2), Type(0), Type(0), + Type(-10), Type(-14), Type(0) + ) + ); + cmp + ( + " Subtract SymmTensor from Tensor = ", + (T - sT), + Tensor + ( + Type(0), Type(0), Type(0), + Type(2), Type(0), Type(0), + Type(10), Type(14), Type(0) + ) + ); + cmp + ( + " Hodge dual of a Tensor = ", + *T, + Vector(T.yz(), -T.xz(), T.xy()) + ); + cmp + ( + " Hodge dual of a Vector = ", + *v, + Tensor + ( + Zero, -v.z(), v.y(), + v.z(), Zero, -v.x(), + -v.y(), v.x(), Zero + ) + ); + /*cmp + ( + " Division of Vector by Tensor = ", + (v/T), + Tensor + ( + Type(-3), Type(1), Type(-0.33333333), + Type(0.75), Type(0.4), Type(-0.16666667), + Type(0.42857143), Type(0.25), Type(-0.11111111) + ) + );*/ + cmp + ( + " Division of Tensor by Type = ", + (T/x), + Tensor + ( + Type(-0.25), Type(0.5), Type(-0.75), + Type(1), Type(1.25), Type(-1.5), + Type(1.75), Type(2), Type(-2.25) + ) + ); + cmp + ( + " Inner-product of Tensor-Tensor = ", + (T & T), + Tensor + ( + Type(-12), Type(-16), Type(18), + Type(-26), Type(-15), Type(12), + Type(-38), Type(-18), Type(12) + ) + ); + cmp + ( + " Inner-product of SpTensor-Tensor = ", + (spT & T), + Tensor + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ) + ); + cmp + ( + " Inner-product of Tensor-SpTensor = ", + (T & spT), + Tensor + ( + Type(-1), Type(2), Type(-3), + Type(4), Type(5), Type(-6), + Type(7), Type(8), Type(-9) + ) + ); + cmp + ( + " Inner-product of SymmTensor-Tensor = ", + (sT & T), + Tensor + ( + Type(-12), Type(-16), Type(18), + Type(-24), Type(-19), Type(18), + Type(-84), Type(-108), Type(126) + ) + ); + cmp + ( + " Inner-product of Tensor-SymmTensor = ", + (T & sT), + Tensor + ( + Type(14), Type(26), Type(18), + Type(24), Type(69), Type(12), + Type(36), Type(108), Type(12) + ) + ); + cmp + ( + " Inner-product of Tensor-Vector = ", + (T & v), + Vector(Type(-2), Type(16), Type(28)) // Column-vector + ); + cmp + ( + " Inner-product of Vector-Tensor = ", + (v & T), + Vector(Type(12), Type(24), Type(-30)) // Row-vector + ); + cmp(" D-inner-product of SpTensor-Tensor = ", (spT && T), Type(-5)); + cmp(" D-inner-product of Tensor-SpTensor = ", (T && spT), Type(-5)); + cmp(" D-inner-product of SymmTensor-Tensor = ", (sT && T), Type(95)); + cmp(" D-inner-product of Tensor-SymmTensor = ", (T && sT), Type(95)); + cmp + ( + " Outer-product of Vector-Vector = ", + (v*v), + Tensor + ( + Type(9), Type(6), Type(3), + Type(6), Type(4), Type(2), + Type(3), Type(2), Type(1) + ) + ); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main() +{ + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "Tensor", + "Tensor", + "Tensor" + }); + + run_tests(types, typeID); + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/test/Tensor2D/Make/files b/applications/test/Tensor2D/Make/files new file mode 100644 index 0000000000..76d3b20002 --- /dev/null +++ b/applications/test/Tensor2D/Make/files @@ -0,0 +1,3 @@ +Test-Tensor2D.C + +EXE = $(FOAM_USER_APPBIN)/Test-Tensor2D diff --git a/applications/test/Tensor2D/Make/options b/applications/test/Tensor2D/Make/options new file mode 100644 index 0000000000..18e6fe47af --- /dev/null +++ b/applications/test/Tensor2D/Make/options @@ -0,0 +1,2 @@ +/* EXE_INC = */ +/* EXE_LIBS = */ diff --git a/applications/test/Tensor2D/Test-Tensor2D.C b/applications/test/Tensor2D/Test-Tensor2D.C new file mode 100644 index 0000000000..4b8950e611 --- /dev/null +++ b/applications/test/Tensor2D/Test-Tensor2D.C @@ -0,0 +1,841 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2014 OpenFOAM Foundation + Copyright (C) 2019-2020 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM 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. + + OpenFOAM 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 + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Application + Test-Tensor2D + +Description + Tests for \c Tensor2D constructors, member functions and operators + using \c floatScalar, \c doubleScalar, and \c complex base types. + + Eigen decomposition tests for \c tensor2D, i.e. Tensor2D. + + Cross-checks were obtained from 'NumPy 1.15.1' and 'SciPy 1.1.0' if no + theoretical cross-check exists (like eigendecomposition relations), and + were hard-coded for elementwise comparisons. + + For \c complex base type, the cross-checks do only involve zero imag part. + +\*---------------------------------------------------------------------------*/ + +#include "vector2DField.H" +#include "tensor2D.H" +#include "symmTensor2D.H" +#include "transform.H" +#include "Random.H" +#include "floatScalar.H" +#include "doubleScalar.H" +#include "complex.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// Total number of unit tests +unsigned nTest_ = 0; + + +// Total number of failed unit tests +unsigned nFail_ = 0; + + +// Create a random tensor2D +tensor2D makeRandomContainer(Random& rnd) +{ + tensor2D A(Zero); + std::generate(A.begin(), A.end(), [&]{ return rnd.GaussNormal(); }); + return A; +} + + +// Compare two floating point types, and print output. +// Do ++nFail_ if values of two objects are not equal within a given tolerance. +// The function is converted from PEP-485. +template +typename std::enable_if +< + std::is_same::value || + std::is_same::value || + std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, // +typename std::enable_if +< + !std::is_same::value && + !std::is_same::value && + !std::is_same::value, + void +>::type cmp +( + const word& msg, + const Type& x, + const Type& y, + const scalar relTol = 1e-8, + const scalar absTol = 0 +) +{ + Info<< msg << x << endl; + + unsigned nFail = 0; + + for (label i = 0; i < pTraits::nComponents; ++i) + { + if (max(absTol, relTol*max(mag(x[i]), mag(y[i]))) < mag(x[i] - y[i])) + { + ++nFail; + } + } + + if (nFail) + { + Info<< nl + << " #### Fail in " << nFail << " comps ####" << nl << endl; + ++nFail_; + } + ++nTest_; +} + + +// Create each constructor of Tensor2D, and print output +template +void test_constructors(Type) +{ + { + Info<< "# Construct initialized to zero:" << nl; + const Tensor2D T(Zero); + Info<< T << endl; + } + { + Info<< "# Construct given VectorSpace:" << nl; + const VectorSpace, Type, 4> V(Zero); + const Tensor2D T(V); + Info<< T << endl; + } + { + Info<< "# Construct given SymmTensor2D:" << nl; + const SymmTensor2D S + ( + Type(1), Type(2), + Type(3) + ); + const Tensor2D T(S); + Info<< T << endl; + } + { + Info<< "# Construct given SphericalTensor2D:" << nl; + const SphericalTensor2D Sp(Type(5)); + const Tensor2D T(Sp); + Info<< T << endl; + } + { + Info<< "# Construct given the two row vectors:" << nl; + const Vector2D x(Type(1), Type(2)); + const Vector2D y(Type(3), Type(4)); + const Tensor2D T(x, y); + Info<< T << endl; + } + { + Info<< "# Construct given the four components:" << nl; + const Tensor2D T + ( + Type(1), Type(2), + Type(3), Type(4) + ); + Info<< T << endl; + } + { + Info<< "# Copy construct:" << nl; + const Tensor2D T(Zero); + const Tensor2D Tcopy(T); + Info<< T << endl; + } +} + + +// Execute each member function of Tensor2D, and print output +template +void test_member_funcs(Type) +{ + Tensor2D T + ( + Type(1), Type(2), + Type(4), Type(5) + ); + Tensor2D Tbak = T; + const Tensor2D cT + ( + Type(-9), Type(8), + Type(-6), Type(5) + ); + + Info<< "# Operand: " << nl + << " Tensor2D = " << T << endl; + + + { + Info<< "# Component access:" << nl; + + Tensor2D cpT + ( + T.xx(), T.xy(), + T.yx(), T.yy() + ); + cmp(" 'Tensor2D' access:", T, cpT); + + const Tensor2D cpcT + ( + cT.xx(), cT.xy(), + cT.yx(), cT.yy() + ); + cmp(" 'const Tensor2D' access:", cT, cpcT); + } + { + Info<< "# Column-vector access:" << nl; + cmp(" cx():", T.cx(), Vector2D(Type(1), Type(4))); + cmp(" cy():", T.cy(), Vector2D(Type(2), Type(5))); + cmp(" col(0):", T.col(0), Vector2D(Type(1), Type(4))); + cmp(" col(1):", T.col(1), Vector2D(Type(2), Type(5))); + cmp + ( + " col<0>:", + T.template col<0>(), + Vector2D(Type(1), Type(4)) + ); + cmp + ( + " col<1>:", + T.template col<1>(), + Vector2D(Type(2), Type(5)) + ); + // Compilation error: Info << " col<2> = " << T.col<2>() << nl; + + + Info<< "# Column-vector manipulation:" << nl; + T.col(1, Vector2D(Type(0), Type(1))); + cmp + ( + " col(1, Vector):", + T.col(1), + Vector2D(Type(0), Type(1)) + ); + + T.cols + ( + Vector2D(Type(1), Type(1)), + Vector2D(Type(-1), Type(1)) + ); + cmp + ( + " cols(Vectors):", + T, + Tensor2D + ( + Type(1), Type(-1), + Type(1), Type(1) + ) + ); + } + { + Info<< "# Row-vector access:" << nl; + T = Tbak; + cmp(" x():", T.x(), Vector2D(Type(1), Type(2))); + cmp(" y():", T.y(), Vector2D(Type(4), Type(5))); + cmp(" row(0):", T.row(0), Vector2D(Type(1), Type(2))); + cmp(" row(1):", T.row(1), Vector2D(Type(4), Type(5))); + cmp + ( + " row<0>:", + T.template row<0>(), + Vector2D(Type(1), Type(2)) + ); + cmp + ( + " row<1>:", + T.template row<1>(), + Vector2D(Type(4), Type(5)) + ); + // Compilation error: Info << " row<2> = " << T.row<2>() << nl; + + + Info<< "# Row-vector manipulation:" << nl; + T.row(1, Vector2D(Type(0), Type(1))); + cmp + ( + " row(1, Vector):", + T.row(1), + Vector2D(Type(0), Type(1)) + ); + + T.rows + ( + Vector2D(Type(1), Type(1)), + Vector2D(Type(-1), Type(1)) + ); + cmp + ( + " rows(Vectors):", + T, + Tensor2D + ( + Type(1), Type(1), + Type(-1), Type(1) + ) + ); + } + { + Info<< "# Diagonal access:" << nl; + + T = Tbak; + cmp + ( + " 'Tensor2D'.diag():", + T.diag(), + Vector2D(Type(1), Type(5)) + ); + cmp + ( + " 'const Tensor2D'.diag():", + cT.diag(), + Vector2D(Type(-9), Type(5)) + ); + + + Info<< "# Diagonal manipulation:" << nl; + + T.diag(Vector2D(Type(-10), Type(-15))); + cmp + ( + " 'Tensor2D'.diag('Vector'):", + T.diag(), + Vector2D(Type(-10), Type(-15)) + ); + } + { + Info<< "# Tensor operations:" << nl; + + T = Tbak; + cmp(" Transpose:", T, (T.T()).T()); + cmp + ( + " Inner-product:", + T.inner(T), + Tensor2D + ( + Type(9), Type(12), + Type(24), Type(33) + ) + ); + cmp + ( + " Schur-product:", + T.schur(T), + Tensor2D + ( + Type(1), Type(4), + Type(16), Type(25) + ) + ); + } + { + Info<< "# Member operators:" << nl; + + T = SphericalTensor2D(Type(5)); + cmp + ( + " Assign to a SphericalTensor2D:", + T, + Tensor2D + ( + Type(5), Zero, + Zero, Type(5) + ) + ); + + T = SymmTensor2D + ( + Type(1), Type(2), + Type(5) + ); + cmp + ( + " Assign to a SymmTensor2D:", + T, + Tensor2D + ( + Type(1), Type(2), + Type(2), Type(5) + ) + ); + } +} + + +// Execute each global function of Tensor2D, and print output +template +void test_global_funcs(Type) +{ + const Tensor2D T + ( + Type(-1), Type(2), + Type(4), Type(5) + ); + + Info<< "# Operand: " << nl + << " Tensor2D = " << T << endl; + + + cmp(" Trace = ", tr(T), Type(4)); + cmp(" Spherical part = ", sph(T), SphericalTensor2D(tr(T)/Type(2))); + cmp + ( + " Symmetric part = ", + symm(T), + SymmTensor2D + ( + Type(-1), Type(3), + Type(5) + ) + ); + cmp + ( + " Twice the symmetric part = ", + twoSymm(T), + SymmTensor2D + ( + Type(-2), Type(6), + Type(10) + ) + ); + cmp + ( + " Skew-symmetric part = ", + skew(T), + Tensor2D + ( + Type(0), Type(-1), + Type(1), Type(0) + ) + ); + cmp + ( + " Deviatoric part = ", + dev(T), + Tensor2D + ( + Type(-3), Type(2), + Type(4), Type(3) + ) + ); + cmp(" Two-third deviatoric part = ", dev2(T), T - 2*sph(T)); + cmp(" Determinant = ", det(T), Type(-13)); + cmp + ( + " Cofactor tensor2D = ", + cof(T), + Tensor2D + ( + Type(5), Type(-4), + Type(-2), Type(-1) + ) + ); + cmp + ( + " Inverse = ", + inv(T, det(T)), + Tensor2D + ( + Type(-0.38461538), Type(0.15384615), + Type(0.30769231), Type(0.07692308) + ), + 1e-6, + 1e-6 + ); + cmp + ( + " Inverse (another) = ", + inv(T), + Tensor2D + ( + Type(-0.38461538), Type(0.15384615), + Type(0.30769231), Type(0.07692308) + ), + 1e-6, + 1e-6 + ); + cmp(" First invariant = ", invariantI(T), Type(4)); + cmp(" Second invariant = ", invariantII(T), Type(-13)); +} + + +// Execute each global operator of Tensor2D, and print output +template +void test_global_opers(Type) +{ + const Tensor2D T + ( + Type(-1), Type(2), + Type(4), Type(5) + ); + const SymmTensor2D sT + ( + Type(1), Type(2), + Type(5) + ); + const SphericalTensor2D spT(Type(1)); + const Vector2D v(Type(3), Type(2)); + const Type x(4); + + Info<< "# Operands:" << nl + << " Tensor2D = " << T << nl + << " SymmTensor2D = " << sT << nl + << " SphericalTensor2D = " << spT << nl + << " Vector2D = " << v << nl + << " Type = " << x << endl; + + + cmp + ( + " Sum of SpTensor2D-Tensor2D = ", + (spT + T), + Tensor2D + ( + Type(0), Type(2), + Type(4), Type(6) + ) + ); + cmp + ( + " Sum of Tensor2D-SpTensor2D = ", + (T + spT), + Tensor2D + ( + Type(0), Type(2), + Type(4), Type(6) + ) + ); + cmp + ( + " Sum of SymmTensor2D-Tensor2D = ", + (sT + T), + Tensor2D + ( + Type(0), Type(4), + Type(6), Type(10) + ) + ); + cmp + ( + " Sum of Tensor2D-SymmTensor2D = ", + (T + sT), + Tensor2D + ( + Type(0), Type(4), + Type(6), Type(10) + ) + ); + cmp + ( + " Subtract Tensor2D from SpTensor2D = ", + (spT - T), + Tensor2D + ( + Type(2), Type(-2), + Type(-4), Type(-4) + ) + ); + cmp + ( + " Subtract SpTensor2D from Tensor2D = ", + (T - spT), + Tensor2D + ( + Type(-2), Type(2), + Type(4), Type(4) + ) + ); + cmp + ( + " Subtract Tensor2D from SymmTensor2D = ", + (sT - T), + Tensor2D + ( + Type(2), Type(0), + Type(-2), Type(0) + ) + ); + cmp + ( + " Subtract SymmTensor2D from Tensor2D = ", + (T - sT), + Tensor2D + ( + Type(-2), Type(0), + Type(2), Type(0) + ) + ); + cmp + ( + " Division of Tensor2D by Type = ", + (T/x), + Tensor2D + ( + Type(-0.25), Type(0.5), + Type(1), Type(1.25) + ) + ); + cmp + ( + " Inner-product of Tensor2D-Tensor2D = ", + (T & T), + Tensor2D + ( + Type(9), Type(8), + Type(16), Type(33) + ) + ); + cmp + ( + " Inner-product of SpTensor2D-Tensor2D = ", + (spT & T), + Tensor2D + ( + Type(-1), Type(2), + Type(4), Type(5) + ) + ); + cmp + ( + " Inner-product of Tensor2D-SpTensor2D = ", + (T & spT), + Tensor2D + ( + Type(-1), Type(2), + Type(4), Type(5) + ) + ); + cmp + ( + " Inner-product of SymmTensor2D-Tensor2D = ", + (sT & T), + Tensor2D + ( + Type(7), Type(12), + Type(18), Type(29) + ) + ); + cmp + ( + " Inner-product of Tensor2D-SymmTensor2D = ", + (T & sT), + Tensor2D + ( + Type(3), Type(8), + Type(14), Type(33) + ) + ); + cmp + ( + " Inner-product of Tensor2D-Vector2D = ", + (T & v), + Vector2D(Type(1), Type(22)) // Column-vector + ); + cmp + ( + " Inner-product of Vector2D-Tensor2D = ", + (v & T), + Vector2D(Type(5), Type(16)) // Row-vector + ); + cmp(" D-inner-product of SpTensor2D-Tensor2D = ", (spT && T), Type(4)); + cmp(" D-inner-product of Tensor2D-SpTensor2D = ", (T && spT), Type(4)); + cmp(" D-inner-product of SymmTensor2D-Tensor2D = ", (sT && T), Type(36)); + cmp(" D-inner-product of Tensor2D-SymmTensor2D = ", (T && sT), Type(36)); + cmp + ( + " Outer-product of Vector2D-Vector2D = ", + (v*v), + Tensor2D + ( + Type(9), Type(6), + Type(6), Type(4) + ) + ); +} + + +// Do compile-time recursion over the given types +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID){} + + +template +inline typename std::enable_if::type +run_tests(const std::tuple& types, const List& typeID) +{ + Info<< nl << " ## Test constructors: "<< typeID[I] <<" ##" << nl; + test_constructors(std::get(types)); + + Info<< nl << " ## Test member functions: "<< typeID[I] <<" ##" << nl; + test_member_funcs(std::get(types)); + + Info<< nl << " ## Test global functions: "<< typeID[I] << " ##" << nl; + test_global_funcs(std::get(types)); + + Info<< nl << " ## Test global operators: "<< typeID[I] <<" ##" << nl; + test_global_opers(std::get(types)); + + run_tests(types, typeID); +} + + +// * * * * * * * * * * * * * * * Main Program * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + + const std::tuple types + ( + std::make_tuple(Zero, Zero, Zero) + ); + + const List typeID + ({ + "Tensor2D", + "Tensor2D", + "Tensor2D" + }); + + run_tests(types, typeID); + + + { + Info<< "# Pre-v2006 tests:" << nl; + vector2D v1(1, 2), v2(3, 4); + tensor2D t3 = v1*v2; + + Info<< v1 << "*" << v2 << " = " << t3 << endl; + + { + Info<< "rows:" << nl; + for (direction i = 0; i < 2; ++i) + { + Info<< " (" << i << ") = " << t3.row(i) << nl; + } + } + + { + Info<< "cols:" << nl; + for (direction i = 0; i < 2; ++i) + { + Info<< " (" << i << ") = " << t3.col(i) << nl; + } + Info<< "col<0> = " << t3.col<0>() << nl; + Info<< "col<1> = " << t3.col<1>() << nl; + // Compilation error: Info << "col<3> = " << t3.col<3>() << nl; + + t3.col<0>({0, 2}); + Info<< "replaced col<0> = " << t3.col<0>() << nl; + Info<< "tensor " << t3 << nl; + + t3.row<1>(Zero); + Info<< "replaced row<1> = " << t3.row<1>() << nl; + Info<< "tensor " << t3 << nl; + } + + + { + vector2DField vfld1(8, Zero); + + forAll(vfld1, i) + { + vfld1[i] = (i + 1) * ((i % 2) ? v1 : v2); + } + + Info<< "vector: " << flatOutput(vfld1) << nl; + + scalarField xvals(8); + scalarField yvals(8); + unzip(vfld1, xvals, yvals); + + Info<< "unzip" << nl + << " x => " << flatOutput(xvals) << nl + << " y => " << flatOutput(yvals) << nl; + + reverse(xvals); + zip(vfld1, xvals, yvals); + + Info<< "rezip (with reversed x)" << nl + << " => " << flatOutput(vfld1) << nl; + } + } + + + if (nFail_) + { + Info<< nl << " #### " + << "Failed in " << nFail_ << " tests " + << "out of total " << nTest_ << " tests " + << "####\n" << endl; + return 1; + } + + Info<< nl << " #### Passed all " << nTest_ <<" tests ####\n" << endl; + return 0; +} + + +// ************************************************************************* // diff --git a/applications/test/tensor/Make/files b/applications/test/tensor/Make/files deleted file mode 100644 index 34d6a3cdac..0000000000 --- a/applications/test/tensor/Make/files +++ /dev/null @@ -1,3 +0,0 @@ -Test-tensor.C - -EXE = $(FOAM_USER_APPBIN)/Test-tensor diff --git a/applications/test/tensor/Test-tensor.C b/applications/test/tensor/Test-tensor.C deleted file mode 100644 index 9371be3439..0000000000 --- a/applications/test/tensor/Test-tensor.C +++ /dev/null @@ -1,235 +0,0 @@ -#include "tensor.H" -#include "triad.H" -#include "symmTensor.H" -#include "transform.H" -#include "stringList.H" -#include "IOstreams.H" - -using namespace Foam; - -int main() -{ - tensor t1(1, 2, 3, 4, 5, 6, 7, 8, 9); - tensor t2(1, 2, 3, 1, 2, 3, 1, 2, 3); - - tensor t3 = t1 + t2; - - Info<< "tensor " << t1 << " + " << t2 << nl - << " = " << t3 << nl << nl; - - { - Info<< "rows:" << nl; - for (direction i=0; i < 3; ++i) - { - Info<< " (" << i << ") = " << t3.row(i) << nl; - } - } - { - Info<< "cols:" << nl; - for (direction i=0; i < 3; ++i) - { - Info<< " (" << i << ") = " << t3.col(i) << nl; - } - Info<< "col<0> = " << t3.col<0>() << nl; - Info<< "col<1> = " << t3.col<1>() << nl; - Info<< "col<2> = " << t3.col<2>() << nl; - // Compilation error: Info << "col<3> = " << t3.col<3>() << nl; - - t3.col<1>({0, 2, 4}); - Info<< "replaced col<1> = " << t3.col<1>() << nl; - Info<< "tensor " << t3 << nl; - - t3.row<2>(Zero); - Info<< "replaced row<2> = " << t3.row<2>() << nl; - Info<< "tensor " << t3 << nl; - - triad tr3(t3); - Info<< "triad " << tr3 << " :: T() " << tr3.T() << nl; - } - Info<< nl; - - - tensor t4(3,-2,1,-2,2,0,1, 0, 4); - - Info<< inv(t4) << endl; - Info<< (inv(t4) & t4) << endl; - - Info<< t1.x() << t1.y() << t1.z() << endl; - - tensor t6(1,0,-4,0,5,4,-4,4,3); - //tensor t6(1,2,0,2,5,0,0,0,0); - - Info<< "tensor " << t6 << endl; - vector e = eigenValues(t6); - Info<< "eigenvalues " << e << endl; - - tensor ev = eigenVectors(t6); - Info<< "eigenvectors " << ev << endl; - - Info<< "Check determinant " << e.x()*e.y()*e.z() << " " << det(t6) << endl; - - Info<< "Check eigenvectors " - << (eigenVectors(t6, e) & t6) << " " - << (e.x() * eigenVectors(t6, e)).x() - << (e.x() * eigenVectors(t6, e)).y() - << (e.z() * eigenVectors(t6, e)).z() - << endl; - - Info<< "Check eigenvalues for symmTensor " - << eigenValues(symm(t6)) - eigenValues(tensor(symm(t6))) << endl; - - Info<< "Check eigenvectors for symmTensor " - << eigenVectors(symm(t6)) - eigenVectors(tensor(symm(t6))) << endl; - - tensor t7(1, 2, 3, 2, 4, 5, 3, 5, 6); - - Info<< "Check transformation " - << (t1 & t7 & t1.T()) << " " << transform(t1, t7) << endl; - - symmTensor st1(1, 2, 3, 4, 5, 6); - symmTensor st2(7, 8, 9, 10, 11, 12); - - Info<< "Check symmetric transformation " - << transform(t1, st1) << endl; - - Info<< "Check dot product of symmetric tensors " - << (st1 & st2) << endl; - - Info<< "Check inner sqr of a symmetric tensor " - << innerSqr(st1) << " " << innerSqr(st1) - (st1 & st1) << endl; - - Info<< "Check symmetric part of dot product of symmetric tensors " - << twoSymm(st1&st2) - ((st1&st2) + (st2&st1)) << endl; - - tensor sk1 = skew(t6); - Info<< "Check dot product of symmetric and skew tensors " - << twoSymm(st1&sk1) - ((st1&sk1) - (sk1&st1)) << endl; - - vector v1(1, 2, 3); - - Info<< sqr(v1) << endl; - Info<< symm(t7) << endl; - Info<< twoSymm(t7) << endl; - Info<< magSqr(st1) << endl; - Info<< mag(st1) << endl; - - Info<< (symm(t7) && t7) - (0.5*(t7 + t7.T()) && t7) << endl; - Info<< (t7 && symm(t7)) - (t7 && 0.5*(t7 + t7.T())) << endl; - - - /* - // Lots of awkward eigenvector tests ... - - tensor T_rand_real - ( - 0.9999996423721313, 0.3330855667591095, 0.6646450161933899, - 0.9745196104049683, 0.0369445420801640, 0.0846728682518005, - 0.6474838852882385, 0.1617118716239929, 0.2041363865137100 - ); - DebugVar(T_rand_real); - vector L_rand_real(eigenValues(T_rand_real)); - DebugVar(L_rand_real); - tensor U_rand_real(eigenVectors(T_rand_real)); - DebugVar(U_rand_real); - - Info << endl << endl; - - tensor T_rand_imag - ( - 0.8668024539947510, 0.1664607226848602, 0.8925783634185791, - 0.9126510620117188, 0.7408077120780945, 0.1499115079641342, - 0.0936608463525772, 0.7615650296211243, 0.8953040242195129 - ); - DebugVar(T_rand_imag); - vector L_rand_imag(eigenValues(T_rand_imag)); - DebugVar(L_rand_imag); - tensor U_rand_imag(eigenVectors(T_rand_imag)); - DebugVar(U_rand_imag); - - Info << endl << endl; - - tensor T_rand_symm - ( - 1.9999992847442627, 1.3076051771640778, 1.3121289014816284, - 1.3076051771640778, 0.0738890841603279, 0.2463847398757935, - 1.3121289014816284, 0.2463847398757935, 0.4082727730274200 - ); - DebugVar(T_rand_symm); - vector L_rand_symm(eigenValues(T_rand_symm)); - DebugVar(L_rand_symm); - tensor U_rand_symm(eigenVectors(T_rand_symm)); - DebugVar(U_rand_symm); - - Info << endl << endl; - - symmTensor T_rand_Symm - ( - 1.9999992847442627, 1.3076051771640778, 1.3121289014816284, - 0.0738890841603279, 0.2463847398757935, - 0.4082727730274200 - ); - DebugVar(T_rand_Symm); - vector L_rand_Symm(eigenValues(T_rand_Symm)); - DebugVar(L_rand_Symm); - tensor U_rand_Symm(eigenVectors(T_rand_Symm)); - DebugVar(U_rand_Symm); - - Info << endl << endl; - - tensor T_rand_diag - ( - 0.8668024539947510, 0, 0, - 0, 0.7408077120780945, 0, - 0, 0, 0.8953040242195129 - ); - DebugVar(T_rand_diag); - vector L_rand_diag(eigenValues(T_rand_diag)); - DebugVar(L_rand_diag); - tensor U_rand_diag(eigenVectors(T_rand_diag)); - DebugVar(U_rand_diag); - - Info << endl << endl; - - tensor T_repeated - ( - 0, 1, 1, - 1, 0, 1, - 1, 1, 0 - ); - DebugVar(T_repeated); - vector L_repeated(eigenValues(T_repeated)); - DebugVar(L_repeated); - tensor U_repeated(eigenVectors(T_repeated)); - DebugVar(U_repeated); - - Info << endl << endl; - - tensor T_repeated_zero - ( - 1, 1, 1, - 1, 1, 1, - 1, 1, 1 - ); - DebugVar(T_repeated_zero); - vector L_repeated_zero(eigenValues(T_repeated_zero)); - DebugVar(L_repeated_zero); - tensor U_repeated_zero(eigenVectors(T_repeated_zero)); - DebugVar(U_repeated_zero); - - Info << endl << endl; - - tensor T_triple - ( - 2, 0, 0, - 0, 2, 0, - 0, 0, 2 - ); - DebugVar(T_triple); - vector L_triple(eigenValues(T_triple)); - DebugVar(L_triple); - tensor U_triple(eigenVectors(T_triple)); - DebugVar(U_triple); - */ - - return 0; -} diff --git a/src/OpenFOAM/primitives/DiagTensor/DiagTensor.H b/src/OpenFOAM/primitives/DiagTensor/DiagTensor.H index 5153303f45..648b6989a4 100644 --- a/src/OpenFOAM/primitives/DiagTensor/DiagTensor.H +++ b/src/OpenFOAM/primitives/DiagTensor/DiagTensor.H @@ -28,11 +28,11 @@ Class Foam::DiagTensor Description - Templated 3D DiagTensor derived from VectorSpace. + A templated (3 x 3) diagonal tensor of objects of \, effectively + containing 3 elements, derived from VectorSpace. - Adding construction from 3 components, element access using xx(), yy() - and zz() member functions and the inner-product (dot-product) and - outer-product operators. +See also + Test-DiagTensor.C SourceFiles DiagTensorI.H diff --git a/src/OpenFOAM/primitives/DiagTensor/DiagTensorI.H b/src/OpenFOAM/primitives/DiagTensor/DiagTensorI.H index a160627795..f34e7e3e4f 100644 --- a/src/OpenFOAM/primitives/DiagTensor/DiagTensorI.H +++ b/src/OpenFOAM/primitives/DiagTensor/DiagTensorI.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -116,6 +117,7 @@ namespace Foam // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * // +//- Sum of a DiagTensor and a Tensor template inline Tensor operator+(const DiagTensor& dt1, const Tensor& t2) @@ -129,6 +131,7 @@ operator+(const DiagTensor& dt1, const Tensor& t2) } +//- Sum of a Tensor and a DiagTensor template inline Tensor operator+(const Tensor& t1, const DiagTensor& dt2) @@ -142,6 +145,7 @@ operator+(const Tensor& t1, const DiagTensor& dt2) } +//- Subtract a Tensor from a DiagTensor template inline Tensor operator-(const DiagTensor& dt1, const Tensor& t2) @@ -155,6 +159,7 @@ operator-(const DiagTensor& dt1, const Tensor& t2) } +//- Subtract a DiagTensor from a Tensor template inline Tensor operator-(const Tensor& t1, const DiagTensor& dt2) @@ -168,7 +173,7 @@ operator-(const Tensor& t1, const DiagTensor& dt2) } -//- Inner-product between two diagonal tensors +//- Inner-product of a DiagTensor and a DiagTensor template inline DiagTensor operator&(const DiagTensor& dt1, const DiagTensor& dt2) @@ -182,7 +187,7 @@ operator&(const DiagTensor& dt1, const DiagTensor& dt2) } -//- Inner-product between a diagonal tensor and a tensor +//- Inner-product of a DiagTensor and a Tensor template inline Tensor operator&(const DiagTensor& dt1, const Tensor& t2) @@ -204,7 +209,7 @@ operator&(const DiagTensor& dt1, const Tensor& t2) } -//- Inner-product between a tensor and a diagonal tensor +//- Inner-product of a Tensor and a DiagTensor template inline Tensor operator&(const Tensor& t1, const DiagTensor& dt2) @@ -226,7 +231,7 @@ operator&(const Tensor& t1, const DiagTensor& dt2) } -//- Inner-product between a diagonal tensor and a vector +//- Inner-product of a DiagTensor and a Vector template inline Vector operator&(const DiagTensor& dt, const Vector& v) @@ -240,7 +245,7 @@ operator&(const DiagTensor& dt, const Vector& v) } -//- Inner-product between a vector and a diagonal tensor +//- Inner-product of a Vector and a DiagTensor template inline Vector operator&(const Vector& v, const DiagTensor& dt) @@ -254,25 +259,67 @@ operator&(const Vector& v, const DiagTensor& dt) } -//- Division of a scalar by a diagonalTensor +//- Division of a Cmpt by a DiagTensor template inline DiagTensor -operator/(const scalar s, const DiagTensor& dt) +operator/(const Cmpt s, const DiagTensor& dt) { + #ifdef FULLDEBUG + if (mag(det(dt)) < VSMALL) + { + FatalErrorInFunction + << "Cmpt = " << s + << " is not divisible by the DiagTensor due to a zero element:" + << "DiagTensor = " << dt + << abort(FatalError); + } + #endif + return DiagTensor(s/dt.xx(), s/dt.yy(), s/dt.zz()); } -//- Division of a vector by a diagonalTensor +//- Division of a DiagTensor by a Cmpt +template +inline DiagTensor +operator/(const DiagTensor& dt, const Cmpt s) +{ + #ifdef FULLDEBUG + if (mag(s) < VSMALL) + { + FatalErrorInFunction + << "DiagTensor = " << dt + << " is not divisible due to a zero value in Cmpt:" + << "Cmpt = " << s + << abort(FatalError); + } + #endif + + return DiagTensor(dt.xx()/s, dt.yy()/s, dt.zz()/s); +} + + +//- Division of a Vector by a DiagTensor template inline Vector operator/(const Vector v, const DiagTensor& dt) { + #ifdef FULLDEBUG + if (mag(det(dt)) < VSMALL) + { + FatalErrorInFunction + << "Vector = " << v + << " is not divisible by the DiagTensor due to a zero element:" + << "DiagTensor = " << dt + << abort(FatalError); + } + #endif + return Vector(v.x()/dt.xx(), v.y()/dt.yy(), v.z()/dt.zz()); } -//- Return the trace of a diagonal tensor +//- Return the trace of a DiagTensor template inline Cmpt tr(const DiagTensor& dt) { @@ -280,42 +327,44 @@ inline Cmpt tr(const DiagTensor& dt) } -//- Return the spherical part of a diagonal tensor +//- Return the spherical part of a DiagTensor as a SphericalTensor template inline SphericalTensor sph(const DiagTensor& dt) { return SphericalTensor ( - 1.0/3.0*tr(dt) + (1.0/3.0)*tr(dt) ); } -//- Return the determinant of a diagonal tensor +//- Return the determinant of a DiagTensor template -inline Cmpt det(const DiagTensor& t) +inline Cmpt det(const DiagTensor& dt) { - return t.xx()*t.yy()*t.zz(); + return dt.xx()*dt.yy()*dt.zz(); } -//- Return the inverse of a diagonal tensor +//- Return the inverse of a DiagTensor as a DiagTensor template inline DiagTensor inv(const DiagTensor& dt) { - return DiagTensor(1.0/dt.xx(), 1.0/dt.yy(), 1.0/dt.zz()); + #ifdef FULLDEBUG + if (mag(det(dt)) < VSMALL) + { + FatalErrorInFunction + << "DiagTensor is not invertible due to the zero determinant:" + << "det(DiagTensor) = " << det(dt) + << abort(FatalError); + } + #endif + + return DiagTensor(1/dt.xx(), 1/dt.yy(), 1/dt.zz()); } -//- Return the diagonal of a symmetric tensor as a diagonal tensor -template -inline DiagTensor diag(const SymmTensor& t) -{ - return DiagTensor(t.xx(), t.yy(), t.zz()); -} - - -//- Return the diagonal of a tensor as a diagonal tensor +//- Return the diagonal of a Tensor as a DiagTensor template inline DiagTensor diag(const Tensor& t) { @@ -323,6 +372,14 @@ inline DiagTensor diag(const Tensor& t) } +//- Return the diagonal of a SymmTensor as a DiagTensor +template +inline DiagTensor diag(const SymmTensor& st) +{ + return DiagTensor(st.xx(), st.yy(), st.zz()); +} + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/OpenFOAM/primitives/DiagTensor/diagTensor/diagTensor.H b/src/OpenFOAM/primitives/DiagTensor/diagTensor/diagTensor.H index 892885f2b8..71ac04855f 100644 --- a/src/OpenFOAM/primitives/DiagTensor/diagTensor/diagTensor.H +++ b/src/OpenFOAM/primitives/DiagTensor/diagTensor/diagTensor.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2019-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -28,7 +28,10 @@ Typedef Foam::diagTensor Description - A scalar version of the templated DiagTensor + DiagTensor of scalars, i.e. DiagTensor. + +See also + Test-DiagTensor.C SourceFiles diagTensor.C diff --git a/src/OpenFOAM/primitives/SphericalTensor/SphericalTensor.H b/src/OpenFOAM/primitives/SphericalTensor/SphericalTensor.H index 5d52404dd0..a203f3019c 100644 --- a/src/OpenFOAM/primitives/SphericalTensor/SphericalTensor.H +++ b/src/OpenFOAM/primitives/SphericalTensor/SphericalTensor.H @@ -28,9 +28,11 @@ Class Foam::SphericalTensor Description - Templated 3D SphericalTensor derived from VectorSpace adding construction - from 1 component, element access using th ii() member function and the - inner-product (dot-product) and outer-product operators. + A templated (3 x 3) diagonal tensor of objects of \, effectively + containing 1 element, derived from VectorSpace. + +See also + Test-SphericalTensor.C SourceFiles SphericalTensorI.H @@ -121,7 +123,7 @@ public: inline Cmpt& ii(); - //- Transpose (no-op) + //- Return non-Hermitian transpose (no-op) inline const SphericalTensor& T() const; }; diff --git a/src/OpenFOAM/primitives/SphericalTensor/SphericalTensorI.H b/src/OpenFOAM/primitives/SphericalTensor/SphericalTensorI.H index 8d0c454303..ff577ee653 100644 --- a/src/OpenFOAM/primitives/SphericalTensor/SphericalTensorI.H +++ b/src/OpenFOAM/primitives/SphericalTensor/SphericalTensorI.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -92,7 +93,7 @@ namespace Foam // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * // -//- Inner-product between two spherical tensors +//- Inner-product of a SphericalTensor and a SphericalTensor template inline SphericalTensor operator&(const SphericalTensor& st1, const SphericalTensor& st2) @@ -101,7 +102,7 @@ operator&(const SphericalTensor& st1, const SphericalTensor& st2) } -//- Inner-product between a spherical tensor and a vector +//- Inner-product of a SphericalTensor and a Vector template inline Vector operator&(const SphericalTensor& st, const Vector& v) @@ -115,7 +116,7 @@ operator&(const SphericalTensor& st, const Vector& v) } -//- Inner-product between a vector and a spherical tensor +//- Inner-product of a Vector and a SphericalTensor template inline Vector operator&(const Vector& v, const SphericalTensor& st) @@ -129,7 +130,7 @@ operator&(const Vector& v, const SphericalTensor& st) } -//- Double-dot-product between a spherical tensor and a spherical tensor +//- Double-inner-product of a SphericalTensor and a SphericalTensor template inline Cmpt operator&&(const SphericalTensor& st1, const SphericalTensor& st2) @@ -138,22 +139,55 @@ operator&&(const SphericalTensor& st1, const SphericalTensor& st2) } -//- Division of a scalar by a sphericalTensor +//- Division of a Cmpt by a SphericalTensor template inline SphericalTensor -operator/(const scalar s, const SphericalTensor& st) +operator/(const Cmpt s, const SphericalTensor& st) { + #ifdef FULLDEBUG + if (mag(st.ii()) < VSMALL) + { + FatalErrorInFunction + << "Cmpt = " << s + << " is not divisible due to a zero element in SphericalTensor:" + << "SphericalTensor = " << st + << abort(FatalError); + } + #endif + return SphericalTensor(s/st.ii()); } +//- Division of a SphericalTensor by a Cmpt template -inline Cmpt magSqr(const SphericalTensor& st) +inline SphericalTensor +operator/(const SphericalTensor& st, const Cmpt s) { - return 3*magSqr(st.ii()); + #ifdef FULLDEBUG + if (mag(s) < VSMALL) + { + FatalErrorInFunction + << "SphericalTensor = " << st + << " is not divisible due to a zero value in Cmpt:" + << "Cmpt = " << s + << abort(FatalError); + } + #endif + + return SphericalTensor(st.ii()/s); } +//- Return the square of Frobenius norm of a SphericalTensor as a Cmpt +template +inline Cmpt magSqr(const SphericalTensor& st) +{ + return Cmpt(3*mag(st.ii()*st.ii())); +} + + +//- Return the max component of a SphericalTensor template inline Cmpt cmptMax(const SphericalTensor& st) { @@ -161,6 +195,7 @@ inline Cmpt cmptMax(const SphericalTensor& st) } +//- Return the min component of a SphericalTensor template inline Cmpt cmptMin(const SphericalTensor& st) { @@ -168,6 +203,7 @@ inline Cmpt cmptMin(const SphericalTensor& st) } +//- Return the sum of components of a SphericalTensor template inline Cmpt cmptSum(const SphericalTensor& st) { @@ -175,6 +211,7 @@ inline Cmpt cmptSum(const SphericalTensor& st) } +//- Return the arithmetic average of components of a SphericalTensor template inline Cmpt cmptAv(const SphericalTensor& st) { @@ -182,7 +219,7 @@ inline Cmpt cmptAv(const SphericalTensor& st) } -//- Return the trace of a spherical tensor +//- Return the trace of a SphericalTensor template inline Cmpt tr(const SphericalTensor& st) { @@ -190,7 +227,7 @@ inline Cmpt tr(const SphericalTensor& st) } -//- Return the spherical part of a spherical tensor, i.e. itself +//- Return the spherical part of a SphericalTensor, i.e. itself template inline SphericalTensor sph(const SphericalTensor& st) { @@ -198,7 +235,7 @@ inline SphericalTensor sph(const SphericalTensor& st) } -//- Return the determinant of a spherical tensor +//- Return the determinant of a SphericalTensor template inline Cmpt det(const SphericalTensor& st) { @@ -206,14 +243,26 @@ inline Cmpt det(const SphericalTensor& st) } -//- Return the inverse of a spherical tensor +//- Return the inverse of a SphericalTensor template inline SphericalTensor inv(const SphericalTensor& st) { - return SphericalTensor(1.0/st.ii()); + #ifdef FULLDEBUG + if (mag(st.ii()) < VSMALL) + { + FatalErrorInFunction + << "SphericalTensor is not invertible due to the zero determinant:" + << "det(SphericalTensor) = " << det(st) + << abort(FatalError); + } + #endif + + return SphericalTensor(1/st.ii()); } +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + template class outerProduct, Cmpt> { diff --git a/src/OpenFOAM/primitives/SphericalTensor/labelSphericalTensor/labelSphericalTensor.H b/src/OpenFOAM/primitives/SphericalTensor/labelSphericalTensor/labelSphericalTensor.H index 038e102aef..3ebc5ecff4 100644 --- a/src/OpenFOAM/primitives/SphericalTensor/labelSphericalTensor/labelSphericalTensor.H +++ b/src/OpenFOAM/primitives/SphericalTensor/labelSphericalTensor/labelSphericalTensor.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2019-2020 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -28,7 +28,10 @@ Typedef Foam::labelSphericalTensor Description - SphericalTensor of labels. + SphericalTensor of labels, i.e. SphericalTensor