122 lines
3.0 KiB
C++
122 lines
3.0 KiB
C++
//This file is part of Bertini 2.
|
|
//
|
|
//bertini2/double_extensions.hpp is free software: you can redistribute it and/or modify
|
|
//it under the terms of the GNU General Public License as published by
|
|
//the Free Software Foundation, either version 3 of the License, or
|
|
//(at your option) any later version.
|
|
//
|
|
//bertini2/double_extensions.hpp is distributed in the hope that it will be useful,
|
|
//but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
//GNU General Public License for more details.
|
|
//
|
|
//You should have received a copy of the GNU General Public License
|
|
//along with bertini2/double_extensions.hpp. If not, see <http://www.gnu.org/licenses/>.
|
|
//
|
|
// Copyright(C) 2015 - 2017 by Bertini2 Development Team
|
|
//
|
|
// See <http://www.gnu.org/licenses/> for a copy of the license,
|
|
// as well as COPYING. Bertini2 is provided with permitted
|
|
// additional terms in the b2/licenses/ directory.
|
|
|
|
// individual authors of this file include:
|
|
// silviana amethyst, University of Wisconsin Eau Claire
|
|
|
|
/**
|
|
\file bertini2/double_extensions.hpp
|
|
|
|
\brief Provides bertini extensions to the double and complex<double> types.
|
|
*/
|
|
|
|
#ifndef BERTINI_DOUBLE_EXTENSIONS_HPP
|
|
#define BERTINI_DOUBLE_EXTENSIONS_HPP
|
|
|
|
#pragma once
|
|
|
|
#include <random>
|
|
#include <complex>
|
|
|
|
|
|
namespace bertini{
|
|
|
|
using dbl = std::complex<double>;
|
|
using dbl_complex = std::complex<double>;
|
|
|
|
/**
|
|
\brief Overload * for unsigned * complex<double>
|
|
*/
|
|
inline
|
|
std::complex<double> operator*(unsigned i, std::complex<double> z)
|
|
{
|
|
z*=i;
|
|
return z;
|
|
}
|
|
|
|
/**
|
|
an overload of isnan, for std::complex<double>
|
|
*/
|
|
inline
|
|
bool isnan(std::complex<double> const& z)
|
|
{
|
|
using std::isnan;
|
|
return isnan(z.real()) || isnan(z.imag());
|
|
}
|
|
|
|
|
|
/**
|
|
\brief Gets you a random real number between -1 and 1, fwiw
|
|
*/
|
|
inline
|
|
double RandReal()
|
|
{
|
|
static std::default_random_engine generator;
|
|
static std::uniform_real_distribution<double> distribution(-1.0,1.0);
|
|
return distribution(generator);
|
|
}
|
|
|
|
namespace{
|
|
using dbl = std::complex<double>;
|
|
}
|
|
|
|
/**
|
|
Compute +,- integral powers of a std::complex<double> number.
|
|
|
|
This function recursively calls itself if the power is negative, by computing the power on the inverse.
|
|
|
|
\note This overload was removed from C++ in C++11, for some insane reason. Here it is, back in black.
|
|
*/
|
|
inline dbl pow(const dbl & z, int power)
|
|
{
|
|
if (power < 0) {
|
|
return pow(1./z, -power);
|
|
}
|
|
else if (power==0)
|
|
return dbl(1,0);
|
|
else if(power==1)
|
|
return z;
|
|
else if(power==2)
|
|
return z*z;
|
|
else if(power==3)
|
|
return z*z*z;
|
|
else
|
|
{
|
|
unsigned int p(power);
|
|
dbl result(1,0), z_to_the_current_power_of_two = z;
|
|
// have copy of p in memory, can freely modify it.
|
|
do {
|
|
if ( (p & 1) == 1 ) { // get the lowest bit of the number
|
|
result *= z_to_the_current_power_of_two;
|
|
}
|
|
z_to_the_current_power_of_two *= z_to_the_current_power_of_two; // square z_to_the_current_power_of_two
|
|
} while (p >>= 1);
|
|
|
|
return result;
|
|
}
|
|
}
|
|
|
|
} // namespace bertini
|
|
|
|
|
|
#endif // include guard
|
|
|