diff --git a/applications/test/tensor/tensorTest.C b/applications/test/tensor/tensorTest.C index 01f00128acd876166c1b6e2f6a598600ddd365f5..4530bc7657919f8ceab844537c225a9931bd55a4 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 4ad04ff2c4beca1a02d312b7656794e1ed34c2e7..a3dfc2f7ea4f647753b0ac80056d3f07de3bc456 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; }