diff --git a/examples/SupernovaeFeedback/SN_feedback.yml b/examples/SupernovaeFeedback/SN_feedback.yml new file mode 100644 index 0000000000000000000000000000000000000000..a59ae302ff4052a3daf0535e93a0c2cd5e9904f5 --- /dev/null +++ b/examples/SupernovaeFeedback/SN_feedback.yml @@ -0,0 +1,44 @@ +# Define the system of units to use internally. +InternalUnitSystem: + UnitMass_in_cgs: 1 # Grams + UnitLength_in_cgs: 1 # Centimeters + UnitVelocity_in_cgs: 1 # Centimeters per second + UnitCurrent_in_cgs: 1 # Amperes + UnitTemp_in_cgs: 1 # Kelvin + +# Values of some physical constants +PhysicalConstants: + G: 0 # (Optional) Overwrite the value of Newton's constant used internally by the code. + +# Parameters governing the time integration +TimeIntegration: + time_begin: 0. # The starting time of the simulation (in internal units). + time_end: 5e-2 # The end time of the simulation (in internal units). + dt_min: 1e-7 # The minimal time-step size of the simulation (in internal units). + dt_max: 1e-4 # The maximal time-step size of the simulation (in internal units). + +# Parameters governing the snapshots +Snapshots: + basename: SN_feedback # Common part of the name of output files + time_first: 0. # Time of the first output (in internal units) + delta_time: 1e-2 # Time difference between consecutive outputs (in internal units) + compression: 1 + +# Parameters governing the conserved quantities statistics +Statistics: + delta_time: 1e-3 # Time between statistics output + +# Parameters for the hydrodynamics scheme +SPH: + resolution_eta: 1.2348 # Target smoothing length in units of the mean inter-particle separation (1.2348 == 48Ngbs with the cubic spline kernel). + CFL_condition: 0.1 # Courant-Friedrich-Levy condition for time integration. + +# Parameters related to the initial conditions +InitialConditions: + file_name: ./SN_feedback.hdf5 + smoothing_length_scaling: 1. + periodic: 1 # Are we running with periodic ICs? + +# Parameters for the stellar models +Stars: + resolution_eta: 1.2348 # Target smoothing length in units of the mean inter-particle separation (1.2348 == 48Ngbs with the cubic spline kernel). diff --git a/examples/SupernovaeFeedback/getGlass.sh b/examples/SupernovaeFeedback/getGlass.sh new file mode 100644 index 0000000000000000000000000000000000000000..d5c5f590ac37c9c9431d626a2ea61b0c12c1513c --- /dev/null +++ b/examples/SupernovaeFeedback/getGlass.sh @@ -0,0 +1,2 @@ +#!/bin/bash +wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/ICs/glassCube_64.hdf5 diff --git a/examples/SupernovaeFeedback/makeIC.py b/examples/SupernovaeFeedback/makeIC.py new file mode 100644 index 0000000000000000000000000000000000000000..8a0fca2bfd21f2cf1d5052660581d97184705480 --- /dev/null +++ b/examples/SupernovaeFeedback/makeIC.py @@ -0,0 +1,117 @@ +############################################################################### + # This file is part of SWIFT. + # Copyright (c) 2016 Matthieu Schaller (matthieu.schaller@durham.ac.uk) + # + # This program is free software: you can redistribute it and/or modify + # it under the terms of the GNU Lesser General Public License as published + # by the Free Software Foundation, either version 3 of the License, or + # (at your option) any later version. + # + # This program is distributed in the hope that it will be useful, + # but WITHOUT ANY WARRANTY; without even the implied warranty of + # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + # GNU General Public License for more details. + # + # You should have received a copy of the GNU Lesser General Public License + # along with this program. If not, see <http://www.gnu.org/licenses/>. + # + ############################################################################## + +import h5py +from numpy import * + +# Generates a swift IC file for the Sedov blast test in a periodic cubic box + +# Parameters +gamma = 5./3. # Gas adiabatic index +rho0 = 1. # Background density +P0 = 1.e-6 # Background pressure +E0= 1. # Energy of the explosion +N_inject = 15 # Number of particles in which to inject energy +fileName = "SN_feedback.hdf5" + +#--------------------------------------------------- +glass = h5py.File("glassCube_64.hdf5", "r") + +# Read particle positions and h from the glass +pos = glass["/PartType0/Coordinates"][:,:] +eps = 1e-6 +pos = (pos - pos.min()) / (pos.max() - pos.min() + eps) +h = glass["/PartType0/SmoothingLength"][:] * 0.3 * 3.3 + +numPart = size(h) +vol = 1. +Boxsize = 1. + +# Generate extra arrays +v = zeros((numPart, 3)) +ids = linspace(1, numPart, numPart) +m = zeros(numPart) +u = zeros(numPart) +r = zeros(numPart) + +r = sqrt((pos[:,0] - 0.5)**2 + (pos[:,1] - 0.5)**2 + (pos[:,2] - 0.5)**2) +m[:] = rho0 * vol / numPart +u[:] = P0 / (rho0 * (gamma - 1)) + +#-------------------------------------------------- + +star_pos = zeros((1, 3)) +star_pos[:,:] = 0.5 * Boxsize + +star_v = zeros((1, 3)) +star_v[:,:] = 0. + +# increase mass to keep it at center +star_m = 1e3 * array([rho0 * vol / numPart]) +star_ids = array([numPart + 1]) +star_h = array([h.max()]) + +#-------------------------------------------------- + +#File +file = h5py.File(fileName, 'w') + +# Header +grp = file.create_group("/Header") +grp.attrs["BoxSize"] = [Boxsize]*3 +grp.attrs["NumPart_Total"] = [numPart, 0, 0, 0, 1, 0] +grp.attrs["NumPart_Total_HighWord"] = [0, 0, 0, 0, 0, 0] +grp.attrs["NumPart_ThisFile"] = [numPart, 0, 0, 0, 1, 0] +grp.attrs["Time"] = 0.0 +grp.attrs["NumFilesPerSnapshot"] = 1 +grp.attrs["MassTable"] = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0] +grp.attrs["Flag_Entropy_ICs"] = 0 +grp.attrs["Dimension"] = 3 + +#Runtime parameters +grp = file.create_group("/RuntimePars") +grp.attrs["PeriodicBoundariesOn"] = 0 + +#Units +grp = file.create_group("/Units") +grp.attrs["Unit length in cgs (U_L)"] = 1. +grp.attrs["Unit mass in cgs (U_M)"] = 1. +grp.attrs["Unit time in cgs (U_t)"] = 1. +grp.attrs["Unit current in cgs (U_I)"] = 1. +grp.attrs["Unit temperature in cgs (U_T)"] = 1. + +#Particle group +grp = file.create_group("/PartType0") +grp.create_dataset('Coordinates', data=pos, dtype='d') +grp.create_dataset('Velocities', data=v, dtype='f') +grp.create_dataset('Masses', data=m, dtype='f') +grp.create_dataset('SmoothingLength', data=h, dtype='f') +grp.create_dataset('InternalEnergy', data=u, dtype='f') +grp.create_dataset('ParticleIDs', data=ids, dtype='L') + +# stellar group +grp = file.create_group("/PartType4") +grp.create_dataset("Coordinates", data=star_pos, dtype="d") +grp.create_dataset('Velocities', data=star_v, dtype='f') +grp.create_dataset('Masses', data=star_m, dtype='f') +grp.create_dataset('SmoothingLength', data=star_h, dtype='f') +grp.create_dataset('ParticleIDs', data=star_ids, dtype='L') + + +file.close() diff --git a/examples/SupernovaeFeedback/run.sh b/examples/SupernovaeFeedback/run.sh new file mode 100644 index 0000000000000000000000000000000000000000..d468bf02dc3228adbaf7c8ef0ee4513f4147e366 --- /dev/null +++ b/examples/SupernovaeFeedback/run.sh @@ -0,0 +1,19 @@ +#!/bin/bash + + # Generate the initial conditions if they are not present. +if [ ! -e glassCube_64.hdf5 ] +then + echo "Fetching initial glass file for the Supernovae feedback example..." + ./getGlass.sh +fi +if [ ! -e SN_feedback.hdf5 ] +then + echo "Generating initial conditions for the Supernovae feedback example..." + python makeIC.py +fi + +# Run SWIFT +../swift -g -b -s -S -t 4 SN_feedback.yml 2>&1 | tee output.log + +# Plot the solution +# TODO