Error executing triple intergration equation. I am getting 'Not a Number (NaN)' as output for the triple integration equation attached below (Image file). I have mentioned the code used also. Kindly help me in knowing where I am wrong.

조회 수: 1 (최근 30일)
CODE:
syms x
syms y
syms phi
N=10;
u=4*pi*(10^-7)
a=0.045;
b=0.05;
fun=@(x,y,phi) (cos(phi).*(x.*y))./(sqrt((x.^2)+(y.^2)-(2.*x.*y.*cos(phi))))
q=integral3(fun,a,b,a,b,0,3.14)
L=q.*(u.*(N.^2))/((b-a).^2)
OUTPUT:
q = NaN
L = NaN
  댓글 수: 8
Walter Roberson
Walter Roberson 2019년 4월 15일
x=y is a singularity only at phi = 0.
x = -y would be a singularity at phi = Pi, but with the range of a and b values were are given, that does not happen.
At all other angles in (0, Pi) exclusive, the division by sqrt(x^2 + y^2 - 2*x*y*cos(theta)) does not lead to a singularity.
David Goodmanson
David Goodmanson 2019년 4월 15일
편집: David Goodmanson 2019년 4월 15일
Hi Torsten,
Not totally obvious, but if xvec and yvec are vectors expressed in polar coordinates (x,phix) and (y,phi) about the origin, this starts out as the integral
cos(phi-phix) / |yvec-xvec|
with area elements
xdx dphix and ydy dphi,
integrated over a circular annulus of inner and outer radii a and b. Toss in some constants and it's the self inductance of the annulus.
After doing the y and phi integrations, by symmetry the answer can't depend on phix. So you can assume that phix = 0 and multiply by 2pi later. This gives the integral in the original question, which has a singularity at phi=0, x=y as you mentioned. But cos(phi) = 1 there so the integrand is basically 1/|yvec-xvec|. Take xvec as the new origin and coordinates r and theta to describe yvec-xvec. For a small circle of radius B centered about xvec, this becomes
Integral (1/r) rdr dtheta = 2pi B
so the integral is bounded. But it takes a new set of coordinates to allow the new area element to cancel out the 1/r behavior. That's why I mentioned separate regions of integration before.

댓글을 달려면 로그인하십시오.

답변 (0개)

카테고리

Help CenterFile Exchange에서 Numeric Solvers에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by