需要帮助修复近似 pi 的算法

我正在尝试为近似的算法编写 C 代码pi。它应该得到一个立方体的体积和该立方体内部球体的体积（球体的半径是立方体边长的 1/2）。然后我应该将立方体的体积除以球体的体积并乘以 6 以获得 pi。

它正在工作，但它在应该获得卷的部分做了一些奇怪的事情。我认为这与我为近似值选择的增量有关。 如果使用边长为 4 的立方体，则体积不是 64，而是 6400。如果是球体，则体积不是 33，而是 3334。等等。

有人能弄清楚吗？这是代码（我评论了相关部分）：

#include <stdio.h>      

int in_esfera(double x, double y, double z, double r_esfera){
    double dist = (x-r_esfera)*(x-r_esfera) + (y-r_esfera)*(y-r_esfera) + (z-r_esfera)*(z-r_esfera);

    return  dist <= (r_esfera)*(r_esfera) ? 1 : 0;   
}   

double get_pi(double l_cubo){   
    double r_esfera = l_cubo/2;   
    double total = 0;
    double esfera = 0;    
//this is delta, for the precision. If I set it to 1E anything less than -1 the program continues endlessly. Is this normal?
    double delta = (1E-1);   

    for(double x = 0; x < l_cubo; x+=delta){
        printf("x => %f; delta => %.6f\n",x,delta);
        for(double y = 0; y <l_cubo; y+=delta){
            printf("y => %f; delta => %.6f\n",y,delta);
            for(double z = 0; z < l_cubo; z+=delta){
                printf("z => %f; delta => %.6f\n",z,delta);
                total+=delta;
                if(in_esfera(x,y,z,r_esfera))
                    esfera+=delta;
            }
        }
    }

    //attempt at fixing this
        //esfera/=delta;
        //total/=delta;
    //

//This printf displays the volumes. Notice how the place of the point is off. If delta isn't a power of 10 the values are completely wrong.   
    printf("v_sphere = %.8f; v_cube = %.8f\n",esfera,total);   

    return (esfera)/(total)*6;
}   

void teste_pi(){        
    double l_cubo = 4;    
    double pi = get_pi(l_cubo);

    printf("%.8f\n",pi);
}   

int main(){   
    teste_pi();
}

total+=delta;
if(in_esfera(x,y,z,r_esfera))
    esfera+=delta;

total and esfera are three-dimensional volumes whereas delta is a one-dimensional length. If you were tracking units you'd have m³ on the left and m on the right. The units are incompatible.

要修复它，立方体delta这样你在概念上就积累了微小的立方体而不是微小的线条。

total+=delta*delta*delta;
if(in_esfera(x,y,z,r_esfera))
    esfera+=delta*delta*delta;

这样做可以修复输出，并且也适用于任何值delta:

v_sphere = 33.37400000; v_cube = 64.00000000
3.12881250

请注意，该算法“适用于”任意delta值，但存在严重的准确性问题。它非常容易出现舍入问题。当delta是二的幂：1/64.0比1/100.0， 例如：

v_sphere = 33.50365448; v_cube = 64.00000000
3.14096761

另外，如果您希望程序运行得更快，请摆脱所有这些打印输出！或者至少是内循环中的那些......