Commits

Sam Skillman committed 48ca7c2

Fixing PressurelessCollapse test problem. A combination of a small Courant
safety factor and relatively large tiny_number led to undesired results that
stop the simulation too soon. This change increases the safety factor, and
returns tiny_number to the default of 1.0e-20. Also cleaned up the parameter
file to more accurately reflect the intended behavior.

Added a make_plots.py, which can be used to make a quick look at the resulting
simulation.

Comments (0)

Files changed (2)

run/Hydro/Hydro-1D/PressurelessCollapse/PressurelessCollapse.enzo

 TopGridDimensions      = 100
 SelfGravity            = 1       // gravity on
 TopGridGravityBoundary = 1       // Isolated BCs
-LeftFaceBoundaryCondition  = 1    // outflow ?
-RightFaceBoundaryCondition = 1    // outflow ?
+LeftFaceBoundaryCondition  = 1   // outflow
+RightFaceBoundaryCondition = 1   // outflow
 PressureFree           = 1       // turn off pressure
 #
 #  set I/O and stop/start parameters
 #  set hydro parameters
 #
 Gamma                  = 1.4
-CourantSafetyNumber    = 0.05    // needs to be lower for pressurefree
-PPMDiffusionParameter  = 0       // diffusion off
+CourantSafetyNumber    = 0.4    // needs to be lower for pressurefree
+PPMDiffusionParameter  = 0      // diffusion off
 #
 #  set grid refinement parameters
 #
-StaticHierarchy           = 1    // dynamic hierarchy
-MaximumRefinementLevel    = 1    // use up to 2 levels
-RefineBy                  = 4    // refinement factor
-MinimumSlopeForRefinement = 0.2  // set this to <= 0.2 to refine CD
+StaticHierarchy           = 1    // static hierarchy
 #
 #  set some global parameters
 #
-SubcycleSafetyFactor   = 2       // 
-tiny_number            = 1.0e-10 // fixes velocity slope problem
-MinimumEfficiency      = 0.4     // better value for 1d than 0.2
-Initialdt = 1e-6
+Initialdt              = 1.0e-6
+

run/Hydro/Hydro-1D/PressurelessCollapse/make_plots.py

+from yt.mods import *
+import os
+import sys
+import pylab
+
+def make_plot(pfname):
+    pf = load(pfname)
+    ### extract an ortho_ray (1D solution vector)
+    ray = pf.h.ortho_ray(0, [0.5, 0.5])
+
+    ### define fields vector
+    fields = ('Density', 'x-velocity', 'TotalEnergy', 'Pressure' )
+
+    ### make plot
+
+    pylab.figure(1, figsize=(8,7))
+
+    # Density Plot
+    a = pylab.axes([0.09, 0.57, 0.38, 0.38])
+    pylab.axhline(0,color='k',linestyle='dotted')
+    pylab.plot(ray['x'],ray['Density'], 'ro', ms=4)
+
+    pylab.xlabel('Position')
+    pylab.ylabel('Density')
+
+    # Velocity Plot
+    a = pylab.axes([0.59, 0.57, 0.38, 0.38])
+    pylab.axhline(0,color='k',linestyle='dotted')
+    pylab.plot(ray['x'],ray['x-velocity'], 'ro', ms=4)
+
+    pylab.xlabel('Position')
+    pylab.ylabel('Velocity')
+
+    # TotalEnergy Plot
+    a = pylab.axes([0.59, 0.07, 0.38, 0.38])
+    pylab.axhline(0,color='k',linestyle='dotted')
+    pylab.plot(ray['x'],ray['TotalEnergy'], 'ro', ms=4)
+
+    pylab.xlabel('Position')
+    pylab.ylabel('Total Energy')
+
+    ### Save plot
+    pylab.savefig('%s.png' % pf)
+    pylab.clf()
+    
+if __name__ == '__main__':
+    for i in range(11):
+        try: 
+            make_plot('DD%04i/data%04i'% (i,i))
+        except:
+            break
+
+    # To make a movie using avconv, uncomment the following 2 lines
+    # os.system('avconv -r 10 -i data%04d.png -threads 8 -pass 1 -an -f webm -b 2000k movie.webm')
+    # os.system('avconv -r 10 -i data%04d.png -threads 8 -pass 2 -an -f webm -b 2000k movie.webm')
+