Definite integral with a parameter in the bound

6 次查看(过去 30 天)
I'm trying to integrate the following function f = ((W-x)^2)*((x^2-1)^(1/2))*x where W is a parameter, from 1 to W with the following code:
syms x W
f = ((W-x)^2)*((x^2-1)^(1/2))*x;
F = int(f,x,1, W)
Now, I don't understand why it doesn't work. Thank you in advance for the help.
  1 个评论
John D'Errico
John D'Errico 2022-10-31
Actually, an analytical solution does exist. It just takes a little effort. Sometimes you need to push to get a symbolic solution.

请先登录,再进行评论。

采纳的回答

John D'Errico
John D'Errico 2022-10-31
编辑:John D'Errico 2022-10-31
syms X W
F = ((W-X)^2)*((X^2-1)^(1/2))*X;
int(F,X,1, W)
ans = 
As you found, int fails to generate a result. That does not mean a closed form solution does not exist. It does suggest that possibility though. It may be just that int is unable to solve the problem. In fact, I checked, and Wolfram Alpha also fails. Still that does not mean nothing can be done.
expand(F)
ans = 
And there we see the integral can be brokn into three terms. Are they all impossible to solve? Clearly not. The W is a constant, so I'll ignore that for now, just computing indefinite integrals.
int(X*sqrt(X^2 - 1),X)
ans = 
int(X^3*sqrt(X^2 - 1),X)
ans = 
It is an integral of this intermediate general form (with an even power of X outside) that seems to be the sticky part, at least for int.
int(X^2*sqrt(X^2 - 1),X)
ans = 
In fact though, Alpha actually was able to integrate a function of that general form. The complete mess, integrating F was just a bit too complicated for Alpha to digest.
A problem is, there is a term in there you won't like. There is an asin(x) term. But if x is GREATER than 1, the arc sine will be complex. In fact, Alpha is actually nice enough to show a plot of the result, and it does indeed have an imaginary part in the indefinite integral. So that may suggest there is no real solution after all.
Is that the case? Well, actually not. The INDEFINTE integral has an imaginary part. But that imaginary part appears to be constant for x above 1. And so the imaginary part will disappear when you form the definite integral, between bounds of 1 and W, as long as W>=1.
Of course, you can be lazy, and just compute a numerical solution.
f = @(x,w) ((w-x).^2).*((x.^2-1).^(1/2)).*x;
intf = @(w) integral(@(x) f(x,w),1,w)
intf = function_handle with value:
@(w)integral(@(x)f(x,w),1,w)
format long g
intf(2)
ans =
0.312068786948633
Anyway, with some effort, we could probably get MATLAB to compute the definite integral in an analytical form, since it appears an analytical solution does exist. Sometimes that requires finding a transformation that solve did not see. The obvious transformation here might be something like x = sin(u), but you will need to be careful, because this will now move into the complex plane, since x is greater than 1. The problem is doable. It just takes a little more effort at this point. Thus, if we make the desired transformation, we will have
x = sin(u)
dx = cos(u)du
int(x^2*sqrt(x^2-1)*dx) = -i*int(sin(u)^2*cos(u)^2*du)
the limits will change of course. Here we will be moving along the imaginary axis for the integral, where the limits of integration become [pi/2,asin(W)]. Remember, when W is greater than 1, asin(W) is complex.
syms u
I2 = -i*int(sin(u)^2*cos(u)^2)
I2 = 
Now we can try to put it all together. I think I've gotten all of my i's dotted at this point.
resultW = int((X^3 + W^2*X)*sqrt(X^2-1),X,[1,W]) - 2*W*(subs(I2,u,asin(W)) - subs(I2,u,pi/2))
resultW = 
Does this result agree with integral, when we did it in a numerical form?
intF = matlabFunction(resultW)
intF = function_handle with value:
@(W)(W.^2-1.0).^(3.0./2.0).*(W.^2.*4.0+1.0).*(2.0./1.5e+1)-W.*(pi.*6.25e-2i-asin(W).*1.25e-1i+W.*sqrt(-W.^2+1.0).*1.25e-1i-W.^3.*sqrt(-W.^2+1.0).*2.5e-1i).*2.0
intF(1) % should be 0
ans =
0
intF(2) % should be 0.3120687...
ans =
0.312068786948632
Yes, it does. That took a little effort, but it seeems to be working. As we hoped, all of the imaginary terms drop out.
  3 个评论
John D'Errico
John D'Errico 2022-10-31
编辑:John D'Errico 2022-10-31
I'm not really a beginner, though not as much of an expert as I want. (I doubt that will ever happen anyway.) But this should be sufficient, at least if int will manage a result.
syms X real
syms W real
assume(W > 1)
F = ((W-X)^2)*((X^2-1)^(1/2))*X;
int(F,X,1, W)
ans = 
And int is still stymied. Honestly, I DID try that before I went to the effort of splitting things apart. :)
The form that Alpha generates is I assume equally valid. The presence of a log is not amazing, as often these things have many ways the same expression can be written.
Antonino Roccaforte
Antonino Roccaforte 2022-10-31
Thank you very much for your time! It seems that in this case doing the calculation by hand is the better choice! I've tried with the trasformation x=cosh(t) in the case 1<x<W and t>0 and apparently it works!

请先登录,再进行评论。

更多回答(1 个)

Antonino Roccaforte
Antonino Roccaforte 2022-10-31
It seems a bit strange, but I've found that a little substitution in the piece that causes the major problem make int to produce the correct result. It is sufficient to put t=x+1;
syms x real
syms W real
F =-2*W*(x-1)^2*sqrt(x^2-2*x);
G=W^2*x*sqrt(x^2-1)+x^3*sqrt(x^2-1);
int(F,x,2,W+1)+int(G,x,1,W)
That's crazy.
  1 个评论
Paul
Paul 2022-11-1
I think this is similar, if not the same, letting Matlab do the work for the substitution.
syms X real
syms W real
assume(W > 1)
F = ((W-X)^2)*((X^2-1)^(1/2))*X;
syms U real
Fu = subs(F,X,U-1);
intF(W) = int(Fu,U,2,W+1)
intF(W) = 
intF(1)
ans = 
0
intF(2)
ans = 
vpa(intF(2))
ans = 
0.31206878694863289560703390535281

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Numbers and Precision 的更多信息

标签

产品


版本

R2022b

Community Treasure Hunt

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

Start Hunting!

Translated by