Matlab Romberg求積公式求積分

2021-10-12 09:20:42 字數 1996 閱讀 6586

用romberg求積法計算積分∫−1

111+

100x2d

x\int_^ \frac} d x

∫−11​1

+100

x21​

dx的近似值,要求誤差不超過0.5×1

0−70.5\times10^

0.5×10

−7。程式1:romberg通用程式

function[result] = romberg(a,b,f,epsilon)

%% 初始賦值

e = inf;

xk = [a,b];%

t = ; % 復化梯形公式

s = ; % 復化simpson公式

c = ; % 復化cotes公式

r = ; % romberg公式

k = 1; % 當前次數

%% 計算

digits(10)

h = b - a;

t1 = vpa( h / 2 * ( f(a) + f(b) ) );

t = [t,t1];

fprintf('k=%d, 2^k=%d, t=%0.8f\n',k,2^k,t1)

while(e>epsilon)

xk2= ;

for i=1:length(xk)-1

xk2 = [xk2, (xk(i) + xk(i+1)) / 2];

end% 計算tn

sum = 0;

for i=1:length(xk2)

sum = vpa(sum + h * f(xk2(i)));

endt_temp = 0.5 * ( t(end) + sum);

t = [t,t_temp];

k = k + 1;

h = h / 2;

xk = [xk, xk2];

xk = sort(xk);

fprintf('k=%d, 2^k=%d, t=%0.8f',k,2^k,t_temp)

% 計算sn

if k>=2

s_temp = 4 / 3 * t(end) - 1 / 3 * t(end-1);

s = [s, s_temp];

fprintf(', s=%0.8f',s_temp)

end% 計算cn

if k>=3

c_temp = 16 / 15 * s(end) - 1 / 15 * s(end-1);

c = [c, c_temp];

fprintf(', c=%0.8f',c_temp)

end% 計算rn

if k>=4

r_temp = 64 / 63 * c(end) - 1 / 63 * c(end-1);

r = [r, r_temp];

fprintf(', r=%0.8f',r_temp)

end% 計算誤差e

if k>=5

e = abs(r(end) - r(end-1));

fprintf(', e=%0.8f',e)

endfprintf('\n')

endfprintf('romberg求積公式的結果為:%0.8f\n',r(end))

result = r(end);

end

程式2:主函式

clc;clear all

%% 初始資料

C 求積分(面積)

下面是兩種辦法 一 函式表示式已知 二 讀入一組資料構成一條曲線,求積分 include using namespace std 方法1,曲線的函式可以用表示式表示的情況 double func1 double x double integral double f double double min...

高斯求積公式 matlab

1.分別用三點和四點gauss chebyshev公式計算積分 並與準確積分值2arctan4比較誤差。若用同樣的三點和四點gauss legendre公式計算,也給出誤差比較結果。2 atan 4 ans 2.6516 gauss chebyshev function i gausscheby f...

利用powerful number求積性函式字首和

好久沒更部落格了,先水一篇再說。其實這個做法應該算是杜教篩的乙個拓展。powerful number的定義是每個質因子次數都 geq 2 的數。首先,leq n 的powerful number個數是 o sqrt 的,這是因為所有powerful number顯然可以表示成 a 2b 3 所以個數...