Home

Staff
     
Consulting
               
Mathematica Tools

Customer Quotes
        
Contact

Older Mathematica Items

 

 

 

 

 

Scientific Arts

Stellar Structure and the Lane-Emden Function: Numerical Solution

Page 1 | Page 2


Numerical Solution

NDSolve

Here we want to numerically integrate Equation 1 and use the result to find the value of [Graphics:../Images/laneemden_gr_43.gif]. Chandrasekhar presents the following table, where the numerical results were obtained from a variety of sources.

[Graphics:../Images/laneemden_gr_44.gif]
Table 1

First try to use NDSolve for [Graphics:../Images/laneemden_gr_45.gif].

[Graphics:../Images/laneemden_gr_46.gif]
[Graphics:../Images/laneemden_gr_47.gif]
[Graphics:../Images/laneemden_gr_48.gif]
[Graphics:../Images/laneemden_gr_49.gif]
[Graphics:../Images/laneemden_gr_50.gif]

Here the numerical integration was not able to leave the starting gate. The same problem happens when [Graphics:../Images/laneemden_gr_51.gif] (as well as in other cases).

[Graphics:../Images/laneemden_gr_52.gif]
[Graphics:../Images/laneemden_gr_53.gif]
[Graphics:../Images/laneemden_gr_54.gif]
[Graphics:../Images/laneemden_gr_55.gif]
[Graphics:../Images/laneemden_gr_56.gif]

But for [Graphics:../Images/laneemden_gr_57.gif] there is the simple analytic solution presented above.

The reason for the difficulty is seen in the expanded form of Equation 1.

[Graphics:../Images/laneemden_gr_58.gif]

From this we can see that, due to the potentially singular term, there might be difficulty with a numerical algorithm that starts at [Graphics:../Images/laneemden_gr_59.gif]. However, if we can move the boundary condition away from the origin, we can use NDSolve from that point on. We can use the differential equation, Equation 1, along with the boundary conditions at the origin, Equation 2, to make a series expansion for [Graphics:../Images/laneemden_gr_60.gif] about [Graphics:../Images/laneemden_gr_61.gif].  

Series Expansion

We will create a general function for doing this, but will develop it by working with a particular case: a test case for [Graphics:../Images/laneemden_gr_62.gif].

[Graphics:../Images/laneemden_gr_63.gif]

First define the series in terms of unknown parameters.

[Graphics:../Images/laneemden_gr_64.gif]
[Graphics:../Images/laneemden_gr_65.gif]

In this expression the boundary conditions were imposed by choosing [Graphics:../Images/laneemden_gr_66.gif] and [Graphics:../Images/laneemden_gr_67.gif]. Now define a table of the unknown coefficients

[Graphics:../Images/laneemden_gr_68.gif]

The differential equation defines an expression for the unknown coefficients.

[Graphics:../Images/laneemden_gr_69.gif]

Collect the coefficients in this expression.

[Graphics:../Images/laneemden_gr_70.gif]

Now solve the equations that result from setting the coefficients to zero (the [Graphics:../Images/laneemden_gr_71.gif] in [Graphics:../Images/laneemden_gr_72.gif] creates the list of equations from the list of coefficients).

[Graphics:../Images/laneemden_gr_73.gif]

Here is the resulting series.

[Graphics:../Images/laneemden_gr_74.gif]
[Graphics:../Images/laneemden_gr_75.gif]

Now we put all of this together into a function. We make the function self-documenting in the conventional Mathematica fashion through the device of a usage message.

[Graphics:../Images/laneemden_gr_76.gif]

The usage message allows us to access the definition of the function.

[Graphics:../Images/laneemden_gr_77.gif]
[Graphics:../Images/laneemden_gr_78.gif]

It also enables the template creation feature the for the function (accessible through [Graphics:../Images/laneemden_gr_79.gif] menu item, or its keyboard equivalent).

Here is an example of the use of LaneEmdenSeries.

[Graphics:../Images/laneemden_gr_80.gif]
[Graphics:../Images/laneemden_gr_81.gif]

(Chandrasekhar shows the result to [Graphics:../Images/laneemden_gr_82.gif]; but, in his case we are sure the he could calculate it out to a significantly higher order: presumably he did. For mortals like us, we can do the same - without error - with Mathematica as a tool.)

Using the Series Solution in NDSolve

Since we designed the function LaneEmdenSeries to automatically incorporate the boundary conditions of Equation 2, it can be used to move the boundary condition to any other point away from the origin where the series is convergent. A prudent approach is to choose a point close to 0, say
[Graphics:../Images/laneemden_gr_83.gif], and explore the result.

First compute the new boundary conditions at [Graphics:../Images/laneemden_gr_84.gif].

[Graphics:../Images/laneemden_gr_85.gif]
[Graphics:../Images/laneemden_gr_86.gif]
[Graphics:../Images/laneemden_gr_87.gif]
[Graphics:../Images/laneemden_gr_88.gif]

Now attempt to use NDSolve with the new boundary conditions. We can experiment with machine arithmetic (the default) versus higher precision arithmetic by using various options to NDSolve. We find that the quoted results (of Table 1) are achieved in the latter case (for the series expansion point and degree that we have picked) if we set the option [Graphics:../Images/laneemden_gr_89.gif].

Here is the result for [Graphics:../Images/laneemden_gr_90.gif] (we use Module to localize the parameters).

[Graphics:../Images/laneemden_gr_91.gif]
[Graphics:../Images/laneemden_gr_92.gif]

Let's take a look at the result.

[Graphics:../Images/laneemden_gr_93.gif]

[Graphics:../Images/laneemden_gr_94.gif]

Find the first real zero using FindRoot (with the option [Graphics:../Images/laneemden_gr_95.gif]).

[Graphics:../Images/laneemden_gr_96.gif]
[Graphics:../Images/laneemden_gr_97.gif]
[Graphics:../Images/laneemden_gr_98.gif]
[Graphics:../Images/laneemden_gr_99.gif]

These numbers agree with the results quoted in Table 1. The open question here is "how many of the digits of this answer are correct?" This can be investigated by changing the various options to NDSolve and FindRoot as well as the degree and starting point of the Lane-Emden series used to move the boundary conditions. The advantage of having access to arbitrary precision arithmetic is that you can investigate the numerical properties of a solution with much greater versatility than when you are under the constraints of machine arithmetic.

The procedure that we just used for the [Graphics:../Images/laneemden_gr_100.gif] case can be put together into a function so that you can experiment with and generate results for any case of interest. In this way you can reconstruct all of the entries of Table 1.


For further information on our services send email to info@scientificarts.com .
Contents of this web site Copyright © 1999-2011, Scientific Arts, LLC.

s