matlab求解复杂不定积分

Question

模拟小球在仅在重力作用下沿轨迹运动
在以水平速度对x轴求不定积分求水平位移关于时间t的关系的过程中，被积函数复杂无解析解，如何解决

clear
clc
syms t_0 %%定义在目标点的角度%%
g = 9.8;
x_0=10*rand(1)+5;
y_0=10*rand(1)+5;%%随机生成目标点
eqn =@(t_0)(t_0-sin(t_0))/(1-cos(t_0))-(x_0/y_0);
t_0 = fzero(eqn, pi);%摆线在目标点t的参数t_0%
r = y_0/(1-cos(t_0));%%解出摆线圆的半径%%
T = 0:0.005:t_0;
a_1 = -rand(1); %%二次函数二次项系数%%
x = r*(T-sin(T));
y_1 = r*(1-cos(T)); %%摆线%%
y_2 = a_1*(x.^2)+(y_0-a_1*x_0^2).*x/x_0;%%抛物线%%
y_3 = y_0*x/x_0;%%直线%%
y_4 = y_0*x./0.5.*(x<=0.5)+y_0*(x>0.5);%%折线%%
%%创建速度在x位置的水平分量对应的余弦值,延拓matlab求导缺少的元素%%
alpha_1 = cos(atan(diff(y_1)));alpha_1 = [alpha_1 alpha_1(end)];
alpha_2 = cos(atan(diff(y_2)));alpha_2 = [alpha_2 alpha_2(end)];
alpha_3 = cos(atan(diff(y_3)));alpha_3 = [alpha_3 alpha_3(end)];
alpha_4 = cos(atan(diff(y_4)));alpha_4 = [alpha_4 alpha_4(end)];
%%水平速度随x的变化%%
v_x1 = alpha_1.*sqrt(2*g*y_1);
v_x2 = alpha_2.*sqrt(2*g*y_2);
v_x3 = alpha_3.*sqrt(2*g*y_3);
v_x4 = alpha_4.*sqrt(2*g*y_4);
%%求出水平位移随时间变化%%
p = polyfit(x,1./v_x1,5)
plot(x,polyval(p,x))
a = polyval(p,x)
a(244)
int(1./v_x2,x)
int(1./v_x3,x)
int(1./v_x4,x)%%在这个地方报错%%

Answer 1

该回答引用于gpt与OKX安生共同编写：

• 该回答引用于gpt与OKX安生共同编写：

这个报错是因为被积函数比较复杂，Matlab 求解不定积分时无法直接给出解析解，需要使用数值积分来求解。

可以使用 MATLAB 的数值积分函数 integral 来求解这个积分。下面是对原代码做出修改，使用数值积分求解的结果：

clear
clc
syms t_0 %%定义在目标点的角度%%
g = 9.8;
x_0=10*rand(1)+5;
y_0=10*rand(1)+5;%%随机生成目标点
eqn =@(t_0)(t_0-sin(t_0))/(1-cos(t_0))-(x_0/y_0);
t_0 = fzero(eqn, pi);%摆线在目标点t的参数t_0%
r = y_0/(1-cos(t_0));%%解出摆线圆的半径%%
T = 0:0.005:t_0;
a_1 = -rand(1); %%二次函数二次项系数%%
x = r*(T-sin(T));
y_1 = r*(1-cos(T)); %%摆线%%
y_2 = a_1*(x.^2)+(y_0-a_1*x_0^2).*x/x_0;%%抛物线%%
y_3 = y_0*x/x_0;%%直线%%
y_4 = y_0*x./0.5.*(x<=0.5)+y_0*(x>0.5);%%折线%%
%%创建速度在x位置的水平分量对应的余弦值,延拓matlab求导缺少的元素%%
alpha_1 = cos(atan(diff(y_1)));alpha_1 = [alpha_1 alpha_1(end)];
alpha_2 = cos(atan(diff(y_2)));alpha_2 = [alpha_2 alpha_2(end)];
alpha_3 = cos(atan(diff(y_3)));alpha_3 = [alpha_3 alpha_3(end)];
alpha_4 = cos(atan(diff(y_4)));alpha_4 = [alpha_4 alpha_4(end)];
%%水平速度随x的变化%%
v_x1 = alpha_1.*sqrt(2*g*y_1);
v_x2 = alpha_2.*sqrt(2*g*y_2);
v_x3 = alpha_3.*sqrt(2*g*y_3);
v_x4 = alpha_4.*sqrt(2*g*y_4);

%%求出水平位移随时间变化%%
p = polyfit(x,1./v_x1,5);
a = @(t) integral(@(t) 1./interp1(x, 1./v_x1, interp1(y_1, T, t)), y_1(1), y_1(end));
a_1 = a(T);

% 画图
figure
plot(T, a_1)
title('水平位移随时间变化')
xlabel('时间')
ylabel('水平位移')

使用 integral 函数来求解复杂积分时，必须注意被积函数不应在被积区间内出现极点和不连续点，否则可能会得到错误的结果。

有用的话，还请采纳一下！

• 如果点集太多，可能会造成拟合效果不佳和计算时间过长，因此可以尝试通过以下几种方法解决：

1. 减少点集数量：可以在程序中将点集进行筛选，去除其中的一些冗余点，从而减少点集数量。

2. 尝试其他拟合方法：如果使用numpy.polyfit函数拟合效果不佳，可以尝试使用其他的拟合方法。例如，可以使用scipy.optimize.curve_fit函数拟合多项式曲线，该函数支持拟合多种函数形式。

3. 调整拟合参数：可以调整numpy.polyfit函数的拟合阶数和权重，以得到更好的拟合效果。例如，可以尝试使用更高阶的多项式拟合、增加权重分配、调整拟合间距等。

当输出向量元素为NAN时，可能是因为拟合过程中存在异常值或数据有误，建议先检查数据集是否存在问题，以及是否运用了适当的数据预处理方法。

关于a_1值只返回一个值的问题，可能是因为代码中只传递了一个参数，或者使用了只返回单个数值的函数。可以检查代码中的函数定义和调用语句，尝试为函数传递需要的参数，并确保函数能够输出正确数量的数值。