曲线下面积的蒙特卡罗积分

Posted 2023-02-16

【中文标题】曲线下面积的蒙特卡罗积分

【英文标题】：Monte carlo integration for area under curve

【发布时间】：2021-06-18 10:34:14

【问题描述】：

我已经编写了一个使用蒙特卡罗方法来估计 pi ​​的程序。现在我试图改为估计平面曲线的面积，特别是四叶草。See image for reference 到目前为止，我一直无法做到这一点。当然，这只涉及对我当前代码的调整？任何帮助或建议将不胜感激。这是我已经拥有的：

#include <math.h>
#include <ctime>
#include <xmemory>
using namespace std;
double pi_(double accuracy)

int n = 0, d = 0;
double x, y, latest_pi = 0;
double Origin_dist = 0;
do

    x = 0;
    y = 0;
    x = rand() % 100;
    y = rand() % 100;
    Origin_dist = sqrt(x * x + y * y);
    if (Origin_dist < 100.0)
    
        d++;
        n++;
    
    else
    
        n++;
    
    latest_pi = 4.0 * (d + 1.0) / (n + 1.0);
 while ((d < 3100) || (4.0 / (n + 1.0) < accuracy));
return latest_pi;

int main()

double accuracy;
srand((int)time(0));
cout << "Enter the accuracy: \n";
cin >> accuracy;
cout << pi_(accuracy) << endl;

【问题讨论】：

你给了一个图像，但你的代码使用距离来估计一个点是否属于光盘。您是否有一个方程可以确定一个点是否在您的四叶草中？

对不起，我没有，也没有提供问题。我只有四叶草的方程。 ((x^2 + y^2)^3) = 4 * x^2 * y^2。我也设法编写了一个程序来估计数字的自然对数。不过我想不通。

这已经是东西了。由于您的表面居中，因此解决问题的一种方法可能是为每个点找到该点与原点形成的线与位于平面同一四分之一的曲线部分之间的交点。如果你发现原点到这个交点的距离优于点到原点的距离，那么这个点就属于曲面。如果不是，或者你根本找不到交点，那么它不属于表面。

【参考方案1】：

首先，我必须明确指出，蒙特卡洛并不是解决这个问题的最佳方法。该解决方案利用了这样一个事实，即对于 x 和 y 都优于 0， ((x^2+y^2)^3 x^2y^2) => (x,y)属于表面。

#include <math.h>
#include <ctime>
#include <iostream>

double pow(double x, int n)
    double r=1.0;
    for (int i=0; i<n; i++)
        r=r*x;
        
    return x;
    

bool belongTo_quadrifolium(double x, double y)
    return pow(pow(x,2) + pow(y,2), 3) - 4 * (pow(x, 2) * pow(y, 2)) < 0; 
    



double montecarlo(double accuracy)
    unsigned int n = 0, d = 0;
    double x, y;
    do
    
        x = 1.0*rand() / (RAND_MAX-1);
        y = 1.0*rand() / (RAND_MAX-1);
        if (belongTo_quadrifolium(x, y))
            d++;
            n++;
            
        else
            n++;
            
         while (d < 3100|| (1.0/n > accuracy));
    return 4.0*d/n;
    

int main()
    double accuracy;
    srand((int)time(0));
    std::cout << "Enter the accuracy: \n";
    std::cin >> accuracy;
    std::cout << montecarlo(accuracy) << std::endl;
    return 0;
    

【讨论】：

非常感谢！这非常有帮助。

我删除了答案代码中的一些行。它们没有任何用处