iT邦幫忙

1

Euler 圓周率函數問題求解

數學大師歐勒(Euler) 找到一個計算圓周率的無窮乘積 :
pi/2 = 3/2 * 5/6 * 7/6 * 11/10 * 13/14 * 17/18 * 19/18 * 23/22 * …

有趣的是, 這公式裡, 所有分子都是大於 2的質數, 分母則是不能被 4整除, 且最接近分子的偶數.

試撰寫一函數 double Euler(int n), 用來估算圓周率的值到第 n項, 並計算 Euler(10), Euler(100), Euler(1000) 與 Euler(10000) 的結果.

我稍微修改了一下題目 想讓程式程式從Euler(1)跑到Euler(n)
不過程式碼有問題 求解惑 謝謝

#include <stdio.h>
#include <stdlib.h>
double Euler(int);
int is_prime(int);
int find_k(int);

int main(void)
{
    int i,n=10;
    for(i=1;i<=n;i++)
    {
    	printf("Euler(%d)=%lf\n",i,Euler(i));
	}
    system("pause");
	return 0;
}

double Euler(int n)
{
	int i;
	double sum=1.0;
	for(i=1;i<=n;i++)
	{
		sum*=is_prime(i)/find_k(i);
	}
	return sum*2;
}

int is_prime(int n)
{
	int k=3,cnt=1;
	while(n>=cnt)
	{
		int i;
		for(i=2;i<k;i++)
		{
			if(k%i==0)
			{
				k++;
				continue;
			}
		}
		cnt++;
	}
	return k;
}

int find_k(int n)
{
	int quo;
    int x=is_prime(n);
	quo=(x-2)/4;
	if(x-quo*4+2 > (quo+1)*4+2-n)	
	{
		return (quo+1)*4+2;
	}
	else
	{
		return quo*4+2;
	}
	
}

2 個回答

1
小魚
iT邦高手 1 級 ‧ 2018-06-15 20:54:29

你的邏輯不大對,
剛好有時間稍微玩了一下,

以下是執行的結果
至少到小數點第2位是對的...
https://ithelp.ithome.com.tw/upload/images/20180615/201056945mGEkDzpxv.jpg

以下是程式碼

#include <iostream>
#include <vector>
using namespace std;

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

vector<int> primeList;
double Euler(int);
void find_prime(int);
bool is_prime(int);
int find_k(int);

int main()
{
    int i,n=10000;

    //先找出n個質數
    find_prime(n);

    for(i=1;i<=n;i++)
    {
        printf("Euler(%d)=%lf\n",i,Euler(i));
    }

    system("pause");
    return 0;
}

double Euler(int n)
{
    double sum = 1;
    for(int i=0;i<n;i++)
    {
        int prime = primeList[i];
        sum*= (double)prime /find_k(prime);
    }
    return sum*2;
}

//找出n個質數
void find_prime(int n)
{
    int count = 0;
    int number = 0;
    primeList.push_back(3);
    primeList.push_back(5);
    primeList.push_back(7);
    primeList.push_back(11);
    count = 4;
    number = 12;

    while(count < n)
    {
        //如果是質數
        if(is_prime(number))
        {
            primeList.push_back(number);
            count++;
        }
        number++;
    }
}

//判斷是否是質數
bool is_prime(int n)
{
    for(int i=2;i<n/2;i++)
    {
        if(n%i==0)
            return false;
    }
    return true;
}

//取得比n小的2的倍數但非4的倍數
int find_k(int n)
{
    int k = n;
	if(k % 4 != 0 && k % 2 == 0)
		return k;
	else if(k % 4 == 1)
		return k+1;
	else if(k % 4 == 3)
		return k-1;
	else
		return 1;
}

看更多先前的回應...收起先前的回應...

謝謝 然後我下面修改後的
程式碼經過測試應該是沒有問題

想問大大 請問 原程式碼 是哪個區塊的邏輯出現問題 謝謝/images/emoticon/emoticon41.gif

謝謝大大
剛剛我稍微做了修正

#include <stdio.h>
#include <stdlib.h>
double Euler(int);
int is_prime(int);
int find_k(int);
int main(void)
{
    int i,n=10000;
    for(i=1;i<=n;i++)
    {
    	printf("Euler(%d)=%lf\n",i,Euler(i));
	}
    system("pause");
	return 0;
}

double Euler(int n)
{
	int k=3,cnt=0;
	double sum=1.0;
	while(cnt<n)
	{
		if(is_prime(k))
		{
			sum*=(double)k/find_k(k);
			cnt++;
		}
		k++;
	}
	return sum*2;
}

int is_prime(int n)
{
	int i;
	for(i=2;i<n;i++)
	{
		if(n%i==0)
			return 0;
	}
	return 1;
}

int find_k(int n)
{
	int quo;
	quo=(n-2)/4;
	if(n-(quo*4+2) > 4*(quo+1)+2-n)	
	{
		return (4*(quo+1)+2);
	}
	else
	{
		return quo*4+2;
	}
	
}

我想我的問題大概是讓is_find 功能複雜化
最初的寫法是想讓is_find 找出3之後的prime 後再傳給find_k找最接近的偶數 還請大大指教 是哪裡的邏輯出了問題

後來把is_prime 想做的事稍微改了一下讓Euler函數去解決
至於find_k跟大大寫法不同 感謝

小魚 iT邦高手 1 級 ‧ 2018-06-15 22:50:55 檢舉

原本的is_prime應該改這樣,
就可以正常執行了,

int is_prime(int n)
{
	int k=3,cnt=1;
	while(n>cnt)
	{
		int i;
		for(i=2;i<k;i++)
		{
			if(k%i==0)
			{
				cnt++;
				break;
			}
		}
		k++;
	}
	return k;
}

另外,質數算到 根號n就可以了,再算下去是多餘的...
不過你每次都要重頭算,
你把n調到10000試看看要跑多久...
不過我的程式跑100000也要跑很久,
我還沒讓它跑完過...

然後is_prime已經算了一次,
find_k裡面就不要再算一次了...
寫程式是要講究效率的...

好的 謝謝您的指教/images/emoticon/emoticon41.gif

小魚 iT邦高手 1 級 ‧ 2018-06-15 23:28:07 檢舉

也對齁,
我只有看前面兩個數而已...
後面沒注意看...

不過最好的方法還是先算好n個質數,
放進列表中,
才不會每次都要算...

0
fysh711426
iT邦研究生 4 級 ‧ 2018-06-17 12:57:29

看到熟悉的題目,前陣子朋友也問過這題,
我的做法和大大類似,不過因為是直接拿算單次的程式來改,
所以就像小魚大大說的 find 函數每次重算,有效能上的問題。
/images/emoticon/emoticon16.gif

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

//判斷是否為質數
bool isPrime(int num)
{
	for (int i = 2; i <= num / 2; i++)
		if (num % i == 0)
			return false;
	return true;
}

//不能被 4 整除,且最接近分子的偶數
int num(int num)
{
	if ((num - 1) % 4 != 0)
	{
		return num - 1;
	}
	return num + 1;
}

//找到第 n 項
double find(int n)
{
	double sum = 1.0;
	int prime = 3;
	int count = 0;

	while (count < n)
	{
		if (isPrime(prime))
		{
			sum *= (double)prime / num(prime);
			count++;
		}
		prime++;
	}
	return sum * 2;
}

int main(void)
{
	int n = 10000;

	for (int i = 1; i <= n; i++)
	{
		printf("Euler(%d)=%lf\n", i, find(i));
	}

	system("pause");
	return 0;
}

我要發表回答

立即登入回答