From bb034f8ba62eb22ded8b5f980669ae26a6e830b2 Mon Sep 17 00:00:00 2001
From: henry <Henry Weller h.weller@opencfd.co.uk>
Date: Wed, 11 Jun 2008 13:32:42 +0100
Subject: [PATCH] Further tinkering and additional test of the eigenvaues of
symmTensors.
---
applications/test/tensor/tensorTest.C | 8 +++++++-
src/OpenFOAM/primitives/tensor/tensor.C | 18 +++++++++---------
2 files changed, 16 insertions(+), 10 deletions(-)
diff --git a/applications/test/tensor/tensorTest.C b/applications/test/tensor/tensorTest.C
index 01f00128acd..4530bc76579 100644
--- a/applications/test/tensor/tensorTest.C
+++ b/applications/test/tensor/tensorTest.C
@@ -38,6 +38,12 @@ int main()
<< (eigenVector(t6, e[2]) & t6) << e[2]*eigenVector(t6, e[2])
<< 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 "
@@ -49,7 +55,7 @@ int main()
<< transform(t1, st1) << endl;
vector v1(1, 2, 3);
-
+
Info<< sqr(v1) << endl;
Info<< symm(t7) << endl;
Info<< twoSymm(t7) << endl;
diff --git a/src/OpenFOAM/primitives/tensor/tensor.C b/src/OpenFOAM/primitives/tensor/tensor.C
index 4ad04ff2c4b..a3dfc2f7ea4 100644
--- a/src/OpenFOAM/primitives/tensor/tensor.C
+++ b/src/OpenFOAM/primitives/tensor/tensor.C
@@ -185,7 +185,7 @@ vector eigenValues(const tensor& t)
vector eigenVector(const tensor& t, const scalar lambda)
{
- if (lambda < SMALL)
+ if (mag(lambda) < SMALL)
{
return vector::zero;
}
@@ -292,22 +292,22 @@ vector eigenValues(const symmTensor& t)
+ t.xy()*t.xy()*t.zz() + t.xx()*t.yz()*t.yz();
// If there is a zero root
- if (mag(c) < SMALL)
+ if (mag(c) < 1.0e-100)
{
scalar disc = sqr(a) - 4*b;
- if (disc > 0)
+ if (disc >= -SMALL)
{
- scalar q = -0.5*(a + sign(a)*sqrt(disc));
+ scalar q = -0.5*sqrt(max(0.0, disc));
i = 0;
- ii = q;
- iii = b/q;
+ ii = -0.5*a + q;
+ iii = -0.5*a - q;
}
else
{
- FatalErrorIn("eigenValues(const symmTensor&)")
- << "zero and complex eigenvalues in symmTensor: " << t
+ FatalErrorIn("eigenValues(const tensor&)")
+ << "zero and complex eigenvalues in tensor: " << t
<< abort(FatalError);
}
}
@@ -380,7 +380,7 @@ vector eigenValues(const symmTensor& t)
vector eigenVector(const symmTensor& t, const scalar lambda)
{
- if (lambda < SMALL)
+ if (mag(lambda) < SMALL)
{
return vector::zero;
}
--
GitLab