Differences

This shows you the differences between two versions of the page.

--- star_destroyers [2014/06/06 09:34]
cook
+++ star_destroyers [2014/06/06 14:04] (current)
cook
@@ Line 9: / Line 9: @@
-Project Considered:  Post-Processing nucleosynthesis
+**Final Presentation**
+The final presentation may be found here, without any videos.
+{{::presentation.pdf|}}
+**Project Considered:  Post-Processing nucleosynthesis**
 Extract a Lagrangian profile from the blast wave simulation and use the network to calculate the nucleosynthesis.
-The plan: First use the Lane-Emden equation to model the density distribution. Use this as input for XNet and evolve at the right temperature for the initial NSE stage.
+**The plan**: First use the Lane-Emden equation to model the density distribution. Use this as input for XNet and evolve at the right temperature for the initial NSE stage.
 Modify VH1 for Lagrangian dynamics, then apply the Sedov blast on the initial density distribution (from the same solution to the Lane-Emden equation). This solution will then be used as input to XNet, using a 150 element network.
-Interesting outputs: Abundance of elements as a function of time, and after the blast, as a function of radius. It is intended that this will be displayed in video form. Another output of interest is the total energy released/required in nucleosynthesis - is this energy significant enough to effect the hydrodynamics? To what extent?
+**Interesting outputs**: Abundance of elements as a function of time, and after the blast, as a function of radius. It is intended that this will be displayed in video form. Another output of interest is the total energy released/required in nucleosynthesis - is this energy significant enough to effect the hydrodynamics? To what extent?
-------Logbook starts here -------
+**Parameters for the Lane-Emden Equation**:
+  * ρ_c = 1.6*10^5
+  * μ = 1/2
+  * n = 3/2
+  * γ = n +1/n
-/6/14 3pm KJC: Lane-Emden Equation: We solved this in Mathematica, because, well, it only takes two lines. Far more lines to actually give a good output for FORTRAN, but that's ok. The notebook may be found on the git [[https://github.com/KJCook/sedov-solution/blob/master/LaneEmden.nb|here]]. Output is three columns - radius (cm), density (g/cm^3) and pressure (dyn/cm^2). A sample output is in the github
+{{::laneemden.png?300|}}
-Parameters for the Lane-Emden Equation:
+**VH1 Parameters**:
-ρ_c = 1.6*10^5
+  * E_blast = 10^49 erg (all placed within innermost cell). Leads to a max temperature of 3 GK in zone 30/400. A full 10^51 erg explosion causes temperatures to exceed XNet's capabilities. This is because the 5/3 cold degenerate equation of state with mu=1/2 is crap in this regime. A star with radially varying composition, and an equation of state to reflect that, is required for additional realism.
-μ = 1/2
+**Making VH1 Starter Lagrangian**
-n = 3/2
+ppmlr.f90 needs to be modified so that once the Lagrange step was taken, the values needed to not be put back onto the Eulerian grid.
+  * remove the "wiggle" part of the file that adds dissipation for carbuncular noise.  Since this is only really a problem in higher dimension, we completely removed this for simplicity.
+  * Instead of defining xf=xf0 define xf0=xf (i.e. make the Eulerian coordinates move to match the Lagrangian ones)
+  * If the final "remap" command is commented out, the code misbehaves, though we don't know why a remap should be needed. Leave the final remap in.
-γ = n +1/n
+**Units**
+VH1 outputs NaN's when the simulation is run in CGS units. Presumably, this is because it is programmed using floats rather than double precision numbers, and the separation in numbers results in problematic subtractive cancellation. To fix this problem, we re-defined our units of mass and radius to be those of solar mass and solar radius. The following led to successful simulations.
-/6/14 4pm KJC: JGL Got VH1 up and running for the Sedov blast as written for the problems done last week. We now have to combine the solution from the Mathematica output into the VH1 scripts. SR is working on this.
+  * M_cgs = M_sun * M_code
+  * R_cgs = R_sun * R_code
+  * rho_cgs = M_sun/R_sun^3 * rho_code
+  * P_cgs = M_sun/R_sun * P_code
+  * T_cgs = M_sun*R_sun^2 * T_code
+  * E_cgs = M_sun*R_sun^2 * E_code
-/6/14 4pm JGL: Have netCDF installed correctly, and compiled serial version of VH1 using a modified makefile (github upload to follow).  The vh1-serial code outputs an *.nc file, also installed ncview to view this file.  With this, we can view the 2D Sod tube solution via the convenient ncview program.  For the starter version, I added a python script written by Nathan Parzuchowski to animate the 1D Sedov blast solution by turning each time step data file into a frame for the animation. We will probably stick with 1D solutions, since we want to do our final simulation with the Lagrangian remap which is supposedly more well behaved in 1D. The python script requires python-numpy and python-matplotlib.
+**XNet Parameters and Inputs**
-/6/14 6pm JGL: Just added the reaction network code to the git repository.  Should be up and running for the use of 150 isotopes now.  We plan to eventually run the reaction network code for both just the alpha network as well as the 150 isotope network to compare the results.  We may also try other networks if we have the time.
+Three separate initial conditions:
+  * pure He
+  * 50% N, 50% P
+  * pure 56Fe
-/6/14 6pm JGL: Today JGL and SR worked on getting the Lagrange remap up and running.  The ppmlr.f90 needed to be modified so that once the Lagrange step was taken, the values needed to not be put back onto the Eulerian grid.  First we removed the "wiggle" part of the file that adds dissipation for carbuncular noise.  Since this is only really a problem in higher dimension, we removed this for simplicity.  Next, we updated the grid using the Lagrange step, however, the code began taking extremely small time steps (on the order of 10^-7) so the simulation got nowhere.  However, when we put back in the final "remap" call, the code seemed to start working as intended.  We are not entirely sure what this call does that fixes everything, but now we have ran the Sedov blast using both Eulerian and Lagrangian coordinates, and both achieve the same physics, but the Lagrangian coordinates resolves the sharp blast wave better.  The exact changes to the code can be found on our github.
+To modify XNet to output 150 nuclei in the plaintext file. (Instead of trying to read in the fortran binary files in python, which is in principle also possible.) You need to change:
-The next problem we are all now working on, is how to put the Lane-Emden (LE) solutions into our code as input.  Using cgs units, the LE input gives densities on the order of (10^5) and pressures of the order of (10^22), which proptly breaks the code after the first times step, as all of the pressures drop to ~10^-15, and all the physics goes away.  We are currently exploring ways to non-dimensionalize our code so that all of the inputs are closer to order unity as this should help.  We are honestly not sure why the code fails in the first place, since even single precision should be able to handle pressures in the 10^22 regime.
+  * The ''control'' file, you must add in the list of all the nuclei you wish to output. Recall each should use 5 characters.
+  * the ''inout'' variable in ''control.f90''.
+  * ''output_nuc'' and
+  * ''Read(lun_control,"(14a5)") output_nuc'' to ''150a5'',
+  * The index in the loop on line 198,
+  * The write statement for ''lun_ev''on lines 352, 465, 473 in the file ''net.f90''. (Do a search for 14 in the files).
-KJC: Implemented a python script for plotting nucleosynthesis on a chart of the isotopes from XNet output. Animation requires animation.FuncAnimation(), a function that proved rather non-trivial to implement. If you want to do this, make sure your first timestep doesn't have identical values. No idea why this matters, but hey. Code is on the github and is called ''plotchart.py''. Here's a screen shot of the animation. We need to write out more isotopes from XNet before this will look really interesting.
+As usual, all the code is in the github.
- {{::screen_shot_2014-06-03_at_6.20.49_pm.png?200|}}
-KJC: Modified XNet to output 150 nuclei in the plaintext file. (Instead of trying to read in the fortran binary files in python, which is in principle also possible.) You need to change:
+------Logbook starts here -------
-The ''control'' file, you must add in the list of all the nuclei you wish to output. Recall each should use 5 characters.
+/6/14 3pm KJC: Lane-Emden Equation: We solved this in Mathematica, because, well, it only takes two lines. Far more lines to actually give a good output for FORTRAN, but that's ok. The notebook may be found on the git [[https://github.com/KJCook/sedov-solution/blob/master/LaneEmden.nb|here]]. Output is three columns - radius (cm), density (g/cm^3) and pressure (dyn/cm^2). A sample output is in the github
-Also, change the ''inout'' variable in ''control.f90''.
+/6/14 4pm KJC: JGL Got VH1 up and running for the Sedov blast as written for the problems done last week. We now have to combine the solution from the Mathematica output into the VH1 scripts. SR is working on this.
-You also need to change:
+/6/14 4pm JGL: Have netCDF installed correctly, and compiled serial version of VH1 using a modified makefile (github upload to follow).  The vh1-serial code outputs an *.nc file, also installed ncview to view this file.  With this, we can view the 2D Sod tube solution via the convenient ncview program.  For the starter version, I added a python script written by Nathan Parzuchowski to animate the 1D Sedov blast solution by turning each time step data file into a frame for the animation. We will probably stick with 1D solutions, since we want to do our final simulation with the Lagrangian remap which is supposedly more well behaved in 1D. The python script requires python-numpy and python-matplotlib.
-''output_nuc'' and
+/6/14 6pm JGL: Just added the reaction network code to the git repository.  Should be up and running for the use of 150 isotopes now.  We plan to eventually run the reaction network code for both just the alpha network as well as the 150 isotope network to compare the results.  We may also try other networks if we have the time.
-''Read(lun_control,"(14a5)") output_nuc'' to ''150a5'',
+/6/14 6pm JGL: Today JGL and SR worked on getting the Lagrange remap up and running.  The ppmlr.f90 needed to be modified so that once the Lagrange step was taken, the values needed to not be put back onto the Eulerian grid.  First we removed the "wiggle" part of the file that adds dissipation for carbuncular noise.  Since this is only really a problem in higher dimension, we removed this for simplicity.  Next, we updated the grid using the Lagrange step, however, the code began taking extremely small time steps (on the order of 10^-7) so the simulation got nowhere.  However, when we put back in the final "remap" call, the code seemed to start working as intended.  We are not entirely sure what this call does that fixes everything, but now we have ran the Sedov blast using both Eulerian and Lagrangian coordinates, and both achieve the same physics, but the Lagrangian coordinates resolves the sharp blast wave better.  The exact changes to the code can be found on our github.
-The index in the loop on line 198,
+/6/2014: KJC: Implemented a python script for plotting nucleosynthesis on a chart of the isotopes from XNet output. Animation requires animation.FuncAnimation(), a function that proved rather non-trivial to implement. If you want to do this, make sure your first timestep doesn't have identical values. No idea why this matters, but hey. Code is on the github and is called ''plotchart.py''. Here's a screen shot of the animation. We need to write out more isotopes from XNet before this will look really interesting.
-The write statement for ''lun_ev''on lines 352, 465, 473.
+ {{::screen_shot_2014-06-03_at_6.20.49_pm.png?200|}}
-in the file ''net.f90''. (Do a search for 14 in the files).
-As usual, all the code is in the github.
 Modified ''plotchart.py'' to allow as many isotopes to be plotted as are in the network you're interested in. The script now outputs an ''.mp4'' file for display.
@@ Line 79: / Line 102: @@
 {{:npfinal.png?600|}}
+Meanwhile, JGL wrote a python script for taking outputs from ''VH1'' and passing them into the ''XNet'' input file. SR continues to work on the units. It seems that we have used the wrong expression for pressure.
+**Conclusions**
+We successfully simulated the nucleosynthesis occurring in a 1 solar mass polytropic star being given a sedov shock. We found that final abundances depend strongly on blast strength and initial abundances. This is hardly surprising, but it is nice to see.
+**Further work**
+Or, What we would have done given some more time:
+Added in Gravity. Our blasts are sufficiently energetic that this shouldn't matter too much, but it would be good to have. Add a more realistic EOS. More realistic initial conditions for hydro. Further variations of initial abundances. Follow the nucleosynthesis with different blast energies. Do nucleosynthesis on all zones. Feed the nucleosynthesis energy into the hydro.
 [[https://github.com/KJCook/sedov-solution.git|GitHub]]

Talent

User Tools

Site Tools

Differences

Page Tools