Numerical Integration of Chandrasekhar's Equation

vmr101
Gold Member
Messages
25
Reaction score
1

Homework Statement


We need to write an integrator for the Chandrasekhars Equation (CE) for White Dwarfs (WD) using python3/NumPy/Matplotlib. We then need to compute the structure of a WD made of our varying elements. We also need to compute and plot the mass-radius relation for WD.

Homework Equations


We are given equations for C1, C2, K, eta, alpha, Radius, Mass and density.
rho = C2*x3=C2*zc3(phi2-1/(zc2))3/2

The Attempt at a Solution


I have never used python3 before so this is a bit of a challenge.
My script sets up all the constants and variables, sets up the arrays.
#then uses the initial conditions and sets an array for the ce
#then sets the integrator as rk4
#then integrates the CE
#calculates rho, alpha, p, temp, r
M = ((4*pi)/c2**2)*(((2*c1)/(pi*G))**(3/2))*((-x**2)*y[:,0])
#Plot r vs rho, r vs temp etc..

I need to figure out where to vary zc to find a certain WD mass
Any suggestions would be much appreciated.
 
Last edited:
Physics news on Phys.org
vmr101 said:

Homework Statement


We need to write an integrator for the Chandrasekhars Equation (CE) for White Dwarfs (WD) using python3/NumPy/Matplotlib. We then need to compute the structure of a WD made of our varying elements. We also need to compute and plot the mass-radius relation for WD.

Homework Equations


We are given equations for C1, C2, K, eta, alpha, Radius, Mass and density.
rho = C2*x3=C2*zc3(phi2-1/(zc2))3/2

The Attempt at a Solution


I have never used python3 before so this is a bit of a challenge.
My script sets up all the constants and variables, sets up the arrays.
ce1 = functools.partial(le.ce,zc=zc) #this uses the initial conditions and sets an array for the ce
a = functools.partial(le.rk4,h=h) #this sets the integrator as rk4
x, y = le.integrate(ce1,x,y,a,t,zc) #this integrates the CE

rho = c2*zc**3*(y[:,0]**2 -1/zc**2)**(3/2)
alpha = ((2*c1)/(pi*G))**(1/2) * 1/(c2*zc)
pr = K*rho**gamma
temp = (pr*mewe)(rho/Rgas) ##error here
r = alpha*x
M = ((4*pi)/c2**2)*(((2*c1)/(pi*G))**(3/2))*((-x**2)*y[:,0])
Mmax = 4*pi*(2.01824)*(K/pi*G)**(3/2)
#Plot r vs rho, r vs temp etc..

Any suggestions would be much appreciated.
Your line where you have "this is where I get my error" needs a phi^2 in place of the ## (1/zc)^2 ##.
 
Charles Link said:
Your line where you have "this is where I get my error" needs a phi^2 in place of the ## (1/zc)^2 ##.
I just noticed and updated this too! Thanks for confirming it for me.

Looks like it produces a density plot, but I need to vary zc to suit the WD mass.
Also my mass isn't plotting anything.
 

Attachments

Last edited:
This was mostly python issues, and not issues with the physics itself.
I managed to work through it.
 
Hi, I had an exam and I completely messed up a problem. Especially one part which was necessary for the rest of the problem. Basically, I have a wormhole metric: $$(ds)^2 = -(dt)^2 + (dr)^2 + (r^2 + b^2)( (d\theta)^2 + sin^2 \theta (d\phi)^2 )$$ Where ##b=1## with an orbit only in the equatorial plane. We also know from the question that the orbit must satisfy this relationship: $$\varepsilon = \frac{1}{2} (\frac{dr}{d\tau})^2 + V_{eff}(r)$$ Ultimately, I was tasked to find the initial...
The value of H equals ## 10^{3}## in natural units, According to : https://en.wikipedia.org/wiki/Natural_units, ## t \sim 10^{-21} sec = 10^{21} Hz ##, and since ## \text{GeV} \sim 10^{24} \text{Hz } ##, ## GeV \sim 10^{24} \times 10^{-21} = 10^3 ## in natural units. So is this conversion correct? Also in the above formula, can I convert H to that natural units , since it’s a constant, while keeping k in Hz ?
Back
Top