Reputation: 75
Working on a Finite element solver for the fun of it, I encountered a math problem where the area of a triangle has to be solved multiple times for different triangular elements on a mesh.
Note: This is all being done in Cython:
For the 2D case, the determinant formula allows us to do this by the least number of arithmetic operations.
A = 0.5|[x_1(y_2-y_3) + x_2(y_3-y_1) + x_3(y_1-y_2)]|
However this will be called once for each element and, in the case of adaptive meshing, for every iteration too. So finding the absolute value becomes a problem. Going with a function that does comparisons and if-else branching seems wrong. So i thought that bit manipulation could be a solution here maybe and went with this approach.
from libc.stdint cimport uint64_t
cdef union DoubleIntUnion:
double double_nr
uint64_t int_nr
cdef double fastABS(double number) nogil:
cdef DoubleIntUnion u_number
u_number.double_nr = number
u_number.int_nr &= <uint64_t>0x7FFFFFFFFFFFFFFF
return u_number.double_nr
this works. But is there any way to do it better in this case or any suggestions ?
Upvotes: 0
Views: 225
Reputation: 136495
If you just call std::abs
, which inlines floating point number binary and
to mask the sign bit, it would be challenging to beat it.
from libcpp.cmath cimport abs as std_abs
cdef double fastABS(double number) noexcept:
return std_abs(number)
Your solution with union cast involves loading a SIMD register into a general purpose register, masking the sign bit, and loading the general purpose register back into SIMD register. Such loads are not cheap and best avoided.
Benchmark your variants to get full understanding.
Upvotes: 3