#動機
数値計算でxの二乗などをするとき、よく#define SQUARE(x) ((x)*(x))
などとし、std::powは(遅いので)使うなと言われる[要出典]
それもあるのだが、単に整数の整数乗(3^5など)を計算したいときも、いちいちdoubleにキャストしてpowして戻り値をintにキャストするのが馬鹿らしい。(c++に累乗演算子などなかった)
#コード(没)
pow(x,n)
で、x,nがコンパイル時に分かっている場合にはconstexpr演算するのはもちろん、nが整数の場合には再帰的に掛け算してxの型のまま計算するコードを書いた。(なお、c++初心者でSFINAEとかc++11 constexprとか扱えないのでif constexprを使ったため、c++17以降でないと使えません)
hnn::pow(x,n)
は、nの型が整数かによって自動的にふさわしいpowを呼んでくれます。
#include <iostream>
#include <type_traits>
#include <cmath>
#if __has_include("sprout/math.hpp")
#include "sprout/math.hpp"
#define HAS_SPROUT_POW
#endif
namespace hnn{
template<class X,class N>
class pow{
X val = 0;
public:
constexpr pow(X x,N n){
//Nが整数なら倍々イディオムでO(logN)で計算
if constexpr (std::is_integral<N>::value){
if (0 > n){
val = 1 / hnn::pow(x, -n);
return;
}
val = 1;
while (0 < n) {
if ((n % 2) == 0) {
x *= x;
n /= 2;
}
else {
val *= x;
--n;
}
}
}
else{
#ifdef HAS_SPROUT_POW
val = sprout::pow(x,n);
#else
val = std::pow(x,n);
#endif
}
}
constexpr operator X() const{
return val;
}
};
}//namespace hnn{
int main(int argc, char* argv[]){
constexpr auto a=hnn::pow(1.0000001,10000000); //自然対数の近似値
auto b=hnn::pow(2.0,3.1);
constexpr auto c=hnn::pow(1.9,60);
constexpr auto d=hnn::pow(1.9,60.0);
using namespace std;
cout<<a<<","<<b<<","<<c<<","<<d<<endl;
if(argc>1){
const int value = std::atoi(argv[1]);
//int型のまま計算
const int valueCubed=hnn::pow(value,3);
cout<<"input cubed="<<valueCubed<<std::endl;
}
return 0;
}
#何が問題か
一見constexpr pow(X x,N n)
のおかげでnがコンパイル時定数ならx*x*x
のように展開してくれると思いましたが、そううまくは行きません。constexpr
指定はx,n両方がコンパイル時定数の時にしか働かないので、xがコンパイル時に分からない場合はnが定数であってもアセンブラでのcall命令が残り、マクロの場合#define SQUARE(x) ((x)*(x))
にパフォーマンスで負けてしまいます。nを再帰テンプレートで扱うときちんと展開されますが、今度はnがコンパイル時に分からない場合にコンパイルエラーになってしまいます・・・ぐぬぬ
#そのpowやばくないですか
xのn乗を計算するhnn::pow(x,n)
で、
- nがコンパイル定数かつ整数なら
x*x*x*x...
に展開 - x,nが共にコンパイル定数ならsprout::powに投げる
- そうでないならstd::powに投げる
このpow、機能多くないですか?そもそも、hnn::pow(x,n)の行を見た段階で、x,nの型を見に行かないと計算量やふるまいが分からないのは可読上良くありませんね。まるで邪悪なc形式キャストみたいだ
#顧客が本当に欲していたpow
元々の問題は#define SQUARE(x) ((x)*(x))
みたいなマクロをインライン関数などに直したい!ってことでした。その際、
・xの型のままで計算する(std::powのように浮動小数へのキャストをしたくない)
・インライン関数でsquare(double x),cubic(double x)...
などと並べるのはいかにも格好悪い(任意の整数乗に対応したい)
・パフォーマンスはマクロの場合を下回ってはならない
ということは、pow<N>(x)
のNに整数をテンプレートで渡して、xはクラステンプレートと引数の型推論を使って任意の型に対応させればよいですね。
#コード(真)
#include <iostream>
namespace hnn{
template <class T,int N>
struct Impl{
static constexpr T calc(T x){
if constexpr(N < 0){
return 1/x*Impl<T,-N>::calc(x);
}
else if constexpr(N&1){
return x*Impl<T,N-1>::calc(x);
}
else if constexpr(N==0){
(void)x;
return 1.0;
}
else{
return Impl<T,N/2>::calc(x*x);
}
}
};
template <int N,class T>
constexpr T pow(T x){
return Impl<T,N>::calc(x);
}
}//namespace hnn{
int main(int argc, char* argv[]){
//コンパイル時自然対数計算
constexpr auto napier=hnn::pow<1234567>(1.0+1.0/1234567);
//コンパイル時に分からない数の3乗
const auto input = argc>1? std::atoi(argv[1]) :0;
const auto cubic=hnn::pow<3>(input);
//確認
std::cout<<napier<<","<<cubic<<std::endl;
return 0;
}
#結論
これで少なくてもgcc -O2ならちゃんとインライン展開してくれました。あとは既存コードでSQUARE(x)
をhnn::pow<2>(x)
に変えるだけですね。c++に累乗演算子と小数剰余演算子はあっていいと思う