Is your CAS up to this integration?
clicl...@freenet.de
I wonder if there is a Computer Algebra System that doesn't choke on
this double integral. Writing this in cylindrical coordinates makes
the inner (radial) integration elementary; still I couldn't proceed to
the outer (angular) integration. The result in dependence on the
parameters a and z should be an elementary expression that takes less
than 40 ASCII characters to write down (I have obtained this by
guessing/fitting results of numerical integration). For the
integration, it is important to restrict the parameter z to positive
values; note that in jxx and jxy, the zeroes of the denominator at
x=y=0 are then compensated by zeroes of the nominator.
jxx(x,y,z) := -x*y*(2*(x^2+y^2+z^2)^(3/2) - 2*z^3 - 3*z*(x^2+y^2)) /
((x^2+y^2+z^2)^(3/2)*(x^2+y^2)^2)
jxy(x,y,z) := ((x^2-y^2)*(x^2+y^2+z^2)*((x^2+y^2+z^2)^(1/2) - z) -
z*x^2*(x^2+y^2)) / ((x^2+y^2+z^2)^(3/2)*(x^2+y^2)^2)
int(int((jxx(r*cos(p)-a/2, r*SIN(p), z) * jxx(r*COS(p)+a/2, r*SIN(p),
z) + jxy(r*cos(p)-a/2, r*SIN(p), z) * jxy(r*COS(p)+a/2, r*SIN(p),
z))*r, r, 0, inf), p, 0, 2*pi)
z :epsilon Real [0, inf)
In cartesian coordinates, this would simply be
int(int(jxx(x-a/2, y, z) * jxx(x+a/2, y, z) + jxy(x-a/2, y, z) * jxy(x
+a/2, y, z), x, -inf, inf), y, -inf, inf)
z :epsilon Real [0, inf)
Martin.
clicl...@freenet.de
cliclic...@freenet.de schrieb:
From the lack of response over the course of one week I tentatively
conclude that, like Derive, none of the large commercial or free
systems can handle this double integral - at least not without an
unreasonable emount of operator assistance.
This is remarkable since this problem (from the area of
electrodynamics or electrical engineering) is of a size and kind where
one can expect computer algebra to be successful. It seems there is
something left to do in the field of automated symbolic integration.
For those wanting to have a check at the labors of their system, here
is my numerically confirmed guess at the solution (for z>0):
pi/a^2 * (8*z*(2*z^2+a^2) / (4*z^2+a^2)^(3/2) - 2)
Thanks for looking at the problem - and I am remaining interested to
learn about systems that can do it.
Martin.
Axel Vogt
Using Maple I was not able to get something in reasonable time.
However I used it to verify you solution for the first Taylor
coefficients in 'a'=0 up to degree 6 for the 'polar task' (it
is symmetric in +-a, hence 0 <= a, and odd degrees are all 0)
[assuming that COS=cos, SIN=sin was meant].
clicl...@freenet.de
Axel Vogt schrieb:
>
> However I used it to verify you solution for the first Taylor
> coefficients in 'a'=0 up to degree 6 for the 'polar task' (it
> is symmetric in +-a, hence 0 <= a, and odd degrees are all 0)
Thanks for this confirmation. I had also done the calculation for a=0
(this was comparatively straightforward in Derive, even using
cartesian coordinates), but I had not tried Taylor expansion.
Your wording suggests that the stumbling block for Maple is the outer
(polar) integration. What is the approximate size of Maple's result
after the inner (radial) integration? (I didn't even get there: Derive
seems to try too many transformations with log and arctan functions to
finish in acceptable time.)
> [assuming that COS=cos, SIN=sin was meant].
Right, SIN=sin and COS=cos. My apologies for the partial
capitalization.
Martin.
Axel Vogt
Then my wording was misleading, I meant the first double integral.
Maple12 did not return an answer for the inner integral in a time
I was willing to wait (say 10 min), even for the indefinite case.
A Taylor sum in a=0 up to degree=6 and the integral w.r.t. r is found
and has ~ 1000 characters as length. Integrating that w.r.t. p over
the interval is quick and easily reduced to the series in a for your
suggested solution.
I also tried to find a differential equation (in 'a') for your
solution (to see, whether it can be proven to hold by the integral),
but Maple did not give me one (a bit astonishing).
The last thing I tried was to write the 'cartesian' double integral
without roots in the denominator (which results in degrees of ~ 16
and splits in squares for the denominator). But it did not look as
if continuing by partial fractions would give me fun and an answer.
So I gave up (again).
clicl...@freenet.de
So for Maple too, the inner (radial) integration is difficult. I had
tried to feed Derive a broken-down inner integrand, but to no avail:
int(256*r*z^2*(16*a^2*r^4*(8*z^2-5*a^2-12*r^2)*cos(p)^4+r^2*(64*~
a^2*z^4+8*z^2*(5*a^4+8*a^2*r^2+16*r^4)-11*a^6+4*r^2*(65*a^4+28*a~
^2*r^2+48*r^4))*cos(p)^2+4*z^4*(a^4-24*a^2*r^2+16*r^4)+4*z^2*(a^~
2+4*r^2)*(a^4-18*a^2*r^2+8*r^4)+(a^2+4*r^2)^2*(a^4-13*a^2*r^2+4*~
r^4))/((4*a*r*cos(p)+4*r^2+4*z^2+a^2)^(3/2)*(-4*a*r*cos(p)+4*r^2~
+4*z^2+a^2)^(3/2)*(4*a*r*cos(p)+4*r^2+a^2)^2*(-4*a*r*cos(p)+4*r^~
2+a^2)^2),r)-int(64*r*z*(32*a^2*r^4*cos(p)^4-16*a*r^3*(3*a^2+4*r~
^2)*cos(p)^3+2*r^2*(16*a^2*z^2+5*a^4+8*a^2*r^2+16*r^4)*cos(p)^2+~
4*a^3*r*(20*r^2-a^2)*cos(p)+2*z^2*(a^4-24*a^2*r^2+16*r^4)+a^6-14~
*a^4*r^2-64*a^2*r^4+32*r^6)/((-4*a*r*cos(p)+4*z^2+a^2+4*r^2)^(3/~
2)*(4*a*r*cos(p)+4*r^2+a^2)^2*(-4*a*r*cos(p)+4*r^2+a^2)^2),r)-in~
t(64*r*z*(32*a^2*r^4*cos(p)^4+16*a*r^3*(3*a^2+4*r^2)*cos(p)^3+2*~
r^2*(16*a^2*z^2+5*a^4+8*a^2*r^2+16*r^4)*cos(p)^2+4*a^3*r*(a^2-20~
*r^2)*cos(p)+2*z^2*(a^4-24*a^2*r^2+16*r^4)+a^6-14*a^4*r^2-64*a^2~
*r^4+32*r^6)/((4*a*r*cos(p)+4*z^2+a^2+4*r^2)^(3/2)*(4*a*r*cos(p)~
+4*r^2+a^2)^2*(-4*a*r*cos(p)+4*r^2+a^2)^2),r)+int(16*r*(16*a^2*r~
^2*cos(p)^2+a^4-24*a^2*r^2+16*r^4)/((4*a*r*cos(p)+4*r^2+a^2)^2*(~
-4*a*r*cos(p)+4*r^2+a^2)^2),r)
The break-down makes it clear, however, that the innner (radial)
antiderivative is elementary, since the first integral can be
transformed from from dr to d(r^2)/2.
Martin.
Axel Vogt
Maple 12 can do these integrals, however the total result is more
600 pages long ... (and I did not succeed in simplifying it much).
My guess is: the integration routine has no good chance to do
algebraic simplifications in between (but I never looked closer
to the according theory & variants around Risch's algorithm).
Would be residual calculus a possible way? Or integral transforms
(as you say it comes from electrical engineering)?
PS: I do not understand the reasoning in your final remark.
clicl...@freenet.de
> 600 pages long ... (and I did not succeed in simplifying it much).
>
> My guess is: the integration routine has no good chance to do
> algebraic simplifications in between (but I never looked closer
> to the according theory & variants around Risch's algorithm).
>
> Would be residual calculus a possible way? Or integral transforms
> (as you say it comes from electrical engineering)?
>
> PS: I do not understand the reasoning in your final remark.
Wow! How easy it is to write a book these days. I had hoped for
something like 10 pages ... and I suppose it is hopeless to try an
angular integration with a 600-page integrand ... and there is the
issue of discontinuities at r=a/2 and perhaps elsewhere.
I will try to explain my remark that the radial antiderivative of the
rewritten integrand is evidently elementary. This fact was the sole
reason why I turned to polar coordinates:
Above I have broken down the integrand into four parts, where in each
the roots appear only in a common denominator. The last part is simply
rational in r; the 2nd and 3rd are of the type Rational(Root(r),r)
where the radicant of Root(r) is quadratic in r; and the first can be
brought into the same form because it is of the type
r*Rational(Root(r^2),r^2), where Root(r^2) is quadratic in r^2. Note
that in the 1st denominator expr1(r)^(3/2) expr2(r)^(3/2) = (expr1(r)
expr2(r))^(3/2), where expr1(r)*expr2(r) is quadratic in r^2.
So, Euler transformations will work for the first three integrands.
Derive can do the first and the last integral, but the 2nd and 3rd
never return.
I know that higher symmetries are at work (the 2-dimensional vector
field [jxx,jxy] has zero divergence, for example), but I thought: Hey,
this is the Age of Computer Algebra, so let's stop thinking: just
trust your CAS and churn it out. No Go.
Residual calculus looks promising here, only I would have to get used
to it again; some Systems should even be able to do it automatically.
But I have no doubts about my final result, which holds to better than
5-digit accuracy at every (a,z) I tried. Yet, as the integral appeared
to be good material to exercise the teeth of one's CAS on, I thought I
should throw it in here.
Yes, there is something left to do in the field of automated symbolic
integration.
Martin.
Axel Vogt
> Axel Vogt schrieb:
>> clicl...@freenet.de wrote:
>>> Axel Vogt schrieb:
>>>> clicl...@freenet.de wrote:
>>>>> Axel Vogt schrieb:
Ok, I should be able to understand that in patients. More or less.
Residual calculus ... for refreshing that stuff to be forgotten
after month again I use Cartan's little book on complex analysis
"Elementare Theorien .." or Jähnich "Einführung i d FktTheorie":
integral = Sum(residues) if vanishing in infinity is known.
Sometimes Maple seems to use it, but may be it is to complicated
as algorithmic recipe in general (brunch cuts & points, limits,
and handling parameters in this case).
Hat Spass gemacht ... Gruss, Axel
Vladimir Bondarenko
> Axel Vogt schrieb:
>
>
>
> > cliclic...@freenet.de wrote:
> > > Axel Vogt schrieb:
M> "this is the Age of Computer Algebra"
No.
In earnest, to bring some reasonable quality rather than
lots of garbage, to start, the Age of Computer Algebra
needs critically a fully automated QA essence which is
able to help building large scale systems.
Like the VM machine is.
As you could guess, there are news coming from us soon.
Best wishes,
Vladimir Bondarenko
Co-founder, CEO, Mathematical Director
http://www.cybertester.com/ Cyber Tester Ltd.
--------------------------------------------------------
"We must understand that technologies
like these are the way of the future."
--------------------------------------------------------
Jacko
My particular failing with CAS is the solution of differential
equation sets.
m'' = -(f g)^2 m
m Pi (x' g)^3 = -2 h x''
g = ((1-(x'/c)^2)(1 - k m/(x c^2)))^(1/2)
Simply I want x and m as functions of t (the w.r.t var)
Mathematica 5.1 complains of overspecification, but I don't think it
is, and each on there own, the second one does not do a full integral,
and returns many secondary variables with no expresssions for them.
Known constants are h, but f is not know to be constant, but probably
is in local space.
cheers
jacko
p.s. Its to do with quantum particle geometry. ref. uncertain geometry
- Simon Jackson, BEng. (c)2008
clicl...@freenet.de
Jacko schrieb:
>
> My particular failing with CAS is the solution of differential
> equation sets.
>
> m'' = -(f g)^2 m
> m Pi (x' g)^3 = -2 h x''
> g = ((1-(x'/c)^2)(1 - k m/(x c^2)))^(1/2)
>
> Simply I want x and m as functions of t (the w.r.t var)
>
> Mathematica 5.1 complains of overspecification, but I don't think it
> is, and each on there own, the second one does not do a full integral,
> and returns many secondary variables with no expresssions for them.
>
One really doesn't have to be a Mathematica developer to respond to
this. I would like to do so with a rhetorical question:
Can you give an example of a system of three independent ordinary
differential equations for two functions that is not overdetermined?
Martin.
Jacko
> Jacko schrieb:
>
>
>
> > My particular failing with CAS is the solution of differential
> > equation sets.
>
> > m'' = -(f g)^2 m
> > m Pi (x' g)^3 = -2 h x''
> > g = ((1-(x'/c)^2)(1 - k m/(x c^2)))^(1/2)
g = ((1-(x'/c)^2)(1 - k j(m,x)))^(1/2)
such that g displaced to origin is gravity at displacement x
> > Simply I want x and m as functions of t (the w.r.t var)
>
> > Mathematica 5.1 complains of overspecification, but I don't think it
> > is, and each on there own, the second one does not do a full integral,
> > and returns many secondary variables with no expresssions for them.
>
> One really doesn't have to be a Mathematica developer to respond to
> this. I would like to do so with a rhetorical question:
>
> Can you give an example of a system of three independent ordinary
> differential equations for two functions that is not overdetermined?
>
> Martin.
Well silly me, but I thought I would work out that g needs to
substitute before eval. But then it still does not give a fully
defined answer for x(t). Speculation on j(m,x)is by use of
differential quadrapole to be including x-x ln[x] to avoid a pole
infinity problem.
Well cheers for that
jacko
clicl...@freenet.de
Axel Vogt schrieb:
>
> Maple 12 can do these integrals, however the total result is more
> 600 pages long ... (and I did not succeed in simplifying it much).
>
I have succeeded in simplifying my integrand somewhat, using
integration by parts (in principle, a CAS could have managed to do
this too ...). In terms of:
hxz(x,y,z) := 3*x*z/(x^2+y^2+z^2)^(5/2)
vxz(x,y,z) := x*(sqrt(x^2+y^2+z^2) - z)/((x^2+y^2)*sqrt(x^2+y^2+z^2))
the integral can be written in polar coordinates as:
int(int(1/2*(hxz(r*cos(p)-a/2, r*sin(p), z) * vxz(r*cos(p)+a/2,
r*sin(p), z) + vxz(r*cos(p)-a/2, r*sin(p), z) * hxz(r*cos(p)+a/2,
r*sin(p), z))*r, r, 0, inf), p, 0, 2*pi)
and in cartesian coordinates as:
int(int(1/2*(hxz(x-a/2, y, z) * vxz(x+a/2, y, z) + vxz(x-a/2, y, z) *
hxz(x+a/2, y, z)), x, -inf, inf), y, -inf, inf)
The radial integration for the simplified version can be performed by
Derive; the integration is again elementary as the new integrand
splits as follows:
int(192*r*z^2*(a^2-4*c^2*r^2)*(16*z^4*(a^2+4*r^2)+8*z^2*(a^2+4*a~
*c*r+4*r^2)*(a^2-4*a*c*r+4*r^2)+(a^2+4*r^2)*(a^2+4*a*c*r+4*r^2)*~
(a^2-4*a*c*r+4*r^2))/((a^2+4*a*c*r+4*r^2)*(a^2-4*a*c*r+4*r^2)*(4~
*z^2+a^2+4*a*c*r+4*r^2)^(5/2)*(4*z^2+a^2-4*a*c*r+4*r^2)^(5/2)),r~
)+int(48*r*z*(4*c^2*r^2-a^2)/((a^2-4*a*c*r+4*r^2)*(4*z^2+a^2+4*a~
*c*r+4*r^2)^(5/2)),r)+int(48*r*z*(4*c^2*r^2-a^2)/((a^2+4*a*c*r+4~
*r^2)*(4*z^2+a^2-4*a*c*r+4*r^2)^(5/2)),r)
where I have abbreviated c := cos(p).
The result of the integration from r=0 to r=inf is
-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(-2*z*sqrt(z^2+~
a^2*c*(c-sqrt(c^2-1)))*(4*z^2+a^2)+8*z^4-2*a^2*z^2*(2*c*sqrt(c^2~
-1)-2*c^2-1)+a^4*c*(sqrt(c^2-1)*(2*c^2-1)-2*c^3+2*c))/(16*(z^2+a~
^2*c*(c-sqrt(c^2-1)))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*~
c^4-c^2-1)*ln(-2*z*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c))*(4*z^2+a^2)+8~
*z^4+2*a^2*z^2*(2*c*sqrt(c^2-1)+2*c^2+1)-a^4*c*(sqrt(c^2-1)*(2*c~
^2-1)+2*c^3-2*c))/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))+3*a^2*z~
*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(-2*sqrt(z^2+a^2*c*(c-s~
qrt(c^2-1)))*sqrt(4*z^2+a^2)+4*z^2-a^2*(c*sqrt(c^2-1)-c^2-1))/(8~
*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^~
2+1)+4*c^4-c^2-1)*ln(-2*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c))*sqrt(4*z~
^2+a^2)+4*z^2+a^2*(c*sqrt(c^2-1)+c^2+1))/(8*(z^2+a^2*c*(sqrt(c^2~
-1)+c))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(~
-2*z*sqrt(z^2+a^2*c*(c-sqrt(c^2-1)))+2*z^2-a^2*c*sqrt(c^2-1))/(1~
6*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c~
^2+1)-4*c^4+c^2+1)*ln(a*(sqrt(c^2-1)-2*c)-2*sqrt(z^2+a^2*c*(c-sq~
rt(c^2-1))))/(8*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*sq~
rt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1)*ln(-2*z*sqrt(z^2+a^2*c*(sqrt(c^~
2-1)+c))+2*z^2+a^2*c*sqrt(c^2-1))/(16*(z^2+a^2*c*(sqrt(c^2-1)+c)~
)^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1)*ln(-2*sqr~
t(z^2+a^2*c*(sqrt(c^2-1)+c))-a*(sqrt(c^2-1)+2*c))/(8*(z^2+a^2*c*~
(sqrt(c^2-1)+c))^(5/2))+3*a^2*z*((z^2+a^2*c*(sqrt(c^2-1)+c))^(5/~
2)*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)-(z^2+a^2*c*(c-sqrt(c^2-~
1)))^(5/2)*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1))*ln(a^2)/(16*(z~
^2+a^2*c*(c-sqrt(c^2-1)))^(5/2)*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2~
))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(a*(c-sqrt(c^~
2-1)))/(8*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))+3*a^2*z*(c*sqrt(c^2~
-1)*(4*c^2+1)+4*c^4-c^2-1)*ln(a*(sqrt(c^2-1)+c))/(8*(z^2+a^2*c*(~
sqrt(c^2-1)+c))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^~
2-1)*ln(-2*c*sqrt(c^2-1)-2*c^2+1)/(16*(z^2+a^2*c*(sqrt(c^2-1)+c)~
)^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(2*c*sq~
rt(c^2-1)-2*c^2+1)/(16*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-z*(sqr~
t(4*z^2+a^2)*(96*c^2*z^11+128*a*c^3*z^10-24*a^2*z^9*(c^4-4*c^2-3~
)-16*a^3*c*z^8*(28*c^4-27*c^2-15)+2*a^4*z^7*(108*c^6-101*c^4-22*~
c^2+15)-4*a^5*c*z^6*(28*c^6+127*c^4-111*c^2-56)+a^6*z^5*(262*c^6~
-249*c^4-64*c^2+3)-a^7*c*z^4*(140*c^6+187*c^4-198*c^2-65)+8*a^8*~
c^2*z^3*(c^2-1)*(14*c^2+3)-2*a^9*c*z^2*(28*c^6+7*c^4-24*c^2-3)+a~
^10*c^2*z*(c^2-1)*(13*c^2+3)+a^11*c^3*(1-c^2)*(7*c^2+5))-256*c^2~
*z^12+64*a^2*z^10*(5*c^4-6*c^2-4)-16*a^4*z^8*(58*c^6-71*c^4-6*c^~
2+11)-4*a^6*z^6*(332*c^6-303*c^4-83*c^2+10)-a^8*z^4*(698*c^6-493~
*c^4-176*c^2+3)-2*a^10*c^2*z^2*(71*c^4-44*c^2-19)+3*a^12*c^2*(1-~
c^2)*(3*c^2+1))/(2*(z^2+a^2*c*(c-sqrt(c^2-1)))^2*(z^2+a^2*c*(sqr~
t(c^2-1)+c))^2*(4*z^2+a^2)^(3/2)*(4*z^2+a^2*(1-c^2))^2)
This is clearly below 600 pages. However, Derive cannot handle the
angular integration of the logarithmic terms. What would be a good
integration technique for those?
Of course, the final result for z>0 should be once more:
pi/a^2 * (8*z*(2*z^2+a^2) / (4*z^2+a^2)^(3/2) - 2)
Martin.
G. A. Edgar
<803a5932-6a21-46ee...@o2g2000yqd.googlegroups.com>,
<clicl...@freenet.de> wrote:
> Axel Vogt schrieb:
> >
> > Maple 12 can do these integrals, however the total result is more
> > 600 pages long ... (and I did not succeed in simplifying it much).
> >
>
> I have succeeded in simplifying my integrand somewhat, using
> integration by parts (in principle, a CAS could have managed to do
> this too ...). In terms of:
>
> hxz(x,y,z) := 3*x*z/(x^2+y^2+z^2)^(5/2)
>
> vxz(x,y,z) := x*(sqrt(x^2+y^2+z^2) - z)/((x^2+y^2)*sqrt(x^2+y^2+z^2))
>
> the integral can be written in polar coordinates as:
>
> int(int(1/2*(hxz(r*cos(p)-a/2, r*sin(p), z) * vxz(r*cos(p)+a/2,
> r*sin(p), z) + vxz(r*cos(p)-a/2, r*sin(p), z) * hxz(r*cos(p)+a/2,
> r*sin(p), z))*r, r, 0, inf), p, 0, 2*pi)
>
> and in cartesian coordinates as:
>
> int(int(1/2*(hxz(x-a/2, y, z) * vxz(x+a/2, y, z) + vxz(x-a/2, y, z) *
> hxz(x+a/2, y, z)), x, -inf, inf), y, -inf, inf)
>
> [...]
> >
> Of course, the final result for z>0 should be once more:
>
> pi/a^2 * (8*z*(2*z^2+a^2) / (4*z^2+a^2)^(3/2) - 2)
>
> Martin.
>
For a simpler problem, consider the limit for a=0 ...
The answer should be pi/(4*z^2), the problem is:
Int(Int(3*x^2*z*((x^2+y^2+z^2)^(1/2)-z)/((x^2+y^2+z^2)^3*(x^2+y^2)),
x = -infinity .. infinity), y = -infinity .. infinity)
--
G. A. Edgar http://www.math.ohio-state.edu/~edgar/
Axel Vogt
For a=0 Maple12 does it, "Int(Int(P, r=0..infinity),p=0..2*Pi);
subs(a=0,%); value(%);" gives 1/4*Pi/z^2, P in the transformed
coordinates
PS: Maple does not see it is the desired solution for a=0, one
has to do l'Hospital by hand for that.
G. A. Edgar
<&nor...@axelvogt.de> wrote:
> G. A. Edgar wrote:
> >
> > For a simpler problem, consider the limit for a=0 ...
> > The answer should be pi/(4*z^2), the problem is:
> >
> > Int(Int(3*x^2*z*((x^2+y^2+z^2)^(1/2)-z)/((x^2+y^2+z^2)^3*(x^2+y^2)),
> > x = -infinity .. infinity), y = -infinity .. infinity)
>
> For a=0 Maple12 does it, "Int(Int(P, r=0..infinity),p=0..2*Pi);
> subs(a=0,%); value(%);" gives 1/4*Pi/z^2, P in the transformed
> coordinates
In Cartesian coordinates, Maple goes a while then says too many levels
...
>
> PS: Maple does not see it is the desired solution for a=0, one
> has to do l'Hospital by hand for that.
I just did limit using the formula given.
Axel Vogt
More general the formula is true for the first Taylor terms
in a=0, cf my first post in that thread, 05 Nov 2008.
Axel Vogt
>
> This is clearly below 600 pages. However, Derive cannot handle the
> angular integration of the logarithmic terms. What would be a good
> integration technique for those?
>
> Of course, the final result for z>0 should be once more:
>
> pi/a^2 * (8*z*(2*z^2+a^2) / (4*z^2+a^2)^(3/2) - 2)
>
> Martin.
For Maple that did not help (some solution 'is' found by Maple
for 0<r and still of incredible length - while your last task
is not solved too).
My only 'idea' for proving was to replace 'a' by 'sqrt(a)', so
your solution solves the differential equation -4*Pi*z^2+2*Pi*a
+ 2*(8*z^4+6*a*z^2+a^2)*a*D(h)(a)+(10*a*z^2+3*a^2+16*z^4)*h(a)
(using 0 < z). But I am unable to see that the double integral
satisfies it, h(0) = Pi/z^2/4.
And a trivial observation: it is symmetric to p = Pi (so you
only need half the interval of the outer integral).
But I have to give up, again.
clicl...@freenet.de
G. A. Edgar schrieb:
> In article <6oogjeF...@mid.individual.net>, Axel Vogt
> <&nor...@axelvogt.de> wrote:
>
> > G. A. Edgar wrote:
> > >
> > > For a simpler problem, consider the limit for a=0 ...
> > > The answer should be pi/(4*z^2), the problem is:
> > >
> > > Int(Int(3*x^2*z*((x^2+y^2+z^2)^(1/2)-z)/((x^2+y^2+z^2)^3*(x^2+y^2)),
> > > x = -infinity .. infinity), y = -infinity .. infinity)
> >
> > For a=0 Maple12 does it, "Int(Int(P, r=0..infinity),p=0..2*Pi);
> > subs(a=0,%); value(%);" gives 1/4*Pi/z^2, P in the transformed
> > coordinates
>
> In Cartesian coordinates, Maple goes a while then says too many levels
> ...
>
If you are still talking about the case a=0, Maple's must be a
particularly bad Risch integrator! With the elementary approach of
Derive, simplification of the double integral
int(int(3*x^2*z*(sqrt(x^2+y^2+z^2) - z) / ((x^2+y^2)*(x^2+y^2+z^2)^3),
x, -inf, inf), y, -inf, inf)
immediately gives pi/(4*z^2).
Martin.
clicl...@freenet.de
cliclic...@freenet.de schrieb:
>
> I have succeeded in simplifying my integrand somewhat, using
> integration by parts (in principle, a CAS could have managed to do
> this too ...). In terms of:
>
> hxz(x,y,z) := 3*x*z/(x^2+y^2+z^2)^(5/2)
>
> vxz(x,y,z) := x*(sqrt(x^2+y^2+z^2) - z)/((x^2+y^2)*sqrt(x^2+y^2+z^2))
>
> the integral can be written in polar coordinates as:
>
> int(int(1/2*(hxz(r*cos(p)-a/2, r*sin(p), z) * vxz(r*cos(p)+a/2,
> r*sin(p), z) + vxz(r*cos(p)-a/2, r*sin(p), z) * hxz(r*cos(p)+a/2,
> r*sin(p), z))*r, r, 0, inf), p, 0, 2*pi)
>
> and in cartesian coordinates as:
>
> int(int(1/2*(hxz(x-a/2, y, z) * vxz(x+a/2, y, z) + vxz(x-a/2, y, z) *
> hxz(x+a/2, y, z)), x, -inf, inf), y, -inf, inf)
>
> The radial integration for the simplified version can be performed by
> Derive; the integration is again elementary as the new integrand
> splits as follows:
>
> int(192*r*z^2*(a^2-4*c^2*r^2)*(16*z^4*(a^2+4*r^2)+8*z^2*(a^2+4*a~
> *c*r+4*r^2)*(a^2-4*a*c*r+4*r^2)+(a^2+4*r^2)*(a^2+4*a*c*r+4*r^2)*~
> (a^2-4*a*c*r+4*r^2))/((a^2+4*a*c*r+4*r^2)*(a^2-4*a*c*r+4*r^2)*(4~
> *z^2+a^2+4*a*c*r+4*r^2)^(5/2)*(4*z^2+a^2-4*a*c*r+4*r^2)^(5/2)),r~
> )+int(48*r*z*(4*c^2*r^2-a^2)/((a^2-4*a*c*r+4*r^2)*(4*z^2+a^2+4*a~
> *c*r+4*r^2)^(5/2)),r)+int(48*r*z*(4*c^2*r^2-a^2)/((a^2+4*a*c*r+4~
> *r^2)*(4*z^2+a^2-4*a*c*r+4*r^2)^(5/2)),r)
>
> where I have abbreviated c := cos(p).
>
> The result of the integration from r=0 to r=inf is
>
This result was not entirely correct: I had done the first and third
integrals separately and then doubled the third to account for the
second one. Instead, I should have replaced c by -c in the result for
the third integral and then added this to account for the second one.
The result after integration from p=0 to p=2*pi will be unaffected by
this. Still, as somebody may want to check the intermediate output of
his CAS, here is the corrected result of the integration from r=0 to
r=inf:
-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(-2*z*sqrt(z^2+~
a^2*c*(c-sqrt(c^2-1)))*(4*z^2+a^2)+8*z^4-2*a^2*z^2*(2*c*sqrt(c^2~
-1)-2*c^2-1)+a^4*c*(sqrt(c^2-1)*(2*c^2-1)-2*c^3+2*c))/(16*(z^2+a~
^2*c*(c-sqrt(c^2-1)))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*~
c^4-c^2-1)*ln(-2*z*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c))*(4*z^2+a^2)+8~
*z^4+2*a^2*z^2*(2*c*sqrt(c^2-1)+2*c^2+1)-a^4*c*(sqrt(c^2-1)*(2*c~
^2-1)+2*c^3-2*c))/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))+3*a^2*z~
*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(-2*sqrt(z^2+a^2*c*(c-s~
qrt(c^2-1)))*sqrt(4*z^2+a^2)+4*z^2-a^2*(c*sqrt(c^2-1)-c^2-1))/(8~
*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^~
2+1)+4*c^4-c^2-1)*ln(-2*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c))*sqrt(4*z~
^2+a^2)+4*z^2+a^2*(c*sqrt(c^2-1)+c^2+1))/(8*(z^2+a^2*c*(sqrt(c^2~
-1)+c))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(~
-2*z*sqrt(z^2+a^2*c*(c-sqrt(c^2-1)))+2*z^2-a^2*c*sqrt(c^2-1))/(1~
6*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c~
^2+1)-4*c^4+c^2+1)*ln(a*(sqrt(c^2-1)-2*c)-2*sqrt(z^2+a^2*c*(c-sq~
rt(c^2-1))))/(16*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-3*a^2*z*(c*s~
qrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(a*(2*c-sqrt(c^2-1))-2*sqrt(~
z^2+a^2*c*(c-sqrt(c^2-1))))/(16*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2~
))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1)*ln(-2*z*sqrt(z^~
2+a^2*c*(sqrt(c^2-1)+c))+2*z^2+a^2*c*sqrt(c^2-1))/(16*(z^2+a^2*c~
*(sqrt(c^2-1)+c))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-~
c^2-1)*ln(a*(sqrt(c^2-1)+2*c)-2*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c)))~
/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(~
4*c^2+1)+4*c^4-c^2-1)*ln(-2*sqrt(z^2+a^2*c*(sqrt(c^2-1)+c))-a*(s~
qrt(c^2-1)+2*c))/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))+3*a^2*z*~
((z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2)*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^~
4+c^2+1)-(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2)*(c*sqrt(c^2-1)*(4*c^2~
+1)+4*c^4-c^2-1))*ln(a^2)/(16*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2)*~
(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2~
+1)-4*c^4+c^2+1)*ln(a*(c-sqrt(c^2-1)))/(16*(z^2+a^2*c*(c-sqrt(c^~
2-1)))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1)*ln(-~
a*(sqrt(c^2-1)+c))/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2))-3*a^2*~
z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(a*(sqrt(c^2-1)-c))/(1~
6*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))+3*a^2*z*(c*sqrt(c^2-1)*(4*c~
^2+1)+4*c^4-c^2-1)*ln(a*(sqrt(c^2-1)+c))/(16*(z^2+a^2*c*(sqrt(c^~
2-1)+c))^(5/2))-3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)+4*c^4-c^2-1)*ln~
(-2*c*sqrt(c^2-1)-2*c^2+1)/(16*(z^2+a^2*c*(sqrt(c^2-1)+c))^(5/2)~
)+3*a^2*z*(c*sqrt(c^2-1)*(4*c^2+1)-4*c^4+c^2+1)*ln(2*c*sqrt(c^2-~
1)-2*c^2+1)/(16*(z^2+a^2*c*(c-sqrt(c^2-1)))^(5/2))-z*(z*sqrt(4*z~
^2+a^2)*(96*c^2*z^10-24*a^2*z^8*(c^4-4*c^2-3)+2*a^4*z^6*(108*c^6~
-101*c^4-22*c^2+15)+a^6*z^4*(262*c^6-249*c^4-64*c^2+3)+8*a^8*c^2~
*z^2*(c^2-1)*(14*c^2+3)+a^10*c^2*(c^2-1)*(13*c^2+3))-256*c^2*z^1~
2+64*a^2*z^10*(5*c^4-6*c^2-4)-16*a^4*z^8*(58*c^6-71*c^4-6*c^2+11~
)-4*a^6*z^6*(332*c^6-303*c^4-83*c^2+10)-a^8*z^4*(698*c^6-493*c^4~
-176*c^2+3)-2*a^10*c^2*z^2*(71*c^4-44*c^2-19)-3*a^12*c^2*(3*c^4-~
2*c^2-1))/(2*(z^2+a^2*c*(c-sqrt(c^2-1)))^2*(z^2+a^2*c*(sqrt(c^2-~
1)+c))^2*(4*z^2+a^2)^(3/2)*(4*z^2+a^2*(1-c^2))^2)
Martin.
G. A. Edgar
> > In Cartesian coordinates, Maple goes a while then says too many levels
> > ...
> >
>
> If you are still talking about the case a=0, Maple's must be a
> particularly bad Risch integrator! With the elementary approach of
> Derive, simplification of the double integral
>
> int(int(3*x^2*z*(sqrt(x^2+y^2+z^2) - z) / ((x^2+y^2)*(x^2+y^2+z^2)^3),
> x, -inf, inf), y, -inf, inf)
>
> immediately gives pi/(4*z^2).
>
> Martin.
Yes, that's the one. In Maple if you evaluate it leaving the outer
integral indefinite, then take the limits at y=infinity
and y=-infinity and subtract, then it works out fine. But just
giving it the expression shown yields: "too may levels of recursion".
So for the definite integral it must be trying something else!
clicl...@freenet.de
> for 0<r and still of incredible length - while your last task
> is not solved too).
>
Current "Risch" integrators appear to be useless for serious
integration work! The fact that USD-5000 systems are of little help
with my problem has spurned me on. It turns out (this wasn't obvious
to me) that the angular integration task is elementary too, and that
the antiderivative can be found with Derive.
To facilitate progress, I have split the logarithms present in the
result of the integration from r=0 to r=inf (the form given earlier
was simply that produced by the Derive integrator, limit extractor,
and simplifier), and expressed c = cos(p), sqrt(c^2-1) = #i*|sin(p)|
for 0 < p < pi in terms of t := exp(#i*p):
3*sqrt(2)*a^2*t^2*z*(t^2+3)/(8*(2*z^2+a^2*(t^2+1))^(5/2))*(-ln((~
sqrt(2)*t*sqrt(2*z^2+a^2*(t^2+1))-2*t*z+a*(t^2+1))/t)-ln((sqrt(2~
)*t*sqrt(2*z^2+a^2*(t^2+1))-2*t*z-a*(t^2+1))/t)-2*ln(sqrt(2)*sqr~
t(4*z^2+a^2)-sqrt(2*z^2+a^2*(t^2+1))-sqrt(2)*z)-2*ln(sqrt(2)*sqr~
t(4*z^2+a^2)-sqrt(2*z^2+a^2*(t^2+1))+sqrt(2)*z)+ln(sqrt(2*z^2+a^~
2*(t^2+1))-sqrt(2)*z-sqrt(2)*a*t)+ln(sqrt(2*z^2+a^2*(t^2+1))-sqr~
t(2)*z+sqrt(2)*a*t)+2*ln(sqrt(2*z^2+a^2*(t^2+1))-sqrt(2)*z)+ln((~
4*t*z+a*(1-t^2))/t)+ln((4*t*z+a*(t^2-1))/t))+(-3*sqrt(2)*a^2*t^2~
*z*(3*t^2+1)*sqrt((2*t^2*z^2+a^2*(t^2+1))/t^2)/(8*(2*t^2*z^2+a^2~
*(t^2+1))^3))*(ln((sqrt(2)*t*sqrt((2*t^2*z^2+a^2*(t^2+1))/t^2)-2~
*t*z+a*(t^2+1))/t)+ln((sqrt(2)*t*sqrt((2*t^2*z^2+a^2*(t^2+1))/t^~
2)-2*t*z-a*(t^2+1))/t)-ln((t*sqrt((2*t^2*z^2+a^2*(t^2+1))/t^2)-s~
qrt(2)*t*z-sqrt(2)*a)/t)-ln((t*sqrt((2*t^2*z^2+a^2*(t^2+1))/t^2)~
-sqrt(2)*t*z+sqrt(2)*a)/t)+2*ln(sqrt(2)*(sqrt(4*z^2+a^2)-z)-sqrt~
((2*t^2*z^2+a^2*(t^2+1))/t^2))+2*ln(sqrt(2)*(sqrt(4*z^2+a^2)+z)-~
sqrt((2*t^2*z^2+a^2*(t^2+1))/t^2))-2*ln(sqrt((2*t^2*z^2+a^2*(t^2~
+1))/t^2)-sqrt(2)*z)-ln((4*t*z+a*(1-t^2))/t)-ln((4*t*z+a*(t^2-1)~
)/t))-2*t^2*z*(z*sqrt(4*z^2+a^2)*(1536*t^4*z^10*(t^2+1)^2-96*a^2~
*t^2*z^8*(t^8-12*t^6-74*t^4-12*t^2+1)+8*a^4*z^6*(27*t^12+61*t^10~
-87*t^8-2*t^6-87*t^4+61*t^2+27)+2*a^6*z^4*(131*t^12+288*t^10-539~
*t^8-1296*t^6-539*t^4+288*t^2+131)+16*a^8*z^2*(t^2+1)^2*(t^4-2*t~
^2+1)*(7*t^4+20*t^2+7)+a^10*(t^2+1)^2*(t^4-2*t^2+1)*(13*t^4+38*t~
^2+13))-4096*t^4*z^12*(t^2+1)^2+256*a^2*t^2*z^10*(5*t^8-4*t^6-82~
*t^4-4*t^2+5)-32*a^4*z^8*(29*t^12+32*t^10-181*t^8-16*t^6-181*t^4~
+32*t^2+29)-16*a^6*z^6*(83*t^12+195*t^10-299*t^8-662*t^6-299*t^4~
+195*t^2+83)-2*a^8*z^4*(349*t^12+1108*t^10-117*t^8-1656*t^6-117*~
t^4+1108*t^2+349)-2*a^10*z^2*(t^2+1)^2*(71*t^8+108*t^6-230*t^4+1~
08*t^2+71)-3*a^12*(t^2+1)^2*(3*t^8+4*t^6-14*t^4+4*t^2+3))/((2*z^~
2+a^2*(t^2+1))^2*(4*z^2+a^2)^(3/2)*(2*t^2*z^2+a^2*(t^2+1))^2*(16~
*t^2*z^2-a^2*(t^4-2*t^2+1))^2)
Care was taken not to introduce discontinuities along the integration
path from p=0 to p=pi. The last term is rational in t and its
integration (dp = -#i/t dt) therefore straightforward. There remain
9+9 logarithmic terms, the second group of which is the complex
conjugate of the first. The antiderivative of the first term, for
example (this and the second one are the hardest for Derive),
int(-3*sqrt(2)*a^2*t^2*z*(t^2+3)/(8*(2*z^2+a^2*(t^2+1))^(5/2))*l~
n((sqrt(2)*t*sqrt(2*z^2+a^2*(t^2+1))-2*t*z+a*(t^2+1))/t)*(-#i/t)~
,t)
is then found to be:
#i*(z*(43*z^6+63*a^2*z^4+33*a^4*z^2+5*a^6)*atan(a*t/sqrt(2*z^2+a~
^2))/(8*a^2*sqrt(2*z^2+a^2)*(3*z^2+a^2)^3)+z*(sqrt(4*z^2+a^2)+z)~
^3*(3*z*sqrt(4*z^2+a^2)-7*z^2-2*a^2)*ln(sqrt(2)*(sqrt(4*z^2+a^2)~
-z)*sqrt(2*z^2+a^2*(t^2+1))+a*t*sqrt(4*z^2+a^2)+2*z^2-2*a*t*z+a^~
2)/(8*a^2*(3*z^2+a^2)^3)+z*(sqrt(4*z^2+a^2)-z)^3*(3*z*sqrt(4*z^2~
+a^2)+7*z^2+2*a^2)*ln(sqrt(2)*(sqrt(4*z^2+a^2)+z)*sqrt(2*z^2+a^2~
*(t^2+1))+a*t*sqrt(4*z^2+a^2)-2*z^2+2*a*t*z-a^2)/(8*a^2*(3*z^2+a~
^2)^3)-sqrt(2)*z*(4*z^2+a^2*(3*t^2+5))*ln((sqrt(2)*t*sqrt(2*z^2+~
a^2*(t^2+1))-2*t*z+a*(t^2+1))/t)/(8*a^2*(2*z^2+a^2*(t^2+1))^(3/2~
))-sqrt(2)*z*(4*z^2+5*a^2)*ln((sqrt(2*z^2+a^2*(t^2+1))-sqrt(2*z^~
2+a^2))/t)/(8*a^2*(2*z^2+a^2)^(3/2))+(2*z^2+a^2)*ln(sqrt(2*z^2+a~
^2*(t^2+1))-sqrt(2)*z)/(8*a^2*z^2)-(5*z^2+a^2)*(8*z^6+9*a^2*z^4+~
6*a^4*z^2+a^6)*ln(2*z^2+a^2*(t^2+1))/(16*a^2*z^2*(3*z^2+a^2)^3)+~
(a^2-z^2)*(sqrt(2)*sqrt(2*z^2+a^2*(t^2+1))*(10*z^6+2*a*t*z^5+16*~
a^2*z^4+7*a^4*z^2+a^6)+12*z^7+6*a*t*z^6+16*a^2*z^5+5*a^3*t*z^4+7~
*a^4*z^3+a^5*t*z^2+a^6*z)/(8*a^2*z*(2*z^2+a^2)*(2*z^2+a^2*(t^2+1~
))*(3*z^2+a^2)^2))
This is also continuous along the integration path. For twice the
integral from p=0 (t=1) to p=pi (t=-1) of this term plus its complex
conjugate fro the second group one then obtains (for z>0) the fairly
simple result:
pi*z/a^2*((4*z^2+5*a^2)/(sqrt(2)*(2*z^2+a^2)^(3/2))-(z^2+2*a^2)/~
(z^2+a^2)^(3/2))
The remaining terms can be done in the same way, but I refrain from
listing the results here.
Remember: don't pay USD 5000 next time unless a CAS vendor can
demonstrate the usefulness of his integrator with this nice elementary
double integral.
Martin.