math.hを使わずにn乗根(ただしnは整数)を計算するためのC言語プログラム

数学の教育用の覚書。累乗根の近似値を具体的に計算する方法を習わずにn乗根を受け入れていたんだなあ - あらきけいすけの雑記帳の補遺。累乗根の計算でmath.hのpow()関数を使う部分が論点先取っぽくて気に入らなかったので、整数べきを計算する関数iPow()に書き換えた。iPow()はaのn乗をO(log(n))回程度の計算で実行する*1。

#include <stdio.h>

double iPow (    // n乗の計算をする関数
    double base  // (base) の (power) 乗の近似値を求める
   ,int    power // べきは整数
) {
    double buff=1;
    do {
        if ( power & 1 ) buff*= base;
        power>>= 1;
        base*= base;
    } while ( power );
    return buff;
}

double calcNthRoot ( // n乗根の計算をする関数
    double base  // (base) の (1/order) 乗の近似値を求める
   ,int    order //
   ,double PREC  // 収束の精度
){
    double root=base ,rootPrev;
    do {
        rootPrev= root;
        root=((order-1)*root)/order + base/(order*iPow(root,order-1));
    } while ( (rootPrev-root)/root > PREC );
    return root;
}

int main ( int argc, char* argv[] ) {
    double base=3; // 3 の70乗根の計算
    int order=70;
    printf("%fの%d乗根は%f\n",base,order,calcNthRoot(base,order,1.e-10));
}

*1:参考文献： http://oshiete.goo.ne.jp/qa/3609014.html [2020.6.23追記]Alexander A. Stepanov, Daniel E. Rose, 株式会社クイープ (訳), 『その数式、プログラムできますか？』, (2015, 翔泳社), §2.1 エジプト乗法, p.9. によればiPow()の算法は"Egyptian multiplication", "Russian peasant algorithm"と呼ばれる。