从一个椭圆积分看matlab求定积分

一个积分，能求出解析解固然好，但是求不出解析解，求一个数值解基本能够满足实际的需求了。现实生活中，我们遇到的积分大部分是解不出或者很难解出解析解，这时候，就需要我们求其数值解。matlab提供了一个求积分的强大算法。下面以椭圆积分为例，我来说一说，matlab如何求定积分。

clc
clear
syms theta;
a=5;
b=3;
c=sqrt(a^2-b^2);
e=c/a;
p=abs(a^2/c-c)
result=int(e*p/(1-e*cos(theta)),0,2*pi)/2/pi
theta_div=0:0.01:2*pi;
r_div=e*p./(1-e*cos(theta_div));
plot(theta_div,r_div)
hold on;
plot(theta_div,result)
% hold off;
% ezplot('x^2/(53^2)+ y^2/(3^2) = 1',[-53 53 -3 3])

很简单啦，就是一个int命令。所求的积分来源于一道数学题，感兴趣的同学可以看看。结果表明，椭圆焦点到周线上点的积分平均值等于短轴长度，而与长轴无关，很有意思。

但是，奇怪的是我利用椭圆的对称性，手工计算出来的结果却是长轴，显然，我做错了。聪明的你，能看出我哪里错了吗？