Please help, I have no idea how to do this integral in matlab

1 次查看（过去 30 天）
Laura
Laura 2023-11-2
评论： Torsten 2023-11-3

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

采纳的回答

Torsten
Torsten 2023-11-2
编辑：Torsten 2023-11-3
The integral seems to exist for those values of theta0 with cos(theta0) ~= 0, thus theta0 ~= (2*k+1)*pi/2 for k in Z.
In this case, the Taylor expansion of sin(theta0)-sin(theta) around theta0 shows a behaviour like 1/sqrt(x) around x=0.
If cos(theta0) = 0, the Taylor expansion of sin(theta0)-sin(theta) around theta0 shows a behaviour like 1/x around x = 0 which means that the integral diverges.
  3 个评论
显示 1更早的评论
Torsten
Torsten 2023-11-2
编辑：Torsten 2023-11-2
Ran in:
Seems to work:
syms theta theta0
theta0 = 30.5*pi/180;
f = 1/sqrt(sin(theta0)-sin(theta));
sol = vpaintegral(f,theta,0,theta0)
sol = 
1.54036

请先登录，再进行评论。

更多回答（1 个）

Walter Roberson
Walter Roberson 2023-11-2
编辑：Walter Roberson 2023-11-3
Ran in:
syms a_0 positive
syms theta real; assumeAlso(theta >= 0);
theta_0 = sym(30.5);
F = 1/(sind(theta_0) - sind(theta))
F = 
a_0 = sym([5:5:25 26 27 28 29]).';
values = arrayfun(@(A) int(F, theta, 0, A), a_0);
values(1) %for example
ans = 
valued = double(values);
plot(a_0, valued)
  4 个评论
显示 2更早的评论
Walter Roberson
Walter Roberson 2023-11-3
Ran in:
@Torsten you are right, I did forget the sqrt() !
It looks like MATLAB is able to integrate even so, but with a more complicated formula.
The imaginary components of the results are so small that they are surely due to round-off error in the numeric calculations.
format long g
syms a_0 positive
syms theta real; assumeAlso(theta >= 0);
theta_0 = sym(30.5);
F = 1/sqrt(sind(theta_0) - sind(theta))
F = 
a_0 = sym([5:5:25 26 27 28 29 30:.1:30.4]).';
values = arrayfun(@(A) int(F, theta, 0, A), a_0);
values(1) %for example
ans = 
valued = double(values);
plot(a_0, real(valued), a_0, imag(valued))
legend({'real', 'imaginary'})
isimag = imag(valued) ~= 0;
a_0(isimag)
ans = 
valued(isimag)
ans =
7.34891328498187 + 1.60441675001006e-67i 15.5050816756634 + 1.20331256250754e-67i 24.7962226985339 - 6.85758772181718e-68i 35.8443076292128 - 3.0831418652267e-67i 50.1801693418716 - 4.50826319501674e-67i 53.7878830561832 - 4.49347594386457e-67i 57.8335596505002 - 4.34745183873693e-67i 62.5233721908558 - 4.02767753257133e-67i 68.307102046422 - 3.52121418060963e-67i 76.7289013342626 - 2.78739684218337e-67i 77.9450026132779 - 2.70237014805841e-67i 79.325691819995 - 2.60994982835737e-67i 80.9638641965888 - 2.5230747278384e-67i 83.0993245089765 - 2.43250281453138e-67i
Torsten
Torsten 2023-11-3
Ran in:
Your values seem to converge towards
syms theta theta0
theta0 = 30.5*pi/180;
f = 1/sqrt(sin(theta0)-sin(theta));
sol = vpaintegral(f,theta,0,theta0)
sol = 
1.54036
sol = sol * 180/pi
sol = 
vpa(sol)
ans = 
88.256285088573850083615798717037

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Particle & Nuclear Physics 的更多信息

标签

产品


版本

R2023a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by