Matlab Romberg求積公式求積分

用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

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 所以個數...