MATLAB计算积分

函数的积分和椭圆的周长

1．正弦函数的积分

[问题]求正弦函数从0到π的积分

y = sin x

当x = 0时，积分为0，画出积分的函数曲线。

[数学模型]

定积分的结果为

ππ00

sin d cos 2S x x x ==-=⎰ 不定积分的结果为

sin d cos I x x x C ==-+⎰

其中C 是积分常量，由初始条件决定。当x = 0时，积分为I = 0，必有C = 1。结果为

I = -cos x + 1

[算法]根据积分的基本概念，将积分区域分为多份，用矩形法求曲线下的近似面积表示积分的近似值

1()n

i i S f x x ==∆∑

矩形法的函数是sum(f)。

用梯形法求曲线下的近似面积表示积分的近似值

1

101[()()]2

n i i i S f x f x x -+==+∆∑

梯形法的函数是trapz(f)。

用数值积分的函数是quad 和quadl ，常用使用格式是

S = quad(f,a,b)

其中，f 表示被积函数，a 表示积分的下限，b 表示积分的下限。

用符号的函数是int ，常用使用格式是

S = int(f,a,b)

[程序]zqy4_1.m 如下。

%正弦函数的积分

clear %清除变量

x=linspace(0,pi); %自变量向量

dx=x(2); %间隔

y=sin(x); %被积函数

s1=sum(y)*dx %矩形法积分

s2=trapz(y)*dx %梯形法积分

f=inline('sin(x)'); %被积的内线函数

s3=quad(f,0,pi) %数值定积分

s4=int('sin(x)',0,pi) %符号积分

sc1=cumsum(y)*dx; %矩形法累积积分(精度稍差)

sc2=cumtrapz(y)*dx; %梯形法累积积分

figure %创建图形窗口

plot(x,-cos(x)+1,x,sc1,'.',x,sc2,'o') %画解析式和矩阵法以及梯形法积分曲线 s=int('sin(x)') %符号积分

sc3=subs(s,'x',x); %替换数值求符号积分的值

C=-sc3(1) %求积分常数

hold on %保持图像

plot(x,sc3+C,'c*') %画符号法积分曲线

grid on %加网格

fs=16; %字体大小

xlabel('\itx','FontSize',fs) %横坐标

ylabel('\intsin\itx\rmd\itx','FontSize',fs)%纵坐标

title('正弦函数的积分','FontSize',fs) %标题

legend('解析解','矩形法','梯形法','符号法')%图例

zqy4.1图 zqy4.2图

2.三角函数和指数的积分

[问题]求如下函数的积分

y = e ax sin bx

其中a = 0.5，b = 2。积分下限为0。画出积分的函数曲线。

[数学模型]

设 11e sin d sin de {e sin e cos d }ax ax ax ax I bx x bx bx b bx x a a

==

=-⎰⎰⎰ 11{e sin cos de }{e sin [e cos e sin d ]}ax ax ax ax ax b b bx bx bx bx b bx x a a a a =-=-+⎰⎰ 因此不定积分为

221e (sin cos )ax I a bx b bx C a b

=-++ 当x = 0时，I 应该为零，所以 22b C a b

=+

因此，从0开始的积分为

221e (sin cos )ax I a bx b bx b a b

=-++ [算法]可用积分的解析式直接画图，也可利用被积函数通过梯形求和指令cumtrapz 求积分，数值积分指令quad 求积分，还可通过符号积分int 求解。

[程序]zqy4_2cumtrapz.m 如下。

%数值积分和符号积分方法

clear %清除变量

a=0.5; %指数的常数

b=2; %正弦函数的常数

dx=0.1; %间隔

xm=6; %上限

x=0:dx:xm; %自变量向量

s1=(exp(a*x).*(-b*cos(b*x)+a*sin(b*x)))/(a^2+b^2);%积分的解析解

C=-s1(1); %求积分常数

y=exp(a*x).*sin(b*x); %被积函数

s2=cumtrapz(y)*dx; %梯形法积分

%s2=cumsum(y)*dx; %矩形法求积分

figure %创建图形窗口

plot(x,s1+C,x,s2,'.') %画积分曲线

s=['exp(',num2str(a),'*x).*sin(',num2str(b),'*x)'];%被积分函数字符串

f=inline(s); %化为内线函数

s3=0; %第1个值

for i=2:length(x) %按自变量循环

s3=[s3,quad(f,0,x(i))]; %连接积分

end %结束循环

s=int('exp(a*x)*sin(b*x)') %对x 进行符号积分

ss=subs(s,{'a','b'},{a,b}); %替换常数

s4=subs(ss,'x',x); %替换向量

hold on %保持图像

plot(x,s3,'or',x,s4-s4(1),'ch','MarkerSize',16)%画数值积分和符号积分曲线 title(['\ity\rm=e^{',num2str(a),'}\it^x\rmsin',num2str(b),'\itx\rm 的积分

'],'FontSize',16)%标题

legend('公式法','梯形法','数值法','符号法',2)%加图例

3.椭圆周长的计算

[问题]椭圆的长半轴为a ，偏心率为e ，求其周长。

[数学模型]

椭圆的方程为

x = a sin φ，y = b cos φ

其微分为

d x = a cos φd φ，y = -b sin φd φ