Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss
Reproducing a sequence from a binary tree...
143 views
Skip to first unread message
Chris M. Thomasson
May 24, 2020, 3:58:36 PM5/24/20
to
For fun, I reproduced the following sequence:
https://oeis.org/A059893
using a crude tree process:
_________________________
#include <iostream>
#include <cmath>
namespace ct_rifc
{
unsigned long store(unsigned long x)
{
unsigned long i = 0;
unsigned long v = 0;
unsigned long tree_val = 0;
for (i = 0; v != x; ++i)
{
unsigned long bit = ((x >> i) % 2);
v += bit * std::pow(2, i);
unsigned long child_left = tree_val * 2 + 1;
unsigned long child_right = tree_val * 2 + 2;
if (bit == 0)
{
// left
//std::cout << "store left\n";
tree_val = child_left;
}
else
{
// right
//std::cout << "store right\n";
tree_val = child_right;
}
//std::cout << "bit = " << bit << "\n";
//std::cout << "v = " << v << "\n";
//std::cout << "tree_val = " << tree_val << "\n\n";
}
return tree_val / 2;
}
unsigned long load(unsigned long x)
{
x = x * 2;
unsigned long x_tmp = x;
unsigned long v = 0;
unsigned long levels = 0;
for (levels = 0; x_tmp != 0; ++levels)
{
x_tmp = (x_tmp - 1) / 2;
}
for (unsigned long i = 0; x != 0; ++i)
{
unsigned long parent = (x - 1) / 2;
unsigned long child_left = parent * 2 + 1;
unsigned long child_right = parent * 2 + 2;
if (x == child_left)
{
// std::cout << "load left\n";
//v += 0 * std::pow(2, levels - i - 1); // mul by zero,
nothing.
}
else
{
//std::cout << "load right\n";
v += 1 * std::pow(2, levels - i - 1);
}
x = parent;
}
return v;
}
void test()
{
for (unsigned long i = 0; i < 256; ++i)
{
unsigned long plaintext = i;
unsigned long ciphertext = store(plaintext);
unsigned long decrypted = load(ciphertext);
//std::cout << "plaintext = " << plaintext << "\n";
//std::cout << "ciphertext = " << ciphertext << "\n";
//std::cout << "decrypted = " << decrypted << "\n";
std::cout << "(pt, ct, dt) = " <<
plaintext << ", " << ciphertext << ", " << decrypted <<
"\n";
//std::cout << "\n";
if (plaintext != decrypted)
{
std::cout << "FAILURE!\n";
break;
}
}
}
}
int main()
{
{
ct_rifc::test();
std::cout << "\n\nFIN!\n";
}
//std::cin.get();
return 0;
}
_________________________
Can you compile and run this thing? Fwiw, here is a sample of the output:
(pt, ct, dt) = 0, 0, 0
(pt, ct, dt) = 1, 1, 1
(pt, ct, dt) = 2, 2, 2
(pt, ct, dt) = 3, 3, 3
(pt, ct, dt) = 4, 4, 4
(pt, ct, dt) = 5, 6, 5
(pt, ct, dt) = 6, 5, 6
(pt, ct, dt) = 7, 7, 7
(pt, ct, dt) = 8, 8, 8
(pt, ct, dt) = 9, 12, 9
(pt, ct, dt) = 10, 10, 10
(pt, ct, dt) = 11, 14, 11
(pt, ct, dt) = 12, 9, 12
(pt, ct, dt) = 13, 13, 13
(pt, ct, dt) = 14, 11, 14
(pt, ct, dt) = 15, 15, 15
(pt, ct, dt) = 16, 16, 16
(pt, ct, dt) = 17, 24, 17
(pt, ct, dt) = 18, 20, 18
(pt, ct, dt) = 19, 28, 19
(pt, ct, dt) = 20, 18, 20
Notice the pattern in the ct (middle) column:
0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
https://oeis.org/A059893
using a crude tree process:
_________________________
#include <iostream>
#include <cmath>
namespace ct_rifc
{
unsigned long store(unsigned long x)
{
unsigned long i = 0;
unsigned long v = 0;
unsigned long tree_val = 0;
for (i = 0; v != x; ++i)
{
unsigned long bit = ((x >> i) % 2);
v += bit * std::pow(2, i);
unsigned long child_left = tree_val * 2 + 1;
unsigned long child_right = tree_val * 2 + 2;
if (bit == 0)
{
// left
//std::cout << "store left\n";
tree_val = child_left;
}
else
{
// right
//std::cout << "store right\n";
tree_val = child_right;
}
//std::cout << "bit = " << bit << "\n";
//std::cout << "v = " << v << "\n";
//std::cout << "tree_val = " << tree_val << "\n\n";
}
return tree_val / 2;
}
unsigned long load(unsigned long x)
{
x = x * 2;
unsigned long x_tmp = x;
unsigned long v = 0;
unsigned long levels = 0;
for (levels = 0; x_tmp != 0; ++levels)
{
x_tmp = (x_tmp - 1) / 2;
}
for (unsigned long i = 0; x != 0; ++i)
{
unsigned long parent = (x - 1) / 2;
unsigned long child_left = parent * 2 + 1;
unsigned long child_right = parent * 2 + 2;
if (x == child_left)
{
// std::cout << "load left\n";
//v += 0 * std::pow(2, levels - i - 1); // mul by zero,
nothing.
}
else
{
//std::cout << "load right\n";
v += 1 * std::pow(2, levels - i - 1);
}
x = parent;
}
return v;
}
void test()
{
for (unsigned long i = 0; i < 256; ++i)
{
unsigned long plaintext = i;
unsigned long ciphertext = store(plaintext);
unsigned long decrypted = load(ciphertext);
//std::cout << "plaintext = " << plaintext << "\n";
//std::cout << "ciphertext = " << ciphertext << "\n";
//std::cout << "decrypted = " << decrypted << "\n";
std::cout << "(pt, ct, dt) = " <<
plaintext << ", " << ciphertext << ", " << decrypted <<
"\n";
//std::cout << "\n";
if (plaintext != decrypted)
{
std::cout << "FAILURE!\n";
break;
}
}
}
}
int main()
{
{
ct_rifc::test();
std::cout << "\n\nFIN!\n";
}
//std::cin.get();
return 0;
}
_________________________
Can you compile and run this thing? Fwiw, here is a sample of the output:
(pt, ct, dt) = 0, 0, 0
(pt, ct, dt) = 1, 1, 1
(pt, ct, dt) = 2, 2, 2
(pt, ct, dt) = 3, 3, 3
(pt, ct, dt) = 4, 4, 4
(pt, ct, dt) = 5, 6, 5
(pt, ct, dt) = 6, 5, 6
(pt, ct, dt) = 7, 7, 7
(pt, ct, dt) = 8, 8, 8
(pt, ct, dt) = 9, 12, 9
(pt, ct, dt) = 10, 10, 10
(pt, ct, dt) = 11, 14, 11
(pt, ct, dt) = 12, 9, 12
(pt, ct, dt) = 13, 13, 13
(pt, ct, dt) = 14, 11, 14
(pt, ct, dt) = 15, 15, 15
(pt, ct, dt) = 16, 16, 16
(pt, ct, dt) = 17, 24, 17
(pt, ct, dt) = 18, 20, 18
(pt, ct, dt) = 19, 28, 19
(pt, ct, dt) = 20, 18, 20
Notice the pattern in the ct (middle) column:
0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
Sam
May 24, 2020, 10:47:35 PM5/24/20
to
Chris M. Thomasson writes:
> For fun, I reproduced the following sequence:
>
> https://oeis.org/A059893
>
> using a crude tree process:
> For fun, I reproduced the following sequence:
>
> https://oeis.org/A059893
>
> using a crude tree process:
> v += bit * std::pow(2, i);
Are you aware that the above line of code forces your computer to sacrifice
innocent electrons for the sake of performing the following task:
1) Convert i from integer to double
2) Compute the natural logarithm of 2.0, which is .69314718055994530941, and
then multiply it with the value calculated in step 1.
3) Calculate what is 2.71828182845904523536 raised to the power of the value
calculated in step 2.
4) Convert the variable "bit" from integer to double
5) Multiply the result of step 3 by step 4.
6) Convert the result of step 5 back to an integer value, perhaps rounding
it in an unexpected fashion. It's unlikely, but rules of C++ math
pendantically make it possible for 1*pow(2, 8), here, to result in 255
getting added to your "v". Surprise!
Was, all of that, your intent here?
What makes this so puzzling is that I see the previous line of code
demonstrated the knowledge and the correct use of the >> operator. So why
not, instead of all that, this can't be done as simply as:
v += bit * (1UL << i);
?
> Can you compile and run this thing?
> Notice the pattern in the ct (middle) column:
>
> 0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
Chris M. Thomasson
May 26, 2020, 7:53:05 PM5/26/20
to
On 5/24/2020 7:47 PM, Sam wrote:
> Chris M. Thomasson writes:
>
>> For fun, I reproduced the following sequence:
>>
>> https://oeis.org/A059893
>>
>> using a crude tree process:
>
>> v += bit * std::pow(2, i);
>
> Are you aware that the above line of code forces your computer to
> sacrifice innocent electrons for the sake of performing the following task:
>
> 1) Convert i from integer to double
>
> 2) Compute the natural logarithm of 2.0, which is .69314718055994530941,
> and then multiply it with the value calculated in step 1.
>
> 3) Calculate what is 2.71828182845904523536 raised to the power of the
> value calculated in step 2.
>
> 4) Convert the variable "bit" from integer to double
>
> 5) Multiply the result of step 3 by step 4.
>
> 6) Convert the result of step 5 back to an integer value, perhaps
> rounding it in an unexpected fashion. It's unlikely, but rules of C++
> math pendantically make it possible for 1*pow(2, 8), here, to result in
> 255 getting added to your "v". Surprise!
>
> Was, all of that, your intent here?
Basically. However, I want to extend this to handle higher bases. Right
> Chris M. Thomasson writes:
>
>> For fun, I reproduced the following sequence:
>>
>> https://oeis.org/A059893
>>
>> using a crude tree process:
>
>> v += bit * std::pow(2, i);
>
> Are you aware that the above line of code forces your computer to
> sacrifice innocent electrons for the sake of performing the following task:
>
> 1) Convert i from integer to double
>
> 2) Compute the natural logarithm of 2.0, which is .69314718055994530941,
> and then multiply it with the value calculated in step 1.
>
> 3) Calculate what is 2.71828182845904523536 raised to the power of the
> value calculated in step 2.
>
> 4) Convert the variable "bit" from integer to double
>
> 5) Multiply the result of step 3 by step 4.
>
> 6) Convert the result of step 5 back to an integer value, perhaps
> rounding it in an unexpected fashion. It's unlikely, but rules of C++
> math pendantically make it possible for 1*pow(2, 8), here, to result in
> 255 getting added to your "v". Surprise!
>
> Was, all of that, your intent here?
how its only a 2-ary tree. It can be extended to n-ary. The code is down
right crude and very basic for now.
> What makes this so puzzling is that I see the previous line of code
> demonstrated the knowledge and the correct use of the >> operator. So
> why not, instead of all that, this can't be done as simply as:
>
> v += bit * (1UL << i);
> ?
>> Can you compile and run this thing?
>
> Yes, I can. And I got the same results as you, so what's the next step?
experimental encryption based on complex numbers into integer land:
https://github.com/ChrisMThomasson/fractal_cipher/blob/master/RIFC/cpp/ct_rifc_sample.cpp
>> Notice the pattern in the ct (middle) column:
>>
>> 0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
>
> Yes, I noticed it. Happy to confirm it.
>
reproduce the sequence, but there it is! Interesting.
The Real Non Homosexual
May 26, 2020, 9:58:41 PM5/26/20
to
I don't have a C++ compiler on my computer anymore otherwise I would have also tried to run it. Uhh....I just hang out here so that I can troll Dan Cross.
Alf P. Steinbach
May 27, 2020, 6:26:48 AM5/27/20
to
namespace impl
{
constexpr inline auto intpow( const double base, const int
exponent )
-> double
{
double result = 1;
double weight = base;
for( int n = exponent; n != 0; weight *= weight ) {
if( is_odd( n ) ) {
result *= weight;
}
n /= 2;
}
return result;
}
} // namespace impl
/// @endcond
/// \brief Efficient *x* to the *n*'th power, when *n* is an integer.
///
/// \param base The *x* in “*x* to the *n*'th”.
/// \param exponent The *n* in “*x* to the *n*'th”.
///
/// Essentially this is Horner's rule adapted to calculating a
power, so that the
/// number of floating point multiplications is at worst O(log2(n)).
constexpr inline auto intpow( const double base, const int exponent )
-> double
{
return (0?0
: exponent == 0? 1.0
: exponent < 0? 1.0/impl::intpow( base, -exponent )
: impl::intpow( base, exponent )
);
}
<url:
https://github.com/alf-p-steinbach/cppx-core-language/blob/master/source/cppx-core-language/calc/floating-point-operations.hpp>
[snip]
- Alf
Bart
May 27, 2020, 7:21:41 AM5/27/20
to
and 1.0)
>> Can you compile and run this thing?
>
> Yes, I can. And I got the same results as you, so what's the next step?
>
>> Notice the pattern in the ct (middle) column:
>>
>> 0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
>
> Yes, I noticed it. Happy to confirm it.
(ie. pow) not shifts.
Ben Bacarisse
May 27, 2020, 8:02:00 AM5/27/20
to
Bart <b...@freeuk.com> writes:
> On 25/05/2020 03:47, Sam wrote:
>> Chris M. Thomasson writes:
>>
>>> For fun, I reproduced the following sequence:
>>>
>>> https://oeis.org/A059893
>>>
>>> using a crude tree process:
>>
>>> v += bit * std::pow(2, i);
>>
>> Are you aware that the above line of code forces your computer to sacrifice innocent electrons for the sake of performing the following task:
>>
>> 1) Convert i from integer to double
>>
>> 2) Compute the natural logarithm of 2.0, which is .69314718055994530941, and then multiply it with the value calculated in step 1.
>>
>> 3) Calculate what is 2.71828182845904523536 raised to the power of the value calculated in step 2.
>>
>> 4) Convert the variable "bit" from integer to double
>>
>> 5) Multiply the result of step 3 by step 4.
>>
>> 6) Convert the result of step 5 back to an integer value, perhaps rounding it in an unexpected fashion. It's unlikely, but rules of C++ math pendantically make it possible for 1*pow(2, 8), here, to result in 255 getting added to your "v". Surprise!
>>
>> Was, all of that, your intent here?
>>
>> What makes this so puzzling is that I see the previous line of code demonstrated the knowledge and the correct use of the >> operator. So why not, instead of all that, this can't be done as simply as:
>>
>> v += bit * (1UL << i);
>> ?
>
> When i = 0, don't (2<<0) and pow(2,0) give different results? (Namely
> 2 and 1.0)
Yes, but the suggestion was to use 1UL << i, not 2 << i. (I might use
> On 25/05/2020 03:47, Sam wrote:
>> Chris M. Thomasson writes:
>>
>>> For fun, I reproduced the following sequence:
>>>
>>> https://oeis.org/A059893
>>>
>>> using a crude tree process:
>>
>>> v += bit * std::pow(2, i);
>>
>> Are you aware that the above line of code forces your computer to sacrifice innocent electrons for the sake of performing the following task:
>>
>> 1) Convert i from integer to double
>>
>> 2) Compute the natural logarithm of 2.0, which is .69314718055994530941, and then multiply it with the value calculated in step 1.
>>
>> 3) Calculate what is 2.71828182845904523536 raised to the power of the value calculated in step 2.
>>
>> 4) Convert the variable "bit" from integer to double
>>
>> 5) Multiply the result of step 3 by step 4.
>>
>> 6) Convert the result of step 5 back to an integer value, perhaps rounding it in an unexpected fashion. It's unlikely, but rules of C++ math pendantically make it possible for 1*pow(2, 8), here, to result in 255 getting added to your "v". Surprise!
>>
>> Was, all of that, your intent here?
>>
>> What makes this so puzzling is that I see the previous line of code demonstrated the knowledge and the correct use of the >> operator. So why not, instead of all that, this can't be done as simply as:
>>
>> v += bit * (1UL << i);
>> ?
>
> When i = 0, don't (2<<0) and pow(2,0) give different results? (Namely
> 2 and 1.0)
1ull to get be sure to get the widest unsigned standard type.)
>>> Can you compile and run this thing?
>>
>> Yes, I can. And I got the same results as you, so what's the next step?
>>
>>> Notice the pattern in the ct (middle) column:
>>>
>>> 0, 1, 2, 3, 4, 6, 5, 7, 8, 12, 10, 14, 9, 13, 11, 15, ...
>>
>> Yes, I noticed it. Happy to confirm it.
>
> So can I, but in a completely different language. But I had to use **
> (ie. pow) not shifts.
provided you can represent the value in the chosen type, you can rely on
shifts to do the job.
--
Ben.
Bart
May 27, 2020, 10:22:35 AM5/27/20
to
On 27/05/2020 13:01, Ben Bacarisse wrote:
> Bart <b...@freeuk.com> writes:
>
>> On 25/05/2020 03:47, Sam wrote:
>>> v += bit * (1UL << i);
>>> ?
>>
>> When i = 0, don't (2<<0) and pow(2,0) give different results? (Namely
>> 2 and 1.0)
>
> Yes, but the suggestion was to use 1UL << i, not 2 << i. (I might use
> 1ull to get be sure to get the widest unsigned standard type.)
OK, that's probably why it didn't work!
> Bart <b...@freeuk.com> writes:
>
>> On 25/05/2020 03:47, Sam wrote:
>>> v += bit * (1UL << i);
>>> ?
>>
>> When i = 0, don't (2<<0) and pow(2,0) give different results? (Namely
>> 2 and 1.0)
>
> Yes, but the suggestion was to use 1UL << i, not 2 << i. (I might use
> 1ull to get be sure to get the widest unsigned standard type.)
Using 1<<i and i<<(level-i-1) worked as far as I could see.
(Calculations all using signed 64 bits in the language I used.)
Chris M. Thomasson
May 29, 2020, 4:44:18 PM5/29/20
to
:^)
Fwiw, here is an expensive root finding algorihtm for complex numbers:
__________________________
ct_complex
ct_root(
ct_complex const& z,
unsigned int r,
unsigned int b
) {
double angle_base = std::arg(z) / b;
double angle_step = CT_PI2 / b;
double length = std::abs(z);
double angle = angle_base + angle_step * r;
double radius = std::pow(length, 1.0 / b);
ct_complex root(
std::cos(angle) * radius,
std::sin(angle) * radius
);
return root;
}
__________________________
r is the root to find in the base b.
Alf P. Steinbach
May 29, 2020, 6:02:20 PM5/29/20
to
as shown (the code you snipped).
You can reduce the number of expensive calls to `cos` and `sin`
* by recognizing special bases like 2, 3 and 4 (possibly you're only
using a very few bases, in which case this covers all of it!),
* for the more general use cases, by returning a root generator object
rather than directly a complex, where the generator object internally
retains the first non-1 b'th root of 1, and calculates the various angle
roots as powers of that, times root of norm.
For the last point the idea is that where code now calls ct_root
repeatedly to obtain all the various angle roots, with b calls to `pow`,
`cos` and `sin, it can instead call ct_root_gen with just 1 call to
`pow`, `cos` and `sin`, to get a generator which it then can call
repeatedly to more fastishly get the b roots.
- Alf
Alf P. Steinbach
May 29, 2020, 6:07:44 PM5/29/20
to
is a corresponding way to calculate real roots fastishly.
- - -
This however holds:
Chris M. Thomasson
Jun 2, 2020, 2:30:54 AM6/2/20
to
The Real Non Homosexual
Jun 2, 2020, 5:24:04 PM6/2/20
to
e can either mean exponent or the euler constant.
Chris M. Thomasson
Jun 2, 2020, 6:18:21 PM6/2/20
to
On 6/2/2020 2:23 PM, The Real Non Homosexual wrote:
> e can either mean exponent or the euler constant.
>
The latter.
> e can either mean exponent or the euler constant.
>
Alf P. Steinbach
Jun 3, 2020, 12:13:39 AM6/3/20
to
Perhaps you though I was just spouting nonsense. To cover that
possibility, here's a trivial implementation of a root generator:
-----------------------------------------------------------------------------
#ifndef _USE_MATH_DEFINES
# error "_USE_MATH_DEFINES not defined. Needed for Posix constant M_PI."
#endif
#include <complex>
#include <cmath>
using Complex = std::complex<double>;
namespace assumed {
using ct_complex = Complex;
const double CT_PI2 = 2*M_PI;
} // assumed
namespace original {
using assumed::ct_complex, assumed::CT_PI2;
ct_complex
ct_root(
ct_complex const& z,
unsigned int r,
unsigned int b
) {
double angle_base = std::arg(z) / b;
double angle_step = CT_PI2 / b;
double length = std::abs(z);
double angle = angle_base + angle_step * r;
double radius = std::pow(length, 1.0 / b);
ct_complex root(
std::cos(angle) * radius,
std::sin(angle) * radius
);
return root;
}
namespace fastish {
struct Roots
{
Complex m_primary;
Complex m_delta;
Complex m_current;
int m_index;
int m_base;
public:
// Further optimization by merging the two calls to ct_root here.
Roots( const Complex& z, const int base ):
m_primary( original::ct_root( z, 0, base ) ),
m_delta( original::ct_root( 1.0, 1, base ) ),
m_current( 1, 0 ),
m_index( 0 ),
m_base( base )
{}
auto more() const -> bool { return m_index < m_base; }
auto current()
-> Complex
{ return m_current*m_primary; }
void advance()
{
m_current *= m_delta;
++m_index;
}
};
} // fastish
#include <iostream>
using namespace std::literals;
using std::pow,
std::cout, std::endl;
auto main() -> int
{
const int n = 5;
const assumed::ct_complex z = pow( Complex( 3, 0.5 ), n );
for( int i = 0; i < n; ++i ) {
cout << original::ct_root( z, i, n ) << endl;
}
cout << endl;
for( fastish::Roots roots( z, n ); roots.more(); roots.advance() ) {
cout << roots.current() << endl;
}
}
-----------------------------------------------------------------------------
And yes, the two loops produce the same output.
- Alf
Chris M. Thomasson
Jun 3, 2020, 2:39:39 AM6/3/20
to
On 6/2/2020 9:13 PM, Alf P. Steinbach wrote:
> On 02.06.2020 08:30, Chris M. Thomasson wrote:
[...]
> On 02.06.2020 08:30, Chris M. Thomasson wrote:
> cout << endl;
> for( fastish::Roots roots( z, n ); roots.more(); roots.advance() ) {
> cout << roots.current() << endl;
> }
> }
> -----------------------------------------------------------------------------
>
>
>
> And yes, the two loops produce the same output.
I am not trying to find all of the roots, just a single targeted one.
> for( fastish::Roots roots( z, n ); roots.more(); roots.advance() ) {
> cout << roots.current() << endl;
> }
> }
> -----------------------------------------------------------------------------
>
>
>
> And yes, the two loops produce the same output.
Say root 6 of base 10. Say, there are ten roots in any complex number z
wrt z^10, I only need the 6'th root.
The sixth root of complex number, say -1+2i in base ten is:
6 digit precision:
(-0.730018,-0.801057)
Raise this to the tenth power and you get -1+2i.
The seventh root:
(-0.119748,-1.07716)
Raised to the tenth power gets -1+2i.
Chris M. Thomasson
Jun 3, 2020, 2:43:44 AM6/3/20
to
Alf P. Steinbach
Jun 3, 2020, 7:37:42 AM6/3/20
to
Up front disclaimer: I haven't looked at the video.
One optimization opportunity is that you call `std::abs`, which involves
a costly square root, only to feed that later as argument to `pow`.
Instead call `std::norm`, which produces the square, and adjust the
exponent argument to `pow`. Like,
const double radius = std::pow( std::norm(z), (1.0/b)/2 );
Another possible optimization for calculation of a single root of
general z is that e.g. sin(x) can be computed as sqrt(1 - cos(x)^2),
plus relevant sign adjustment. But I'm not sure if it's faster. After
all the trig functions are done in hardware now.
And maybe you could compute the integral real root faster via
Newton/Raphson approximation rather than via general `pow`, just as
`sqrt(x)` is probably faster than `pow(x, 0.5)`. I haven't tried though,
and it precludes the optimization of using `norm`, and again, since the
basic math ops are done in hardware now it's questionable whether a
software Newton/Raphson loop would be faster than hardware `pow`. Still
one could try, and measure.
Cheers!,
- Alf
Alf P. Steinbach
Jun 3, 2020, 3:12:28 PM6/3/20
to
myself again, for the second time this thread.
> and again, since the
> basic math ops are done in hardware now it's questionable whether a
> software Newton/Raphson loop would be faster than hardware `pow`. Still
> one could try, and measure.
Alf P. Steinbach
Jun 4, 2020, 4:16:57 PM6/4/20
to
the ones mentioned above, is to use `norm` instead of `abs`.
- - -
For the express-sin-vis-cos possibility I measured these functions:
struct Stdlib
{
static auto func( const double v )
-> double
{ return sin( v )*cos( v ); }
};
struct Pythagoras
{
static auto func( const double v )
-> double
{
const double cos_value = cos( v );
const double abs_sin_value = sqrt( 1 - cos_value*cos_value );
return (sgn( v )*abs_sin_value)*cos_value;
}
};
... with the following results with 64-bit compilers in Windows 10,
where a.exe is produced by MinG g++ 9.2, and b.exe by MSVC 2019:
[H:\forums\clc++\058 complex roots]
> a
Stdlib : 0.315015 in 4.17709e-10 seconds (average over
262144 iterations).
Pythagoras : 0.315015 in 4.18091e-10 seconds (average over
262144 iterations).
[H:\forums\clc++\058 complex roots]
> b
Stdlib : 0.315015 in 4.17709e-10 seconds (average over
262144 iterations).
Pythagoras : 0.315015 in 2.12463e-08 seconds (average over
16384 iterations).
With g++ I used option `-O3`, and with MSVC I used option `/O2`.
So with MinGW g++ it's a tie, and with MSVC the standard library wins
handily.
- - -
For the compute-root-via-Newton/Raphson idea I measured these functions:
struct Stdlib
{
static auto func( const double power, const int base )
-> double
{ return pow( power, 1.0/base ); }
};
struct Newton
{
static auto func( const double power, const int base )
-> double
{
// x^base = power => f(x) = x^base - power, f'(x) =
base*x^(base - 1)
double xp = power; // previous x
for( ;; ) {
const double p = intpow( xp, base - 1 );
const double fx = p*xp - power;
const double fdx = base*p;
const double new_guess = xp - fx/fdx;
if( abs( new_guess - xp )/xp < 1e-14 ) {
return new_guess;
}
xp = new_guess;
}
}
};
... with results
[H:\forums\clc++\058 complex roots]
> a
Stdlib : 5 in 4.17709e-10 seconds (average over 262144
iterations).
Newton : 5 in 7.95996e-07 seconds (average over 2048
iterations).
[H:\forums\clc++\058 complex roots]
> b
Stdlib : 5 in 4.20666e-10 seconds (average over 1048576
iterations).
Newton : 5 in 4.17709e-10 seconds (average over 262144
iterations).
So for these functions with MinGW g++ the standard library wins handily,
while with MSVC it's a tie, opposite of the previous function pair.
- Alf
Chris M. Thomasson
Jun 4, 2020, 10:07:04 PM6/4/20
to
On 6/4/2020 1:16 PM, Alf P. Steinbach wrote:
> On 03.06.2020 21:12, Alf P. Steinbach wrote:
>> On 03.06.2020 13:37, Alf P. Steinbach wrote:
>>> On 03.06.2020 08:43, Chris M. Thomasson wrote:
>>>> On 6/2/2020 11:39 PM, Chris M. Thomasson wrote:
>>>>> On 6/2/2020 9:13 PM, Alf P. Steinbach wrote:
>>>>>> On 02.06.2020 08:30, Chris M. Thomasson wrote:
[...]
Wrt fast, check this out wrt obtaining the hypotenuse:
> On 03.06.2020 21:12, Alf P. Steinbach wrote:
>> On 03.06.2020 13:37, Alf P. Steinbach wrote:
>>> On 03.06.2020 08:43, Chris M. Thomasson wrote:
>>>> On 6/2/2020 11:39 PM, Chris M. Thomasson wrote:
>>>>> On 6/2/2020 9:13 PM, Alf P. Steinbach wrote:
>>>>>> On 02.06.2020 08:30, Chris M. Thomasson wrote:
[...]
https://forums.parallax.com/discussion/147522/dog-leg-hypotenuse-approximation
Using that for the roots might be good for generating a faster rendering
of a fractal, however, it might not be such a good thing to use in my
fractal root storage experiment thing here:
https://github.com/ChrisMThomasson/fractal_cipher/blob/master/RIFC/cpp/ct_rifc_sample.cpp
Btw, can you run this?
It really should be using an arbitrary precision floating point package.
Or perhaps, rationals might work as well.
Ralf Goertz
Jun 5, 2020, 4:01:12 AM6/5/20
to
Am Thu, 4 Jun 2020 22:16:41 +0200
schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
> For the express-sin-vis-cos possibility I measured these functions:
>
> struct Stdlib
> {
> static auto func( const double v )
> -> double
> { return sin( v )*cos( v ); }
> };
>
> struct Pythagoras
> {
> static auto func( const double v )
> -> double
> {
> const double cos_value = cos( v );
> const double abs_sin_value = sqrt( 1 -
> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
> }
> };
Did you know that with g++ you have `void sincos(double x, double *sin,
double *cos)` available? It is much faster than computing sin and cos
separately.
>
> For the compute-root-via-Newton/Raphson idea I measured these
> functions:
>
> struct Stdlib
> {
> static auto func( const double power, const int base )
> -> double
> { return pow( power, 1.0/base ); }
> };
I find your naming convention a bit confusing. In b^e=p I know b as the
base, e as the exponent and p the power.
schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
> For the express-sin-vis-cos possibility I measured these functions:
>
> struct Stdlib
> {
> static auto func( const double v )
> -> double
> { return sin( v )*cos( v ); }
> };
>
> struct Pythagoras
> {
> static auto func( const double v )
> -> double
> {
> const double cos_value = cos( v );
> const double abs_sin_value = sqrt( 1 -
> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
> }
> };
double *cos)` available? It is much faster than computing sin and cos
separately.
>
> For the compute-root-via-Newton/Raphson idea I measured these
> functions:
>
> struct Stdlib
> {
> static auto func( const double power, const int base )
> -> double
> { return pow( power, 1.0/base ); }
> };
base, e as the exponent and p the power.
Alf P. Steinbach
Jun 5, 2020, 4:52:10 AM6/5/20
to
On 05.06.2020 10:01, Ralf Goertz wrote:
> Am Thu, 4 Jun 2020 22:16:41 +0200
> schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
>
>> For the express-sin-vis-cos possibility I measured these functions:
>>
>> struct Stdlib
>> {
>> static auto func( const double v )
>> -> double
>> { return sin( v )*cos( v ); }
>> };
>>
>> struct Pythagoras
>> {
>> static auto func( const double v )
>> -> double
>> {
>> const double cos_value = cos( v );
>> const double abs_sin_value = sqrt( 1 -
>> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
>> }
>> };
>
> Did you know that with g++ you have `void sincos(double x, double *sin,
> double *cos)` available? It is much faster than computing sin and cos
> separately.
Thanks, I didn't know, but should have guessed.
> Am Thu, 4 Jun 2020 22:16:41 +0200
> schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
>
>> For the express-sin-vis-cos possibility I measured these functions:
>>
>> struct Stdlib
>> {
>> static auto func( const double v )
>> -> double
>> { return sin( v )*cos( v ); }
>> };
>>
>> struct Pythagoras
>> {
>> static auto func( const double v )
>> -> double
>> {
>> const double cos_value = cos( v );
>> const double abs_sin_value = sqrt( 1 -
>> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
>> }
>> };
>
> Did you know that with g++ you have `void sincos(double x, double *sin,
> double *cos)` available? It is much faster than computing sin and cos
> separately.
There is also an x86 instruction, fsincos, but an SO answer claims that
"As CPUs become faster, old hardware trig instructions like fsin have
not kept pace. With some CPUs, a software function, using a polynomial
approximation for sine or another trig function, is now faster than a
hardware instruction. In short, fsincos is too slow."
<url: https://stackoverflow.com/a/24470801>
>> For the compute-root-via-Newton/Raphson idea I measured these
>> functions:
>>
>> struct Stdlib
>> {
>> static auto func( const double power, const int base )
>> -> double
>> { return pow( power, 1.0/base ); }
>> };
>
> I find your naming convention a bit confusing. In b^e=p I know b as the
> base, e as the exponent and p the power.
the root to find in the base b".
I believe the standard term for the n in n'th root is "degree". The OP's
use of "base" for the n in a root can be viewed as consistent with use
of "base" for the fixed parameter of a logarithm function; a view where
"base" refers to the fixed parameter. However, the standard use of
"base" in this context refers to the base of an exponentiation, e.g. x
is the "base" in both x^n=p and n=logₓp.
The function name `func` supports general templated measuring code.
- Alf
David Brown
Jun 5, 2020, 5:32:38 AM6/5/20
to
On 05/06/2020 10:01, Ralf Goertz wrote:
> Am Thu, 4 Jun 2020 22:16:41 +0200
> schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
>
>> For the express-sin-vis-cos possibility I measured these functions:
>>
>> struct Stdlib
>> {
>> static auto func( const double v )
>> -> double
>> { return sin( v )*cos( v ); }
>> };
>>
>> struct Pythagoras
>> {
>> static auto func( const double v )
>> -> double
>> {
>> const double cos_value = cos( v );
>> const double abs_sin_value = sqrt( 1 -
>> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
>> }
>> };
>
> Did you know that with g++ you have `void sincos(double x, double *sin,
> double *cos)` available? It is much faster than computing sin and cos
> separately.
>
No, "g++" does not have such a function - it is part of glibc, not the
> Am Thu, 4 Jun 2020 22:16:41 +0200
> schrieb "Alf P. Steinbach" <alf.p.stein...@gmail.com>:
>
>> For the express-sin-vis-cos possibility I measured these functions:
>>
>> struct Stdlib
>> {
>> static auto func( const double v )
>> -> double
>> { return sin( v )*cos( v ); }
>> };
>>
>> struct Pythagoras
>> {
>> static auto func( const double v )
>> -> double
>> {
>> const double cos_value = cos( v );
>> const double abs_sin_value = sqrt( 1 -
>> cos_value*cos_value ); return (sgn( v )*abs_sin_value)*cos_value;
>> }
>> };
>
> Did you know that with g++ you have `void sincos(double x, double *sin,
> double *cos)` available? It is much faster than computing sin and cos
> separately.
>
compiler, and it is a non-standard extension in <math.h>. That means
you need to define the macro "_GNU_SOURCE" before including <math.h> to
get access to it.
And it is /not/ faster than calculating sin and cos separately, because
gcc is smart enough to combine calls to sin and cos with the same
operand to a call to sincos.
You can see it here (with test code that is suitable for C and C++, in
case there was a difference):
#define _GNU_SOURCE
#include <math.h>
typedef struct SC {
double s;
double c;
} SC;
SC test1(double x) {
SC sc;
sc.s = sin(x);
sc.c = cos(x);
return sc;
}
SC test2(double x) {
SC sc;
sincos(x, &sc.s, &sc.c);
return sc;
}
With -O2, gcc (via <https://godbolt.org>) gives the same code for both
functions:
test1:
sub rsp, 24
mov rsi, rsp
lea rdi, [rsp+8]
call sincos
movsd xmm0, QWORD PTR [rsp+8]
movsd xmm1, QWORD PTR [rsp]
add rsp, 24
ret
(clang and Intel icc do the same optimisation.)
"sincos" is like the "div" function in <stdlib.h> - a relic from the
past when compilers were weaker, that was useful at the time but has
lost its purpose. (Unless, of course, you feel its use makes your code
clearer in some way - in which case, great.)
Ralf Goertz
Jun 5, 2020, 5:51:14 AM6/5/20
to
Am Fri, 5 Jun 2020 11:32:28 +0200
schrieb David Brown <david...@hesbynett.no>:
> On 05/06/2020 10:01, Ralf Goertz wrote:
> >
> > Did you know that with g++ you have `void sincos(double x, double
> > *sin, double *cos)` available? It is much faster than computing sin
> > and cos separately.
> >
>
> No, "g++" does not have such a function - it is part of glibc, not
> the compiler, and it is a non-standard extension in <math.h>. That
> means you need to define the macro "_GNU_SOURCE" before including
> <math.h> to get access to it.
Well, that's why I wrote that the function was available with g++. As
Alf mentioned MinGW he'd be using glibc in that case.
> And it is /not/ faster than calculating sin and cos separately,
> because gcc is smart enough to combine calls to sin and cos with the
> same operand to a call to sincos.
That might very well be. I remember to have played around with this when
it still mattered.
schrieb David Brown <david...@hesbynett.no>:
> On 05/06/2020 10:01, Ralf Goertz wrote:
> >
> > Did you know that with g++ you have `void sincos(double x, double
> > *sin, double *cos)` available? It is much faster than computing sin
> > and cos separately.
> >
>
> No, "g++" does not have such a function - it is part of glibc, not
> the compiler, and it is a non-standard extension in <math.h>. That
> means you need to define the macro "_GNU_SOURCE" before including
> <math.h> to get access to it.
Alf mentioned MinGW he'd be using glibc in that case.
> And it is /not/ faster than calculating sin and cos separately,
> because gcc is smart enough to combine calls to sin and cos with the
> same operand to a call to sincos.
it still mattered.
David Brown
Jun 5, 2020, 6:35:39 AM6/5/20
to
On 05/06/2020 11:51, Ralf Goertz wrote:
> Am Fri, 5 Jun 2020 11:32:28 +0200
> schrieb David Brown <david...@hesbynett.no>:
>
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>
>>> Did you know that with g++ you have `void sincos(double x, double
>>> *sin, double *cos)` available? It is much faster than computing sin
>>> and cos separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not
>> the compiler, and it is a non-standard extension in <math.h>. That
>> means you need to define the macro "_GNU_SOURCE" before including
>> <math.h> to get access to it.
>
> Well, that's why I wrote that the function was available with g++. As
> Alf mentioned MinGW he'd be using glibc in that case.
For his case, perhaps. In almost all /my/ uses of g++, I don't have
> Am Fri, 5 Jun 2020 11:32:28 +0200
> schrieb David Brown <david...@hesbynett.no>:
>
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>
>>> Did you know that with g++ you have `void sincos(double x, double
>>> *sin, double *cos)` available? It is much faster than computing sin
>>> and cos separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not
>> the compiler, and it is a non-standard extension in <math.h>. That
>> means you need to define the macro "_GNU_SOURCE" before including
>> <math.h> to get access to it.
>
> Well, that's why I wrote that the function was available with g++. As
> Alf mentioned MinGW he'd be using glibc in that case.
glibc. It is not always easy to know if someone is talking about a very
specific implementation, or in more general terms.
>
>> And it is /not/ faster than calculating sin and cos separately,
>> because gcc is smart enough to combine calls to sin and cos with the
>> same operand to a call to sincos.
>
> That might very well be. I remember to have played around with this when
> it still mattered.
>
are usually best aiming at writing the code in a clear and simple way,
with high-level algorithmic efficiency in mind, and letting the compiler
figure out the low-level efficiency details.
And if you want fast floating point in gcc, "-ffast-math" can make a
/lot/ of difference to some code.
Manfred
Jun 5, 2020, 7:50:28 AM6/5/20
to
On 6/5/2020 11:32 AM, David Brown wrote:
> On 05/06/2020 10:01, Ralf Goertz wrote:
>>
>> Did you know that with g++ you have `void sincos(double x, double *sin,
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>>
>
> No, "g++" does not have such a function - it is part of glibc, not the
> compiler, and it is a non-standard extension in <math.h>. That means
> you need to define the macro "_GNU_SOURCE" before including <math.h> to
> get access to it.
>
> And it is /not/ faster than calculating sin and cos separately, because
> gcc is smart enough to combine calls to sin and cos with the same
> operand to a call to sincos.
Yes, it /is/ faster than calculating sin and cos separately, that is why
> On 05/06/2020 10:01, Ralf Goertz wrote:
>>
>> Did you know that with g++ you have `void sincos(double x, double *sin,
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>>
>
> No, "g++" does not have such a function - it is part of glibc, not the
> compiler, and it is a non-standard extension in <math.h>. That means
> you need to define the macro "_GNU_SOURCE" before including <math.h> to
> get access to it.
>
> And it is /not/ faster than calculating sin and cos separately, because
> gcc is smart enough to combine calls to sin and cos with the same
> operand to a call to sincos.
the gcc optimizer replaces separate calls to sin and cos into a single
sincos call.
Which, by the way, only works when the optimizer can detect that such
calls are effectively part of the same operation.
in pairs is very common, so it is a very useful function.
David Brown
Jun 5, 2020, 8:08:14 AM6/5/20
to
On 05/06/2020 13:50, Manfred wrote:
> On 6/5/2020 11:32 AM, David Brown wrote:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not the
>> compiler, and it is a non-standard extension in <math.h>. That means
>> you need to define the macro "_GNU_SOURCE" before including <math.h>
>> to get access to it.
>>
>> And it is /not/ faster than calculating sin and cos separately,
>> because gcc is smart enough to combine calls to sin and cos with the
>> same operand to a call to sincos.
>
> Yes, it /is/ faster than calculating sin and cos separately, that is why
> the gcc optimizer replaces separate calls to sin and cos into a single
> sincos call.
I don't care how the compiler and/or library implements sin, cos or
> On 6/5/2020 11:32 AM, David Brown wrote:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not the
>> compiler, and it is a non-standard extension in <math.h>. That means
>> you need to define the macro "_GNU_SOURCE" before including <math.h>
>> to get access to it.
>>
>> And it is /not/ faster than calculating sin and cos separately,
>> because gcc is smart enough to combine calls to sin and cos with the
>> same operand to a call to sincos.
>
> Yes, it /is/ faster than calculating sin and cos separately, that is why
> the gcc optimizer replaces separate calls to sin and cos into a single
> sincos call.
sincos. It could use polynomials, or FPU instructions, or a Cordic
algorithm, or a circle and a ruler. It doesn't matter if the compiler
and/or library calculates sincos in a single cycle, or if it calculates
sin first then cos. The speed of the user code is /identical/ when you
write the sin and cos separately, or when you write a call to sincos,
because the code generated by the compiler is identical.
> Which, by the way, only works when the optimizer can detect that such
> calls are effectively part of the same operation.
>
expect sincos to be an alternative.
>>
>> "sincos" is like the "div" function in <stdlib.h> - a relic from the
>> past when compilers were weaker, that was useful at the time but has
>> lost its purpose. (Unless, of course, you feel its use makes your
>> code clearer in some way - in which case, great.)
>>
>>
>
> I wouldn't call it a relic - in fact the need to calculate sin and cos
> in pairs is very common, so it is a very useful function.
sin and cos efficiently in pairs, and will typically lead to uglier code.
From Alf's code (with more conventional syntax):
#include <math.h>
double foo1(double v) {
return sin(v) * cos(v);
}
Using sincos:
#define _GNU_SOURCE
#include <math.h>
double foo2(double v) {
double s, c;
sincos(v, &s, &c);
return s * c;
}
Improvements in the compiler mean you no longer need to write awkward
and convoluted code like foo2 - you can write it in a neater manner like
foo1 and get /identical/ performance. That is why sincos is a relic.
(Of course the library function sincos itself is still important to get
the efficiency - but that is an implementation detail. It can be hidden
from the user - it could be called __sin_and_cos, or implemented as
inline assembly or a builtin function. None of that matters to the user.)
Melzzzzz
Jun 5, 2020, 8:20:44 AM6/5/20
to
On 2020-06-05, David Brown <david...@hesbynett.no> wrote:
> On 05/06/2020 13:50, Manfred wrote:
>> On 6/5/2020 11:32 AM, David Brown wrote:
>>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>>
>>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>>> double *cos)` available? It is much faster than computing sin and cos
>>>> separately.
>>>>
>>>
>>> No, "g++" does not have such a function - it is part of glibc, not the
>>> compiler, and it is a non-standard extension in <math.h>. That means
>>> you need to define the macro "_GNU_SOURCE" before including <math.h>
>>> to get access to it.
>>>
>>> And it is /not/ faster than calculating sin and cos separately,
>>> because gcc is smart enough to combine calls to sin and cos with the
>>> same operand to a call to sincos.
>>
>> Yes, it /is/ faster than calculating sin and cos separately, that is why
>> the gcc optimizer replaces separate calls to sin and cos into a single
>> sincos call.
>
> I don't care how the compiler and/or library implements sin, cos or
> sincos. It could use polynomials, or FPU instructions, or a Cordic
> algorithm, or a circle and a ruler. It doesn't matter if the compiler
> and/or library calculates sincos in a single cycle, or if it calculates
> sin first then cos. The speed of the user code is /identical/ when you
> write the sin and cos separately, or when you write a call to sincos,
> because the code generated by the compiler is identical.
So compiler uses sincos when appropriate?
> On 05/06/2020 13:50, Manfred wrote:
>> On 6/5/2020 11:32 AM, David Brown wrote:
>>> On 05/06/2020 10:01, Ralf Goertz wrote:
>>>>
>>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>>> double *cos)` available? It is much faster than computing sin and cos
>>>> separately.
>>>>
>>>
>>> No, "g++" does not have such a function - it is part of glibc, not the
>>> compiler, and it is a non-standard extension in <math.h>. That means
>>> you need to define the macro "_GNU_SOURCE" before including <math.h>
>>> to get access to it.
>>>
>>> And it is /not/ faster than calculating sin and cos separately,
>>> because gcc is smart enough to combine calls to sin and cos with the
>>> same operand to a call to sincos.
>>
>> Yes, it /is/ faster than calculating sin and cos separately, that is why
>> the gcc optimizer replaces separate calls to sin and cos into a single
>> sincos call.
>
> I don't care how the compiler and/or library implements sin, cos or
> sincos. It could use polynomials, or FPU instructions, or a Cordic
> algorithm, or a circle and a ruler. It doesn't matter if the compiler
> and/or library calculates sincos in a single cycle, or if it calculates
> sin first then cos. The speed of the user code is /identical/ when you
> write the sin and cos separately, or when you write a call to sincos,
> because the code generated by the compiler is identical.
> foo1 and get /identical/ performance. That is why sincos is a relic.
calculating only one.
>
> (Of course the library function sincos itself is still important to get
> the efficiency - but that is an implementation detail. It can be hidden
> from the user - it could be called __sin_and_cos, or implemented as
> inline assembly or a builtin function. None of that matters to the user.)
expect from compiler to figure out?
--
current job title: senior software engineer
skills: c++,c,rust,go,nim,haskell...
press any key to continue or any other to quit...
U ničemu ja ne uživam kao u svom statusu INVALIDA -- Zli Zec
Svi smo svedoci - oko 3 godine intenzivne propagande je dovoljno da jedan narod poludi -- Zli Zec
Na divljem zapadu i nije bilo tako puno nasilja, upravo zato jer su svi
bili naoruzani. -- Mladen Gogala
David Brown
Jun 5, 2020, 8:55:54 AM6/5/20
to
>> foo1 and get /identical/ performance. That is why sincos is a relic.
>
> calculating both sincos is much cheaper as you get both results when
> calculating only one.
>
>>
>> (Of course the library function sincos itself is still important to get
>> the efficiency - but that is an implementation detail. It can be hidden
>> from the user - it could be called __sin_and_cos, or implemented as
>> inline assembly or a builtin function. None of that matters to the user.)
>
> So instead of using sincos when needed we call both sin and cos and
> expect from compiler to figure out?
>
figure out the low-level details.
Manfred
Jun 5, 2020, 9:11:51 AM6/5/20
to
On 6/5/2020 10:52 AM, Alf P. Steinbach wrote:
> On 05.06.2020 10:01, Ralf Goertz wrote:
<...>
> On 05.06.2020 10:01, Ralf Goertz wrote:
>> Did you know that with g++ you have `void sincos(double x, double *sin,
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>
> Thanks, I didn't know, but should have guessed.
>
> There is also an x86 instruction, fsincos, but an SO answer claims that
>
> "As CPUs become faster, old hardware trig instructions like fsin have
> not kept pace. With some CPUs, a software function, using a polynomial
> approximation for sine or another trig function, is now faster than a
> hardware instruction. In short, fsincos is too slow."
> <url: https://stackoverflow.com/a/24470801>
>
There's more to that: "All 64-bit SIMD integer instructions use MMX
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>
> Thanks, I didn't know, but should have guessed.
>
> There is also an x86 instruction, fsincos, but an SO answer claims that
>
> "As CPUs become faster, old hardware trig instructions like fsin have
> not kept pace. With some CPUs, a software function, using a polynomial
> approximation for sine or another trig function, is now faster than a
> hardware instruction. In short, fsincos is too slow."
> <url: https://stackoverflow.com/a/24470801>
>
registers, which share register state with the x87 floating point
stack." (from e.g. the Intel® 64 and IA-32 Architectures Optimization
Reference Manual)
This means that in order to use any of the x87 instructions with SIMD
code (most common on x64) special care must be taken by the compiler.
Speed in floating point mode is a tricky subject, for example the x87
added the hardware parallel pipeline for FPU instructions, while SIMD
instructions follow the main CPU pipeline, so it really depends on
details of the code generated by the compiler.
As far as I know on x64 the x87 instructions have lost interest mostly
because of the spread of SIMD code, part of which is incompatible with
(or at least unpractical to mix with) x87 instructions.
However, most algorithms that calculate trigonometric functions do
produce the combined sin/cos pair much faster than the separate ones.
Melzzzzz
Jun 5, 2020, 9:18:19 AM6/5/20
to
On 2020-06-05, Manfred <non...@add.invalid> wrote:
> On 6/5/2020 10:52 AM, Alf P. Steinbach wrote:
>> On 05.06.2020 10:01, Ralf Goertz wrote:
><...>
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>
>> Thanks, I didn't know, but should have guessed.
>>
>> There is also an x86 instruction, fsincos, but an SO answer claims that
>>
>> "As CPUs become faster, old hardware trig instructions like fsin have
>> not kept pace. With some CPUs, a software function, using a polynomial
>> approximation for sine or another trig function, is now faster than a
>> hardware instruction. In short, fsincos is too slow."
>> <url: https://stackoverflow.com/a/24470801>
>>
>
> There's more to that: "All 64-bit SIMD integer instructions use MMX
> registers, which share register state with the x87 floating point
> stack." (from e.g. the Intel® 64 and IA-32 Architectures Optimization
> Reference Manual)
That's for MMX, not sor SSE...
> On 6/5/2020 10:52 AM, Alf P. Steinbach wrote:
>> On 05.06.2020 10:01, Ralf Goertz wrote:
><...>
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>
>> Thanks, I didn't know, but should have guessed.
>>
>> There is also an x86 instruction, fsincos, but an SO answer claims that
>>
>> "As CPUs become faster, old hardware trig instructions like fsin have
>> not kept pace. With some CPUs, a software function, using a polynomial
>> approximation for sine or another trig function, is now faster than a
>> hardware instruction. In short, fsincos is too slow."
>> <url: https://stackoverflow.com/a/24470801>
>>
>
> There's more to that: "All 64-bit SIMD integer instructions use MMX
> registers, which share register state with the x87 floating point
> stack." (from e.g. the Intel® 64 and IA-32 Architectures Optimization
> Reference Manual)
Manfred
Jun 5, 2020, 10:53:09 AM6/5/20
to
On 6/5/2020 3:18 PM, Melzzzzz wrote:
> On 2020-06-05, Manfred <non...@add.invalid> wrote:
>> On 6/5/2020 10:52 AM, Alf P. Steinbach wrote:
>>> On 05.06.2020 10:01, Ralf Goertz wrote:
>> <...>
>>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>>> double *cos)` available? It is much faster than computing sin and cos
>>>> separately.
>>>
>>> Thanks, I didn't know, but should have guessed.
>>>
>>> There is also an x86 instruction, fsincos, but an SO answer claims that
>>>
>>> "As CPUs become faster, old hardware trig instructions like fsin have
>>> not kept pace. With some CPUs, a software function, using a polynomial
>>> approximation for sine or another trig function, is now faster than a
>>> hardware instruction. In short, fsincos is too slow."
>>> <url: https://stackoverflow.com/a/24470801>
>>>
>>
>> There's more to that: "All 64-bit SIMD integer instructions use MMX
>> registers, which share register state with the x87 floating point
>> stack." (from e.g. the Intel® 64 and IA-32 Architectures Optimization
>> Reference Manual)
>
> That's for MMX, not sor SSE...
>
>
I believe my quote was clear enough: "All 64-bit SIMD integer
> On 2020-06-05, Manfred <non...@add.invalid> wrote:
>> On 6/5/2020 10:52 AM, Alf P. Steinbach wrote:
>>> On 05.06.2020 10:01, Ralf Goertz wrote:
>> <...>
>>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>>> double *cos)` available? It is much faster than computing sin and cos
>>>> separately.
>>>
>>> Thanks, I didn't know, but should have guessed.
>>>
>>> There is also an x86 instruction, fsincos, but an SO answer claims that
>>>
>>> "As CPUs become faster, old hardware trig instructions like fsin have
>>> not kept pace. With some CPUs, a software function, using a polynomial
>>> approximation for sine or another trig function, is now faster than a
>>> hardware instruction. In short, fsincos is too slow."
>>> <url: https://stackoverflow.com/a/24470801>
>>>
>>
>> There's more to that: "All 64-bit SIMD integer instructions use MMX
>> registers, which share register state with the x87 floating point
>> stack." (from e.g. the Intel® 64 and IA-32 Architectures Optimization
>> Reference Manual)
>
> That's for MMX, not sor SSE...
>
>
instructions...", which are part of SSE.
James Kuyper
Jun 5, 2020, 11:03:52 AM6/5/20
to
On 6/5/20 5:32 AM, David Brown wrote:
> On 05/06/2020 10:01, Ralf Goertz wrote:
...
> On 05/06/2020 10:01, Ralf Goertz wrote:
>> Did you know that with g++ you have `void sincos(double x, double *sin,
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>>
>
> No, "g++" does not have such a function - it is part of glibc, not the
> compiler, and it is a non-standard extension in <math.h>. That means
> you need to define the macro "_GNU_SOURCE" before including <math.h> to
> get access to it.
>
> And it is /not/ faster than calculating sin and cos separately, because
> gcc is smart enough to combine calls to sin and cos with the same
> operand to a call to sincos.
Actually, it's a little too smart about that - the way it's implemented
>> double *cos)` available? It is much faster than computing sin and cos
>> separately.
>>
>
> No, "g++" does not have such a function - it is part of glibc, not the
> compiler, and it is a non-standard extension in <math.h>. That means
> you need to define the macro "_GNU_SOURCE" before including <math.h> to
> get access to it.
>
> And it is /not/ faster than calculating sin and cos separately, because
> gcc is smart enough to combine calls to sin and cos with the same
> operand to a call to sincos.
can break well-formed code even when compiling with every option turned
on to make g++ conform as accurately as possible. See
<https://gcc.gnu.org/bugzilla/show_bug.cgi?id=46926>. It's still an open
bug nearly a decade later - no one has been assigned to work on it yet.
That was a bug report against gcc about it's handling of C code.
However, I just rewrote the code using C++, and it fails the same way
when compiled using g++. g++ --version says
(Ubuntu 7.5.0-3ubuntu1~18.04) 7.5.0/
Oddly enough, I can no longer reproduce this bug using gcc.
Summary:
A third-party library I was linking to defined it's own function named
sincos(). The definition of that function didn't use any of the tricks
that glibc's sincos() does to make it significantly faster than separate
calls to sin() and cos(). It simply made separate calls to sin() and
cos(). That's dumb, but the code was strictly conforming C (and my
modified version is well-formed C++). sincos is not a reserved
identifier in either language, so there's nothing wrong with naming a
function that way.
My own code had a pair of sin() and cos() function calls with the same
argument.
In both cases, gcc noticed the sin() and cos() calls, and converted them
into a single call to sincos(). This made the third-party library's
sincos() function infinitely recursive, and caused my function to call
that one - the program ended up aborting when it ran out of stack.
gcc was invoked in a fully-conforming mode which means this was not a
valid optimization. It should either have converted the sin() and cos()
calls to call to a function with a name reserved to the implementation,
such as _Sincos(), or it shouldn't have fiddled with the sin() and cos()
calls at all.
David Brown
Jun 5, 2020, 11:37:00 AM6/5/20
to
On 05/06/2020 17:03, James Kuyper wrote:
> On 6/5/20 5:32 AM, David Brown wrote:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
> ...
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not the
>> compiler, and it is a non-standard extension in <math.h>. That means
>> you need to define the macro "_GNU_SOURCE" before including <math.h> to
>> get access to it.
>>
>> And it is /not/ faster than calculating sin and cos separately, because
>> gcc is smart enough to combine calls to sin and cos with the same
>> operand to a call to sincos.
>
> Actually, it's a little too smart about that - the way it's implemented
> can break well-formed code even when compiling with every option turned
> on to make g++ conform as accurately as possible.
Yes, I saw that bug - and I agree it appears to be a bug. I've seen
> On 6/5/20 5:32 AM, David Brown wrote:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
> ...
>>> Did you know that with g++ you have `void sincos(double x, double *sin,
>>> double *cos)` available? It is much faster than computing sin and cos
>>> separately.
>>>
>>
>> No, "g++" does not have such a function - it is part of glibc, not the
>> compiler, and it is a non-standard extension in <math.h>. That means
>> you need to define the macro "_GNU_SOURCE" before including <math.h> to
>> get access to it.
>>
>> And it is /not/ faster than calculating sin and cos separately, because
>> gcc is smart enough to combine calls to sin and cos with the same
>> operand to a call to sincos.
>
> Actually, it's a little too smart about that - the way it's implemented
> can break well-formed code even when compiling with every option turned
> on to make g++ conform as accurately as possible.
similar recursive calls generated from home-made memcpy functions (the
compiler sees the contents of the function, recognizes it as a pattern
matching memcpy, so turns the loop into a call to memcpy...). But in
that case, the compiler has the excuse that memcpy is a standard
function, and not one you should be implementing.
I think the way to solve these things would be to have implementation
functions in the library use internal (reserved identifier) names, and
make the user-visible names be weak aliases. That way you could define
your own sincos() function, and if the compiler figured out the best
implementation was to call __library_internal_sincos(), there would be
no problem.
Keith Thompson
Jun 5, 2020, 2:59:06 PM6/5/20
to
Ralf Goertz <m...@myprovider.invalid> writes:
> Am Fri, 5 Jun 2020 11:32:28 +0200
> schrieb David Brown <david...@hesbynett.no>:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>> > Did you know that with g++ you have `void sincos(double x, double
>> > *sin, double *cos)` available? It is much faster than computing sin
>> > and cos separately.
>>
>> No, "g++" does not have such a function - it is part of glibc, not
>> the compiler, and it is a non-standard extension in <math.h>. That
>> means you need to define the macro "_GNU_SOURCE" before including
>> <math.h> to get access to it.
>
> Well, that's why I wrote that the function was available with g++. As
> Alf mentioned MinGW he'd be using glibc in that case.
I thought MinGW used Microsoft's C runtime, not glibc.
> Am Fri, 5 Jun 2020 11:32:28 +0200
> schrieb David Brown <david...@hesbynett.no>:
>> On 05/06/2020 10:01, Ralf Goertz wrote:
>> > Did you know that with g++ you have `void sincos(double x, double
>> > *sin, double *cos)` available? It is much faster than computing sin
>> > and cos separately.
>>
>> No, "g++" does not have such a function - it is part of glibc, not
>> the compiler, and it is a non-standard extension in <math.h>. That
>> means you need to define the macro "_GNU_SOURCE" before including
>> <math.h> to get access to it.
>
> Well, that's why I wrote that the function was available with g++. As
> Alf mentioned MinGW he'd be using glibc in that case.
In any case, the sincos() function is provided by the runtime
library, not by the compiler. It's available if you use clang++
with glibc, and it's provided by other library implementations as
well (musl and newlib both provide it).
--
Keith Thompson (The_Other_Keith) Keith.S.T...@gmail.com
Working, but not speaking, for Philips Healthcare
void Void(void) { Void(); } /* The recursive call of the void */
Chris M. Thomasson
Jun 5, 2020, 4:42:18 PM6/5/20
to
On 6/4/2020 7:06 PM, Chris M. Thomasson wrote:
> On 6/4/2020 1:16 PM, Alf P. Steinbach wrote:
>> On 03.06.2020 21:12, Alf P. Steinbach wrote:
>>> On 03.06.2020 13:37, Alf P. Steinbach wrote:
>>>> On 03.06.2020 08:43, Chris M. Thomasson wrote:
>>>>> On 6/2/2020 11:39 PM, Chris M. Thomasson wrote:
>>>>>> On 6/2/2020 9:13 PM, Alf P. Steinbach wrote:
>>>>>>> On 02.06.2020 08:30, Chris M. Thomasson wrote:
> [...]
>
> Wrt fast, check this out wrt obtaining the hypotenuse:
>
> https://forums.parallax.com/discussion/147522/dog-leg-hypotenuse-approximation
Humm... Fwiw, imvho, it would be fun to use this interesting
> On 6/4/2020 1:16 PM, Alf P. Steinbach wrote:
>> On 03.06.2020 21:12, Alf P. Steinbach wrote:
>>> On 03.06.2020 13:37, Alf P. Steinbach wrote:
>>>> On 03.06.2020 08:43, Chris M. Thomasson wrote:
>>>>> On 6/2/2020 11:39 PM, Chris M. Thomasson wrote:
>>>>>> On 6/2/2020 9:13 PM, Alf P. Steinbach wrote:
>>>>>>> On 02.06.2020 08:30, Chris M. Thomasson wrote:
> [...]
>
> Wrt fast, check this out wrt obtaining the hypotenuse:
>
> https://forums.parallax.com/discussion/147522/dog-leg-hypotenuse-approximation
approximation for normalization of vectors in my experimental field.
sqrt is used to gain the length of a vector wrt normalization, well,
that makes a triangle, that can be used with the approx...
For instance the following animation uses the test vector field with
sqrt. It uses cos sin as well to gain the spirals wrt interpolating from
0...pi/2 where 0 is neutral and pi/2 is the equipotential in a 2d field.
However, I want to focus on the sqrt approx for now. Actually, I do not
need cos and sin to gain the perpendicular field in 2d. I use them to be
able to smoothly animate an interpolation of angles from 0...pi/2.
https://youtu.be/JIM-QioOhdY
Wonder what the approximation would look like visually when compared to
the original... Humm...
Fwiw, this can go 3d, humm... Need to add a z parameter the the dog leg
approx... In 2d they are alwyas zero. Here is an example in 3d using sqrt:
https://skfb.ly/6QU86
Here is a field where the lines are drawn from a 3d grid:
https://skfb.ly/6St9p
Chris M. Thomasson
Jun 5, 2020, 4:45:47 PM6/5/20
to
Ike Naar
Jun 6, 2020, 1:45:17 AM6/6/20
to
On 2020-06-04, Alf P. Steinbach <alf.p.stein...@gmail.com> wrote:
> For the express-sin-vis-cos possibility I measured these functions:
>
> struct Stdlib
> {
> static auto func( const double v )
> -> double
> { return sin( v )*cos( v ); }
> };
>
> struct Pythagoras
> {
> static auto func( const double v )
> -> double
> {
> const double cos_value = cos( v );
> const double abs_sin_value = sqrt( 1 - cos_value*cos_value );
> return (sgn( v )*abs_sin_value)*cos_value;
> }
> };
Have you also the identity
> For the express-sin-vis-cos possibility I measured these functions:
>
> struct Stdlib
> {
> static auto func( const double v )
> -> double
> { return sin( v )*cos( v ); }
> };
>
> struct Pythagoras
> {
> static auto func( const double v )
> -> double
> {
> const double cos_value = cos( v );
> const double abs_sin_value = sqrt( 1 - cos_value*cos_value );
> return (sgn( v )*abs_sin_value)*cos_value;
> }
> };
sin(v) * cos(v) == 0.5 * sin(2.0 * v)
? I haven't measured, but I would not be surprised if computing
0.5 * sin(2.0 * v)
would run faster than either Stdlib or Pythagoras.
Alf P. Steinbach
Jun 6, 2020, 5:29:07 AM6/6/20
to
sin*cos, but to compute both sin and cos. The sin*cos function result is
just a dependency that, when e.g. stored in a `volatile` variable,
prevents the compiler from optimizing away the calls.
The sin and cos values are used by the OP's code to find a point on a
circle in the complex number plane, specifically one of the n n'th roots
of a complex number.
- Alf
Chris M. Thomasson
Jun 6, 2020, 2:49:57 PM6/6/20
to
On 6/6/2020 2:28 AM, Alf P. Steinbach wrote:
> On 06.06.2020 07:45, Ike Naar wrote:
>> On 2020-06-04, Alf P. Steinbach <alf.p.stein...@gmail.com> wrote:
[...]
> On 06.06.2020 07:45, Ike Naar wrote:
>> On 2020-06-04, Alf P. Steinbach <alf.p.stein...@gmail.com> wrote:
> The sin and cos values are used by the OP's code to find a point on a
> circle in the complex number plane, specifically one of the n n'th roots
> of a complex number.
Correct. Now, it would be nice to keep everything in polar form.
> circle in the complex number plane, specifically one of the n n'th roots
> of a complex number.
However, I need to do two things with it. Wrt my highly experimental
fractal root storage thing, I need to add two complex numbers together.
This is not ideal in polar form. So, I convert it into rectangular form
via cos and sin. Also, the need to actually plot these damn points is in
order. Therefore rectangular form is "needed" to get at the x and y
components anyway.
Fwiw, I am storing information in the roots of these numbers in the
context of a fractal, basically a Julia set. Here is some context on
visually plotting these roots:
https://en.wikipedia.org/wiki/Julia_set#Using_backwards_(inverse)_iteration_(IIM)
Take a fractal, say:
z = z^7+c
Well, this actually creates a seven ary tree to store data. I am doing
this for fun, and it works. An interesting part to me is the the fact
that any message can be encoded in a single complex number using this
method. The problem is precision... An arbitrary precision package can
be used. Its not all that practical to store a large message in a single
complex number for the damn real and imaginary parts can contain
hundreds of digits each...
The main idea is instead of randomly choosing a root, let a stream of
real data determine what root to follow during iteration.
Btw, is there a technique to add two complex numbers together in polar
form that is faster than rectangular form? So far, I am having some
trouble finding any.
Alf P. Steinbach
Jun 6, 2020, 3:07:13 PM6/6/20
to
On 06.06.2020 20:49, Chris M. Thomasson wrote:
>
> [snip]
PForm( v1, m1 ) + PForm( v2, m2 ) = PForm( seriously inefficient calc )
The one exception is when the angles v1 and v2 are identical or
opposite, in which case the magnitudes m1 and m2 can just be added or
subtracted, respectively.
On the other hand /multiplication/ is more efficient in polar form than
in rectangular form, namely 1 add + 1 mul versus more for rectangular.
Disclaimer: my math knowledge is from college, early 1980's.
- Alf
>
> [snip]
> Btw, is there a technique to add two complex numbers together in polar
> form that is faster than rectangular form? So far, I am having some
> trouble finding any.
Rform( re1, im1 ) + Rform( re2, im2 ) = Rform( re1 + re2, im1 + im2 )
> form that is faster than rectangular form? So far, I am having some
> trouble finding any.
PForm( v1, m1 ) + PForm( v2, m2 ) = PForm( seriously inefficient calc )
The one exception is when the angles v1 and v2 are identical or
opposite, in which case the magnitudes m1 and m2 can just be added or
subtracted, respectively.
On the other hand /multiplication/ is more efficient in polar form than
in rectangular form, namely 1 add + 1 mul versus more for rectangular.
Disclaimer: my math knowledge is from college, early 1980's.
- Alf
0 new messages