Convert solution of EllipticE and Inverse Function to C?

Questions
, , ,
1

Hello;

I have this below result as an analytic solution to my differential equation. Can you please help me to understand :

  1. How to Ccode this result?
  2. What is ellipticE[#1/2,…]? What is this dash sign?
  3. what is & sign used for?
  4. what is C1?

I am trying to convert this result to C language to write it in microprocessor. Can someone guide me on this?

DiffDenkSolutionWolfram.JPG
DiffDenkSolutionWolfram.JPG705×137 16.4 KB

2

Maybe this helps a little:

Internal`InheritedBlock[{Solve}, Unprotect[Solve];
Solve[x___] := 
Block[{$guard = True}, Print["Solve called : ", HoldForm[Solve[x]]];
Solve[x]] /; ! TrueQ[$guard];
DSolve[{-x^2 + (0.232866 - 0.232166*Cos[y[x]])*y'[x]^2 == 0}, y[x], x]];

(*Solve called : Solve[10 Sqrt[14] EllipticE[y[x]/2,-(116083/175)]==-((5714200963 x^2)/16162201)+Subscript[\[ConstantC], 1],y[x]]*)

We see that solution is implicit form:

10 Sqrt[14] EllipticE[y[x]/2,-(116083/175)]==-((5714200963 x^2)/16162201)+C[1]

where: C1 is a integration constant.

For: EllipticE,#and & see in Documentation Center

3

It would be easier and probably more reliable to get the numeric solution on some interval and try to convery the interpolating function to C code.

4

Hello Daniel;

Like you said, I tried to have the numeric solution as an interpolating function. Inside mathematica, its possible to use this function , but I couldnt find a way to export this function out or look at the parameters of the function. How its done? Please see attached.

Also the numbers in the equation (0.232866 and -0.232166) are going to be parametric, say A and B. Is it possible to receive the interpolating function parameters interms of A and B outside of Mathematica in this case?

InterpolatingFunction.JPG
InterpolatingFunction.JPG837×186 18.9 KB

5

Hello Mariusz;

Thanks for this great help. Now from the documentation I see that

Elliptic1.JPG

where

Elliptic2.JPG

so in your case

y(x)/2 is Phi (and it means it will go inside Sine term. )

-(116083/175) is m

Is this correct?

6

Yes you are correct. We can equation write:

HoldForm[Integrate[(1 + (116083/175)*Sin[t]^2)^(1/2), {t, 0, 
    y[x]/2}] == (-((5714200963 x^2)/16162201) + C[1])/(10 Sqrt[14])] //  TeXForm
\int_0^{\frac{y(x)}{2}} \sqrt{1+\frac{116083 \sin ^2(t)}{175}} \, dt=\frac{-\frac{5714200963 x^2}{16162201}+c_1}{10 \sqrt{14}}
7

Caner,

Building on Daniels’s suggestion, If your goal is to code this all in C, I would suggest considering using a polynomial approximation of the Elliptic integral. For example you can find it in Abramowitz and stegun p 592 You can use Mathematica to reformulate the result and insert the approximation and then simplify/solve it into a compact form for coding in C. Also, you can check the result against examples in MMA.

I hope this helps,

Regards,

Neil

8

Hello Mariusz;

I am trying to apply Neil’s approximation advice to your solution.
In your code, when I enter initial conditions, nothing changes in the answer. Integration constant remains. How do I add initial condition to your code to have an effect on the answer?

Internal`InheritedBlock[{Solve}, Unprotect[Solve];
  Solve[x___] := 
   Block[{$guard = True}, Print["Solve called : ", HoldForm[Solve[x]]];
     Solve[x]] /; ! TrueQ[$guard];
  DSolve[{-x^2 + (0.232866 - 0.232166*Cos[y[x]])*y'[x]^2 == 0, 
    y[2] == 1.57}, y[x], x]];
9

Try this:

DSolve[{-x^2 + (0.232866 - 0.232166*Cos[y[x]])*y'[x]^2 == 0, y[2] == 1.57}, y[x], x]
(*{{y[x] -> 
   InverseFunction[10 Sqrt[14] EllipticE[#1/2, -(116083/175)] &][
    1699.28 - 353.553 x^2]}, {y[x] -> 
   InverseFunction[
     10 Sqrt[14] EllipticE[#1/2, -(116083/175)] &][-1129.14 + 
     353.553 x^2]}}*)

Then we extracting from outputs two solution:

  {10 Sqrt[14] EllipticE[y[x]/2, -(116083/175)] == 
    1699.2844439496437` - 353.5533905932738` x^2, 
   10 Sqrt[14]
      EllipticE[y[x]/2, -(116083/175)] == -1129.1426807965465` + 
     353.5533905932738` x^2}
10

Hello Mariusz;

Many thanks for your kind support and giving me a hand. I clearly understand that I need to have an analytical solution of Inverse ellipctical integral of second kind. I researched the internet and did not find any approximation for the inverse.

So , in the end, I need to solve the equation EllipticE[y[x]/2 , -k]=a+b*x^2.

How can I analytically solve this equation, by letting y(x) alone? Can mathematica help me on this?

Caner

11

Caner,

My thought was to use the polynomial approximation for the elliptic integral and solve the polynomial in MMA. If the polynomial is too complicated for an analytical solution, I would resort to a simple numerical Solve. You can implement that easily in C.

Regards,

Neil

12

No,It’s mathematicaly impossible,because you have a transcendental equation to solve.
With Mathematica is possible to plot solution from Dsolve:

 p1 = ContourPlot[{10 Sqrt[14] EllipticE[y/2, -(116083/175)] == 
     1699.2844439496437` - 353.5533905932738` x^2, 10 Sqrt[14]
      EllipticE[y/2, -(116083/175)] == -1129.1426807965465` + 
      353.5533905932738` x^2}, {x, -3, 3}, {y, -4, 4}, 
   FrameLabel -> Automatic]

4834Bez tytułu.png
4834Bez tytułu.png432×422 6.53 KB

13

Hello Neil;

You are absolutely right, numerical solution would be a good way to go. But…
The problem with numerical solution is that, i need to equate 2 nonlinear differential equations at a Point X to find X. Inorder to do this I have to have 2 analytical solutions at hand to equate. With numerical solution, I dont have any equation at hand.
Unfortunately the analytic solutions of these differential equations involve inverse Elliptic Integrals. A polynomial approximation for the inverse elliptical integral (should be inverse) would work perfect for me. In your manual, it only gives direct approximation for the elliptic integral not the inverse approximation.