UpsilonNumworks/Upsilon
1
0
Fork 0
mirror of https://github.com/UpsilonNumworks/Upsilon.git synced 2026-01-19 00:37:25 +01:00
Code Issues Packages Projects Releases Wiki Activity

[poincare] ComplexCartesian: fix approximation

Browse Source 
The assertions that the real and imaginary parts are real can be false
due to double lack of precision. We escape these cases by returning an
undefined approximation.
This commit is contained in:
Émilie Feral
2019-01-16 12:13:20 +01:00
committed by LeaNumworks
parent 8f0246b9d6
commit 6948de2ad3
1 changed files with 12 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files

12
poincare/src/complex_cartesian.cpp
View File

@@ -35,6 +35,18 @@ Complex<T> ComplexCartesianNode::templatedApproximate(Context& context, Preferen
assert(imagEvalution.type() == EvaluationNode<T>::Type::Complex);
std::complex<T> a = static_cast<Complex<T> &>(realEvaluation).stdComplex();
std::complex<T> b = static_cast<Complex<T> &>(imagEvalution).stdComplex();
if ((a.imag() != 0.0 && !std::isnan(a.imag())) || (b.imag() != 0.0 && !std::isnan(b.imag()))) {
/* a and b are supposed to be real (if they are not undefined). However,
* due to double precision limit, the approximation of the real part or the
* imaginary part can lead to complex values.
* For instance, let the real part be
* sqrt(2*sqrt(5E23+1)-1E12*sqrt(2)) ~ 1.1892E-6. Due to std::sqrt(2.0)
* unprecision, 2*sqrt(5E23+1)-1E12*sqrt(2) < 0 which leads to
* sqrt(2*sqrt(5E23+1)-1E12*sqrt(2)) being a complex number.
* In this case, we return an undefined complex because the approximation
* is very likely to be false. */
return Complex<T>::Undefined();
}
assert(a.imag() == 0.0 || std::isnan(a.imag()));
assert(b.imag() == 0.0 || std::isnan(b.imag()));
return Complex<T>(a.real(), b.real());