Drexel dragonThe Math ForumDonate to the Math Forum



Search All of the Math Forum:

Views expressed in these public forums are not endorsed by Drexel University or The Math Forum.


Math Forum » Discussions » Software » comp.soft-sys.matlab

Topic: ellipse points regularly spaced, or inverse of the incomplete
elliptic integral of the second (not first) kind

Replies: 8   Last Post: Jun 1, 2010 3:45 PM

Advanced Search

Back to Topic List Back to Topic List Jump to Tree View Jump to Tree View   Messages: [ Previous | Next ]
Bruno Luong

Posts: 8,813
Registered: 7/26/08
Re: ellipse points regularly spaced, or inverse of the incomplete elliptic integral of the second (not first) kind
Posted: Jun 1, 2010 3:45 PM
  Click to see the message monospaced in plain text Plain Text   Click to reply to this topic Reply

"David W. Cantrell" <DWCantrell@sigmaxi.net> wrote in message <htmffh$jnc$1@news.acm.uiuc.edu>...
>
> [By the way, the second link mentioned there, to the sci.math.research
> thread "Inverting elliptic integrals", seems to be broken. However, that
> thread can now be found at the Math Forum:
> <http://mathforum.org/kb/message.jspa?messageID=1672061>.]
>
> Regards,
> David W. Cantrell


Incidentally, I bump to the same problem today, and I try David's method. Unless if I make an error somewhere, the power series seem to have poor accuracy for m ~ 1 and z >> 0.

I hope I make a mistake somewhere:

function E = elliptic12_powerseries(z, m)
% Compute EllipticE by power series in z
% Provided by David Cantrell
% http://mathforum.org/kb/message.jspa?messageID=1672061

% E = z + ...
% (m*z^3)/6 + ...
% (1/120)*(-4*m + 13*m^2)*z^5 + ...
% ((16*m - 284*m^2 + 493*m^3)*z^7)/5040 + ...
% ((-64*m + 4944*m^2 - 31224*m^3 + 37369*m^4)*z^9)/362880 + ...
% ((256*m - 81088*m^2 + 1406832*m^3 - 5165224*m^4 +
% 4732249*m^5)*z^11)/39916800;

P = {1;
(1/6)*[1 0];
(1/120)*[13 -4 0];
(1/5040)*[493 -284 16 0];
(1/362880)*[37369 -31224 4944 -64 0];
(1/39916800)*[4732249 -5165224 1406832 -81088 256 0] };
P = cellfun(@(p) polyval(p,m), P);
E = z.*polyval(P(end:-1:1),z.^2);

end % elliptic12_powerseries

%%
m = 0.9;
theta=linspace(0,pi/2);
Epowerseries = elliptic12_powerseries(theta, m);
% http://www.mathworks.com/matlabcentral/fileexchange/8805
[trash Eanalytic] = elliptic12(theta, m);

plot(theta, Eanalytic, 'r', theta, Epowerseries, 'b')
title(sprintf('Ellipictic integral 2nd kind, m = %g', m));
xlabel('z (rad)')
legend('Analytic', 'power series','Location','NorthWest');

Any comment?

Bruno



Point your RSS reader here for a feed of the latest messages in this topic.

[Privacy Policy] [Terms of Use]

© Drexel University 1994-2014. All Rights Reserved.
The Math Forum is a research and educational enterprise of the Drexel University School of Education.