Numerically integrating a symbolic equation in a for loop?

Question

Harvey Rael 2018-6-25

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/407308-numerically-integrating-a-symbolic-equation-in-a-for-loop

评论： Torsten 2018-6-26

Hi! I'll try to explain qualitatively first what I am trying to do before then giving code. I have an equation that depends on a variable (theta) that changes (though how it changes doesn't matter). However, theta's end point depends on another variable N, and N is significant. I have a function that needs to be numerically integrated in terms of theta, to its end point which is pi/2+2piN. I want to change N with a for loop and store the outputs to be plotted.

Here's my code:

clc;clear;
% Define constants
AU = astroConstants(2);
mu = astroConstants(4);
% Define start and finish parameters for the exponential sinusoid.
r1  = AU;                 % Initial radius
psi = pi/2;               % Final polar angle of Mars/finish transfer
phi = pi/2;
r2 = 1.5*AU;              %Final radius
k2=0.2;
N =  0:1:5;
k1 = zeros(size(N));
y = zeros(size(N));
      for i = 1:length(N)
      k1(i) = sqrt( ( (log(r1/r2) + sin(k2*(psi + 2*pi.*N(i)))*tan(0)/k2) / (1-cos(k2*(psi+2*pi.*N(i)))) )^2 + 
      tan(0)^2/k2^2 );
      k0 = r1./exp(k1(i)*sin(phi));       %The parameters of the exponential sinusoid
      R = @(theta) k0*exp(k1(i)*sin(k2.*theta + phi));       %Equation of the exponential sinusoid
      theta_dot = @(theta) sqrt((mu./(R(theta).^3))./((tan(0))^2 + k1(i)*k2^2*sin(k2*theta + phi) + 1));
      z = @(theta)1./theta_dot(theta);        %change of variable name for ease of use
      y(i) = integral(z,0,(psi+2.*pi.*N(i))) - 94608000; %function to be plotted
end
The error message I keep getting is :

"Error in integral (line 88) Q = integralCalc(fun,a,b,opstruct);

Error in Untitled2 (line 23) y(i) = integral(theta_dot,0,(psi+2.*pi.*N(i))) - 94608000; %function to be plotted"

Anyone know where I am going wrong? When I have N as a single number, everything works fine.

EDIT : As @Torsten pointed out I forgot add in a couple of (i)'s in my code - so I have added them in. Error message persists

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Harvey Rael 2018-6-25

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/407308-numerically-integrating-a-symbolic-equation-in-a-for-loop#answer_326214

Fixed it - simply needed to flick the 'ArrayValued' setting in integral to true.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Torsten 2018-6-26

I don't see the reason why this should be necessary in your code.

Best wishes

Torsten.

请先登录，再进行评论。

Answer 2

Torsten 2018-6-25

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/407308-numerically-integrating-a-symbolic-equation-in-a-for-loop#answer_326187

In the definition of "R" and "thetadot", you use "k1" instead of "k1(i)".

6 个评论
显示 4更早的评论隐藏 4更早的评论

Harvey Rael 2018-6-25

Post your code then? Is y a vector or a single point solution? Because I contiually get errors

Torsten 2018-6-26

在 MATLAB Online 中打开

function main
AU = 2.0;
mu = 0.5;
% Define start and finish parameters for the exponential sinusoid.
r1  = AU;                 % Initial radius
psi = pi/2;               % Final polar angle of Mars/finish transfer
phi = pi/2;
r2 = 1.5*AU;              %Final radius
k2 = 0.2;
N =  0:1:5;
k1 = zeros(size(N));
y = zeros(size(N));
for i = 1:length(N)
      k1(i) = sqrt( ( (log(r1/r2) + sin(k2*(psi + 2*pi.*N(i)))*tan(0)/k2) / (1-cos(k2*(psi+2*pi.*N(i)))) )^2 +... 
      tan(0)^2/k2^2 );
      k0 = r1./exp(k1(i)*sin(phi));       %The parameters of the exponential sinusoid
      R = @(theta) k0*exp(k1(i)*sin(k2.*theta + phi));       %Equation of the exponential sinusoid
      theta_dot = @(theta) sqrt((mu./(R(theta).^3))./((tan(0))^2 + k1(i)*k2^2*sin(k2*theta + phi) + 1));
      z = @(theta)1./theta_dot(theta);        %change of variable name for ease of use
      y(i) = integral(z,0,(psi+2.*pi.*N(i))) - 94608000; %function to be plotted
end
end