立方根の計算 in C DouKaku?

xは0以上1000未満の実数です。 y * y * y = xになるような実数y(立方根)を小数点以下12桁以上の正確さで求める関数cube_rootを作って下さい。

ただし、このお題の趣旨は実数区間での探索なので、立方根関数があっても使ってはいけません。指数関数と対数関数も禁止します。

Pythonで表現した入出力の例：

>>> cube_root(10.0)
2.1544346900318834
>>> _ ** 3
9.9999999999999947
>>> cube_root(100.0)
4.6415888336127793
>>> _ ** 3
100.00000000000003

% ./cube_root 10
2.1544346900318895

% ./cube_root 100
4.6415888336127781

#include <stdio.h>
#include <stdlib.h>

long double cube_root(long double x, long double min, long double max){
    long double m = (min + max) / 2;
    long double d = m * m * m - x;
    if(d > -1E-13 && d < 1E-13) return m;
    if(d > 0)
        return cube_root(x, min, m);
    else
        return cube_root(x, m, max);
}

int main(int argc, char *argv[])
{
    long double x, y;
    if(argc < 2) return EXIT_FAILURE;
    sscanf(argv[1], "%llf", &x);
    y = cube_root(x, 0, 1000);
    printf("%.16llf\n", y);
    return EXIT_SUCCESS;
}

相対誤差も絶対誤差も 1e-12 以下になるようにしてみました。浮動小数点演算の誤差は相対誤差を基本に考えた方が幸せになれるような気がしたので。

#include<stdio.h>
#include<stdlib.h>
#include<time.h>
#include<math.h>
#include<assert.h>

#define EPS 1e-12

double cubicrt(double x)
{
    double y,oldy,z=1.0;
    if(x==0.0){
        return 0.0;
    }
    while(fabs(x)>8){
        x/=8;
        z*=2;
    }
    while(fabs(x)<0.125){
        x*=8;
        z/=2;
    }
    y=x;
    do{
        oldy=y;
        y=2*oldy/3+x/(3*oldy*oldy);
    }while(fabs(y-oldy)>fabs(EPS*oldy) || fabs(y-oldy)>EPS);

    return y*z;
}

void test(double x,double y,double truey,double eps)
{
    assert(fabs(y-truey)<=fabs(eps*y) && fabs(y-truey)<=eps);
    return;
}

int main()
{
    int i=0;
    double x,y;
    srand(time(NULL));

    /* テスト */
    for(i=0;i<1000;i++){
        x=rand()*2000.0/RAND_MAX-500;
        y=cubicrt(x);
        test(x,y*y*y,x,EPS);
    }
    for(i=0;i<1000;i++){
        x=rand()*2e-310/RAND_MAX-1e-310;
        y=cubicrt(x);
        test(x,y*y*y,x,EPS);
    }
    return 0;
}

DouKaku?

立方根の計算

Posted feedbacks - C

Index

Feed

Other

Link

Pathtraq