關(guān)于嚴(yán)謹(jǐn)性的聲明:
在用C語(yǔ)言進(jìn)行定積分的計(jì)算之前,我需要聲明以下幾點(diǎn):
一、我們所進(jìn)行定積分計(jì)算的函數(shù)都是應(yīng)當(dāng)是黎曼可積的,這保證了我們即使均勻地分割區(qū)間也保證了積分的收斂性。
二、我們同時(shí)還應(yīng)該認(rèn)識(shí)到,鑒于某些函數(shù)不一定是在每一點(diǎn)上都是可導(dǎo)的,因此我們將不使用所謂的利用到導(dǎo)數(shù)來(lái)提高計(jì)算精度的這種算法。
三、考慮到某些函數(shù)并不是連續(xù)的,如果涉及到這些計(jì)算,那么我們計(jì)算到的結(jié)果將會(huì)是一個(gè)瑕積分,如果這個(gè)瑕積分是收斂的,那么我們得到的就是收斂值。如果不收斂,這個(gè)瑕積分的值并不是積分主值,因?yàn)槲覀儧](méi)有辦法保證瑕點(diǎn)附近的兩個(gè)點(diǎn)的趨近于瑕點(diǎn)的趨勢(shì)是一樣的。
四、如果你計(jì)算的是一個(gè)無(wú)窮積分,如果這個(gè)無(wú)窮積分是收斂的,那么我們得到的就是收斂值。如果不收斂,那么我們得到的值將會(huì)是它的積分主值。
另外,考慮到C語(yǔ)言精度的限制,我們無(wú)需寫一個(gè)專門的程序用來(lái)計(jì)算瑕積分,因?yàn)楫?dāng)我們用C語(yǔ)言計(jì)算這個(gè)瑕積分的時(shí)候,就相當(dāng)于我們采用了兩邊趨近瑕點(diǎn)取極限的辦法來(lái)計(jì)算瑕積分,當(dāng)然了,由于我們采用了均勻分割的方法,所以當(dāng)我們兩邊趨近的時(shí)候并不能使得其以任意方式趨近瑕點(diǎn)。除此之外,你應(yīng)當(dāng)保證瑕點(diǎn)不在積分區(qū)間的邊界上。對(duì)于無(wú)窮積分,不論我們?cè)趺辞懈顓^(qū)間都不能使得劃分的區(qū)間足夠小,除非我們劃分區(qū)間的時(shí)候?qū)⑺麆潪榱藷o(wú)窮份,當(dāng)然了這個(gè)劃分的份數(shù)應(yīng)當(dāng)是無(wú)窮積分區(qū)間大小的高階無(wú)窮大。事實(shí)上,這是做不到的,因此我們將采用函數(shù)變換的方法將無(wú)窮積分轉(zhuǎn)化為有限積分。
最后,鑒于C語(yǔ)言的精度限制,如果我們只是計(jì)算性質(zhì)較為良好的函數(shù),那么我們將會(huì)獲得非常好的精度,甚至是僅僅只是可測(cè)的函數(shù),也是可以計(jì)算的。但如果我們計(jì)算的是有著許多無(wú)窮間斷點(diǎn),或有著許多瑕點(diǎn)的函數(shù),或者是震蕩很厲害的函數(shù),那么它的積分值誤差將會(huì)較大。最后我還想聲明一點(diǎn),那就是我們不應(yīng)該用該程序去計(jì)算本身就是發(fā)散的積分因?yàn)槟闵踔量赡艿玫降膬H僅只是個(gè)有限值。
1. 有限積分
根據(jù)我們?cè)谥八f(shuō)的,我們可以輕易得到:
?其中:
?因此,我們不難寫出以下代碼:
//定義一個(gè)C語(yǔ)言函數(shù),用于進(jìn)行微積分的數(shù)值運(yùn)算
long double integrate(long double a, long double b, long int n)
{
long double s = 0;
long double dx = (b - a) / n;
for (long int i = 1; i < n; i++)
{
s += f(a + i * dx) * dx;
}
return s;
}
其中函數(shù) long double integrate(long double a, long double b, long int n) 就是用來(lái)計(jì)算定積分的函數(shù)。值得注意的是,在調(diào)用該函數(shù)之前,我們應(yīng)當(dāng)提前定義被積函數(shù)f(x)。其中a,b,n分別為積分下限,積分上限和區(qū)間分割數(shù)。
以上是矩形法,數(shù)值計(jì)算的精度并不是非常高,我們改進(jìn)一下,使用梯形法,得到:
??其中:
??因此,我們不難寫出以下代碼:
long double integral(long double x0, long double xt, long long n)
{
long double sum = 0;
long double dx = (xt - x0) / n;
for (long long k = 0; k < n; k++)
{
sum += ((f(x0 + k * dx) + f(x0 + (k + 1) * dx)) / 2) * dx;
}
return sum;
}
相較于矩形法,梯形法不論是精度還是速度都有著不小的提升。
2. 無(wú)窮積分
我們已經(jīng)討論過(guò)無(wú)窮積分采用直接計(jì)算的困難的地方,因此我們將采用函數(shù)變換的方法。
我們很容易注意到:
?文章來(lái)源地址http://www.zghlxwxcb.cn/news/detail-636919.html
這兩個(gè)式子是等價(jià)的,因此我們可以將所有的無(wú)窮積分轉(zhuǎn)化為無(wú)限積分。只需要利用積分區(qū)間的可疊加性即可。因此不多敘述。
?文章來(lái)源:http://www.zghlxwxcb.cn/news/detail-636919.html
有了相關(guān)的理論鋪墊后,我們添加一點(diǎn)點(diǎn)細(xì)節(jié)。以下直接給出關(guān)于定積分的計(jì)算代碼:(詳細(xì)信息在注釋里)
#include <stdio.h>
#include <math.h>
#define inf INFINITY
// 定義用于被積分的函數(shù)
long double f(long double x)
{
return exp(x);
}
// 定義被積函數(shù)的特殊形式用于計(jì)算無(wú)窮積分
long double g(long double x)
{
long double y = 1.0 / x;
long double result = f(y) * (1.0 / (x * x));
return result;
}
// 一個(gè)用于計(jì)算定積分的C語(yǔ)言子程序
//值得注意的是我們還將數(shù)學(xué)函數(shù)作為參數(shù)傳給了這個(gè)函數(shù),這可以讓我們有多個(gè)數(shù)學(xué)函數(shù)時(shí),可以指定數(shù)學(xué)函數(shù)計(jì)算
long double integral(long double x0, long double xt, long long n, long double (*f)(long double))
{
//對(duì)于黎曼可積的函數(shù),這樣的積分值為零,直接返回。
if ((x0 == xt)&&((x0 != INFINITY && x0 != -INFINITY)&&(xt != INFINITY && xt != -INFINITY)))
return 0;
else
{
long double sum = 0;
long double dx = (xt - x0) / n;
for (long long k = 0; k < n; k++)
{
sum += ((f(x0 + k * dx) + f(x0 + (k + 1) * dx)) / 2) * dx;
}
return sum;
}
}
// 定義絕對(duì)值函數(shù)
static long double lfabs(long double x)
{
if (x >= 0)
return x;
else
return -x;
}
//用于計(jì)算指定精度的定積分的函數(shù)
//值得注意的是我們還將數(shù)學(xué)函數(shù)作為參數(shù)傳給了這個(gè)函數(shù),這可以讓我們有多個(gè)數(shù)學(xué)函數(shù)時(shí),可以指定數(shù)學(xué)函數(shù)計(jì)算
long double definite_integral(long double x0, long double xt, long double esp, long double (*f)(long double))
{
//如果積分限不是無(wú)窮則進(jìn)行經(jīng)典數(shù)值積分計(jì)算
if ((x0 != INFINITY && x0 != -INFINITY) && (xt != INFINITY && xt != -INFINITY))
{
long long i = 10;
long double d = 0;
long double s = 0;
do {
long double s1 = integral(x0, xt, i, f);
long double s2 = integral(x0, xt, 10 * i, f);
i *= 10;
s = s2;
d = lfabs(s1 - s2);
} while (d > esp);
return s;
}
//如果是無(wú)窮積分,則采用特殊方法進(jìn)行數(shù)值計(jì)算
else
{
//積分下限有限大小
//值得注意的是,當(dāng)我們使用g函數(shù)進(jìn)行積分的時(shí)候,有一個(gè)積分限迎丹是0,然而事實(shí)上我們以一個(gè)很小的小數(shù)來(lái)替代
if (x0 != INFINITY && x0 != -INFINITY)
{
if (xt == INFINITY)
{
long double s1 = definite_integral(x0, 1, esp, f);
long double s2 = definite_integral(esp/1000000.0, 1, esp, g);
return (s1 + s2);
}
else if (xt == -INFINITY)
{
long double s1 = definite_integral(-1, x0, esp, f);
long double s2 = definite_integral(-1, -esp/1000000.0, esp, g);
return -(s1 + s2);
}
}
//積分上限有限大小
else if (xt != INFINITY && xt != -INFINITY)
{
if (x0 == INFINITY)
{
long double s1 = definite_integral(xt, 1, esp, f);
long double s2 = definite_integral(esp/1000000.0, 1, esp, g);
return -(s1 + s2);
}
else if (x0 == -INFINITY)
{
long double s1 = definite_integral(-1, xt, esp, f);
long double s2 = definite_integral(-1, -esp/1000000.0, esp, g);
return (s1 + s2);
}
}
//積分上下限都為無(wú)窮
else if (x0 == -INFINITY && xt == INFINITY)
{
long double s1 = definite_integral(-1, 1, esp, f);
long double s2 = definite_integral(-1, -esp/1000000.0, esp, g);
long double s3 = definite_integral(esp/1000000.0, 1, esp, g);
return (s1 + s2 + s3);
}
else if (x0 == INFINITY && xt == -INFINITY)
{
long double s1 = definite_integral(-1, 1, esp, f);
long double s2 = definite_integral(-1, -esp/1000000.0, esp, g);
long double s3 = definite_integral(esp/1000000.0, 1, esp, g);
return -(s1 + s2 + s3);
}
//如果積分是從正無(wú)窮積分到正無(wú)窮或從負(fù)無(wú)窮積分到負(fù)無(wú)窮則需要另外判斷
else
{
if (x0 == INFINITY && xt == INFINITY)
{
if (f(INFINITY) == 0)
return 0;
else
return NAN;
}
else if (x0 == -INFINITY && xt == -INFINITY)
{
if (f(-INFINITY) == 0)
return 0;
else
return NAN;
}
}
return NAN;
}
}
int main()
{
long double x0, xt, esp;
printf("如果要輸入無(wú)窮大或無(wú)窮小,請(qǐng)輸入inf或-inf\n");
printf("請(qǐng)輸入積分區(qū)域(x0,xt),和最大容忍誤差esp:>>");
scanf("%lf %lf %lf", &x0, &xt, &esp);
long double result = definite_integral(x0, xt, esp, &f);
printf("%.16lf\n", result);
return 0;
}
在計(jì)算指定定積分精度的C函數(shù)中,對(duì)于無(wú)窮積分我們通過(guò)轉(zhuǎn)化為有限積分后再調(diào)用自身的方式有效地減傷了大量大代碼量。這也是為什么要把我們定義的數(shù)學(xué)函數(shù)作為參數(shù)傳給函數(shù)的原因。
如下:
?
代碼的運(yùn)行:
?我們指定數(shù)學(xué)函數(shù)為exp(x),可以看到我們很好地計(jì)算出來(lái)該無(wú)窮積分。
?
?
?
到了這里,關(guān)于用C語(yǔ)言實(shí)現(xiàn)定積分計(jì)算(包括無(wú)窮積分/可自定義精度)的文章就介紹完了。如果您還想了解更多內(nèi)容,請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章,希望大家以后多多支持TOY模板網(wǎng)!