Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
SWIFT
SWIFTsim
Commits
cb952932
Commit
cb952932
authored
Oct 07, 2018
by
Matthieu Schaller
Browse files
Added a new example: A small cosmological volume with some cooling on.
parent
c8f45163
Changes
10
Hide whitespace changes
Inline
Side-by-side
examples/SmallCosmoVolume/small_cosmo_volume.yml
View file @
cb952932
...
...
@@ -56,8 +56,3 @@ InitialConditions:
cleanup_velocity_factors
:
1
generate_gas_in_ics
:
1
# Generate gas particles from the DM-only ICs
cleanup_smoothing_lengths
:
1
# Since we generate gas, make use of the (expensive) cleaning-up procedure.
# Constant lambda cooling function
LambdaCooling
:
lambda_cgs
:
1e-22
# Cooling rate (in cgs units)
cooling_tstep_mult
:
1.0
# Dimensionless pre-factor for the time-step condition
examples/SmallCosmoVolume_cooling/README
0 → 100644
View file @
cb952932
Small LCDM cosmological simulation generated by C. Power. Cosmology
is WMAP9 and the box is 100Mpc/h in size with 64^3 particles.
We use a softening length of 1/25th of the mean inter-particle separation.
The ICs have been generated to run with Gadget-2 so we need to switch
on the options to cancel the h-factors and a-factors at reading time.
We generate gas from the ICs using SWIFT's internal mechanism and set the
temperature to the expected gas temperature at this redshift.
The 'plotTempEvolution.py' plots the temperature evolution of the gas
in the simulated volume.
MD5 checksum of the ICs:
08736c3101fd738e22f5159f78e6022b small_cosmo_volume.hdf5
examples/SmallCosmoVolume_cooling/getIC.sh
0 → 100755
View file @
cb952932
#!/bin/bash
wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/ICs/small_cosmo_volume.hdf5
examples/SmallCosmoVolume_cooling/plotTempEvolution.py
0 → 100644
View file @
cb952932
################################################################################
# This file is part of SWIFT.
# Copyright (c) 2018 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/>.
#
################################################################################
# Computes the temperature evolution of the gas in a cosmological box
# Physical constants needed for internal energy to temperature conversion
k_in_J_K
=
1.38064852e-23
mH_in_kg
=
1.6737236e-27
# Number of snapshots generated
n_snapshots
=
160
import
matplotlib
matplotlib
.
use
(
"Agg"
)
from
pylab
import
*
import
h5py
import
os.path
# Plot parameters
params
=
{
'axes.labelsize'
:
10
,
'axes.titlesize'
:
10
,
'font.size'
:
9
,
'legend.fontsize'
:
9
,
'xtick.labelsize'
:
10
,
'ytick.labelsize'
:
10
,
'text.usetex'
:
True
,
'figure.figsize'
:
(
3.15
,
3.15
),
'figure.subplot.left'
:
0.14
,
'figure.subplot.right'
:
0.99
,
'figure.subplot.bottom'
:
0.12
,
'figure.subplot.top'
:
0.99
,
'figure.subplot.wspace'
:
0.15
,
'figure.subplot.hspace'
:
0.12
,
'lines.markersize'
:
6
,
'lines.linewidth'
:
2.
,
'text.latex.unicode'
:
True
}
rcParams
.
update
(
params
)
rc
(
'font'
,
**
{
'family'
:
'sans-serif'
,
'sans-serif'
:[
'Times'
]})
# Read the simulation data
sim
=
h5py
.
File
(
"snap_0000.hdf5"
,
"r"
)
boxSize
=
sim
[
"/Header"
].
attrs
[
"BoxSize"
][
0
]
time
=
sim
[
"/Header"
].
attrs
[
"Time"
][
0
]
scheme
=
sim
[
"/HydroScheme"
].
attrs
[
"Scheme"
][
0
]
kernel
=
sim
[
"/HydroScheme"
].
attrs
[
"Kernel function"
][
0
]
neighbours
=
sim
[
"/HydroScheme"
].
attrs
[
"Kernel target N_ngb"
][
0
]
eta
=
sim
[
"/HydroScheme"
].
attrs
[
"Kernel eta"
][
0
]
alpha
=
sim
[
"/HydroScheme"
].
attrs
[
"Alpha viscosity"
][
0
]
H_mass_fraction
=
sim
[
"/HydroScheme"
].
attrs
[
"Hydrogen mass fraction"
][
0
]
H_transition_temp
=
sim
[
"/HydroScheme"
].
attrs
[
"Hydrogen ionization transition temperature"
][
0
]
T_initial
=
sim
[
"/HydroScheme"
].
attrs
[
"Initial temperature"
][
0
]
T_minimal
=
sim
[
"/HydroScheme"
].
attrs
[
"Minimal temperature"
][
0
]
git
=
sim
[
"Code"
].
attrs
[
"Git Revision"
]
cooling_model
=
sim
[
"/SubgridScheme"
].
attrs
[
"Cooling Model"
]
if
cooling_model
==
"Constant Lambda"
:
Lambda
=
sim
[
"/SubgridScheme"
].
attrs
[
"Lambda [cgs]"
][
0
]
# Cosmological parameters
H_0
=
sim
[
"/Cosmology"
].
attrs
[
"H0 [internal units]"
][
0
]
gas_gamma
=
sim
[
"/HydroScheme"
].
attrs
[
"Adiabatic index"
][
0
]
unit_length_in_cgs
=
sim
[
"/Units"
].
attrs
[
"Unit length in cgs (U_L)"
]
unit_mass_in_cgs
=
sim
[
"/Units"
].
attrs
[
"Unit mass in cgs (U_M)"
]
unit_time_in_cgs
=
sim
[
"/Units"
].
attrs
[
"Unit time in cgs (U_t)"
]
unit_length_in_si
=
0.01
*
unit_length_in_cgs
unit_mass_in_si
=
0.001
*
unit_mass_in_cgs
unit_time_in_si
=
unit_time_in_cgs
# Primoridal ean molecular weight as a function of temperature
def
mu
(
T
,
H_frac
=
H_mass_fraction
,
T_trans
=
H_transition_temp
):
if
T
>
T_trans
:
return
4.
/
(
8.
-
5.
*
(
1.
-
H_frac
))
else
:
return
4.
/
(
1.
+
3.
*
H_frac
)
# Temperature of some primoridal gas with a given internal energy
def
T
(
u
,
H_frac
=
H_mass_fraction
,
T_trans
=
H_transition_temp
):
T_over_mu
=
(
gas_gamma
-
1.
)
*
u
*
mH_in_kg
/
k_in_J_K
ret
=
np
.
ones
(
np
.
size
(
u
))
*
T_trans
# Enough energy to be ionized?
mask_ionized
=
(
T_over_mu
>
(
T_trans
+
1
)
/
mu
(
T_trans
+
1
,
H_frac
,
T_trans
))
if
np
.
sum
(
mask_ionized
)
>
0
:
ret
[
mask_ionized
]
=
T_over_mu
[
mask_ionized
]
*
mu
(
T_trans
*
10
,
H_frac
,
T_trans
)
# Neutral gas?
mask_neutral
=
(
T_over_mu
<
(
T_trans
-
1
)
/
mu
((
T_trans
-
1
),
H_frac
,
T_trans
))
if
np
.
sum
(
mask_neutral
)
>
0
:
ret
[
mask_neutral
]
=
T_over_mu
[
mask_neutral
]
*
mu
(
0
,
H_frac
,
T_trans
)
return
ret
z
=
np
.
zeros
(
n_snapshots
)
a
=
np
.
zeros
(
n_snapshots
)
T_mean
=
np
.
zeros
(
n_snapshots
)
T_std
=
np
.
zeros
(
n_snapshots
)
T_log_mean
=
np
.
zeros
(
n_snapshots
)
T_log_std
=
np
.
zeros
(
n_snapshots
)
T_median
=
np
.
zeros
(
n_snapshots
)
T_min
=
np
.
zeros
(
n_snapshots
)
T_max
=
np
.
zeros
(
n_snapshots
)
# Loop over all the snapshots
for
i
in
range
(
n_snapshots
):
sim
=
h5py
.
File
(
"snap_%04d.hdf5"
%
i
,
"r"
)
z
[
i
]
=
sim
[
"/Cosmology"
].
attrs
[
"Redshift"
][
0
]
a
[
i
]
=
sim
[
"/Cosmology"
].
attrs
[
"Scale-factor"
][
0
]
u
=
sim
[
"/PartType0/InternalEnergy"
][:]
# Compute the temperature
u
*=
(
unit_length_in_si
**
2
/
unit_time_in_si
**
2
)
u
/=
a
[
i
]
**
(
3
*
(
gas_gamma
-
1.
))
Temp
=
T
(
u
)
# Gather statistics
T_median
[
i
]
=
np
.
median
(
Temp
)
T_mean
[
i
]
=
Temp
.
mean
()
T_std
[
i
]
=
Temp
.
std
()
T_log_mean
[
i
]
=
np
.
log10
(
Temp
).
mean
()
T_log_std
[
i
]
=
np
.
log10
(
Temp
).
std
()
T_min
[
i
]
=
Temp
.
min
()
T_max
[
i
]
=
Temp
.
max
()
# CMB evolution
a_evol
=
np
.
logspace
(
-
3
,
0
,
60
)
T_cmb
=
(
1.
/
a_evol
)
**
2
*
2.72
# Plot the interesting quantities
figure
()
subplot
(
111
,
xscale
=
"log"
,
yscale
=
"log"
)
fill_between
(
a
,
T_mean
-
T_std
,
T_mean
+
T_std
,
color
=
'C0'
,
alpha
=
0.1
)
plot
(
a
,
T_max
,
ls
=
'-.'
,
color
=
'C0'
,
lw
=
1.
,
label
=
"${
\\
rm max}~T$"
)
plot
(
a
,
T_min
,
ls
=
':'
,
color
=
'C0'
,
lw
=
1.
,
label
=
"${
\\
rm min}~T$"
)
plot
(
a
,
T_mean
,
color
=
'C0'
,
label
=
"${
\\
rm mean}~T$"
,
lw
=
1.5
)
fill_between
(
a
,
10
**
(
T_log_mean
-
T_log_std
),
10
**
(
T_log_mean
+
T_log_std
),
color
=
'C1'
,
alpha
=
0.1
)
plot
(
a
,
10
**
T_log_mean
,
color
=
'C1'
,
label
=
"${
\\
rm mean}~{
\\
rm log} T$"
,
lw
=
1.5
)
plot
(
a
,
T_median
,
color
=
'C2'
,
label
=
"${
\\
rm median}~T$"
,
lw
=
1.5
)
legend
(
loc
=
"upper left"
,
frameon
=
False
,
handlelength
=
1.5
)
# Cooling model
if
cooling_model
==
"Constant Lambda"
:
text
(
1e-2
,
6e4
,
"$\Lambda_{
\\
rm const} = %.1f
\\
times10^{%d}~[
\\
rm{cgs}]$"
%
(
Lambda
/
10.
**
(
int
(
log10
(
Lambda
))),
log10
(
Lambda
)),
fontsize
=
8
)
else
:
text
(
1e-2
,
6e4
,
"No cooling"
)
# Expected lines
plot
([
1e-10
,
1e10
],
[
H_transition_temp
,
H_transition_temp
],
'k--'
,
lw
=
0.5
,
alpha
=
0.7
)
text
(
2.5e-2
,
H_transition_temp
*
1.07
,
"$T_{
\\
rm HII
\\
rightarrow HI}$"
,
va
=
"bottom"
,
alpha
=
0.7
,
fontsize
=
8
)
plot
([
1e-10
,
1e10
],
[
T_minimal
,
T_minimal
],
'k--'
,
lw
=
0.5
,
alpha
=
0.7
)
text
(
1e-2
,
T_minimal
*
0.8
,
"$T_{
\\
rm min}$"
,
va
=
"top"
,
alpha
=
0.7
,
fontsize
=
8
)
plot
(
a_evol
,
T_cmb
,
'k--'
,
lw
=
0.5
,
alpha
=
0.7
)
text
(
a_evol
[
20
],
T_cmb
[
20
]
*
0.55
,
"$(1+z)^2
\\
times T_{
\\
rm CMB,0}$"
,
rotation
=-
34
,
alpha
=
0.7
,
fontsize
=
8
,
va
=
"top"
,
bbox
=
dict
(
facecolor
=
'w'
,
edgecolor
=
'none'
,
pad
=
1.0
,
alpha
=
0.9
))
redshift_ticks
=
np
.
array
([
0.
,
1.
,
2.
,
5.
,
10.
,
20.
,
50.
,
100.
])
redshift_labels
=
[
"$0$"
,
"$1$"
,
"$2$"
,
"$5$"
,
"$10$"
,
"$20$"
,
"$50$"
,
"$100$"
]
a_ticks
=
1.
/
(
redshift_ticks
+
1.
)
xticks
(
a_ticks
,
redshift_labels
)
minorticks_off
()
xlabel
(
"${
\\
rm Redshift}~z$"
,
labelpad
=
0
)
ylabel
(
"${
\\
rm Temperature}~T~[{
\\
rm K}]$"
,
labelpad
=
0
)
xlim
(
9e-3
,
1.1
)
ylim
(
20
,
2.5e7
)
savefig
(
"Temperature_evolution.png"
,
dpi
=
200
)
examples/SmallCosmoVolume_cooling/run.sh
0 → 100755
View file @
cb952932
#!/bin/bash
# Generate the initial conditions if they are not present.
if
[
!
-e
small_cosmo_volume.hdf5
]
then
echo
"Fetching initial conditions for the small cosmological volume example..."
./getIC.sh
fi
# Run SWIFT
../swift
-c
-s
-G
-C
-t
8 small_cosmo_volume.yml 2>&1 |
tee
output.log
# Plot the temperature evolution
python plotTempEvolution.py
examples/SmallCosmoVolume_cooling/small_cosmo_volume.yml
0 → 100644
View file @
cb952932
# Define the system of units to use internally.
InternalUnitSystem
:
UnitMass_in_cgs
:
1.98848e43
# 10^10 M_sun
UnitLength_in_cgs
:
3.08567758e24
# 1 Mpc
UnitVelocity_in_cgs
:
1e5
# 1 km/s
UnitCurrent_in_cgs
:
1
# Amperes
UnitTemp_in_cgs
:
1
# Kelvin
Cosmology
:
# WMAP9 cosmology
Omega_m
:
0.276
Omega_lambda
:
0.724
Omega_b
:
0.0455
h
:
0.703
a_begin
:
0.019607843
# z_ini = 50.
a_end
:
1.0
# z_end = 0.
# Parameters governing the time integration
TimeIntegration
:
dt_min
:
1e-6
dt_max
:
1e-2
# Parameters for the self-gravity scheme
Gravity
:
eta
:
0.025
theta
:
0.3
comoving_softening
:
0.0889
# 1/25th of the mean inter-particle separation: 88.9 kpc
max_physical_softening
:
0.0889
# 1/25th of the mean inter-particle separation: 88.9 kpc
mesh_side_length
:
64
# Parameters of the hydro scheme
SPH
:
resolution_eta
:
1.2348
# "48 Ngb" with the cubic spline kernel
CFL_condition
:
0.1
initial_temperature
:
7075.
# (1 + z_ini)^2 * 2.72K
minimal_temperature
:
100.
# Parameters governing the snapshots
Snapshots
:
basename
:
snap
delta_time
:
1.02
scale_factor_first
:
0.02
# Parameters governing the conserved quantities statistics
Statistics
:
delta_time
:
1.02
scale_factor_first
:
0.02
Scheduler
:
max_top_level_cells
:
8
cell_split_size
:
50
# Parameters related to the initial conditions
InitialConditions
:
file_name
:
small_cosmo_volume.hdf5
cleanup_h_factors
:
1
cleanup_velocity_factors
:
1
generate_gas_in_ics
:
1
# Generate gas particles from the DM-only ICs
cleanup_smoothing_lengths
:
1
# Since we generate gas, make use of the (expensive) cleaning-up procedure.
# Constant lambda cooling function
LambdaCooling
:
lambda_cgs
:
1e-30
# Cooling rate (in cgs units)
cooling_tstep_mult
:
1.0
# Dimensionless pre-factor for the time-step condition
src/cooling/const_lambda/cooling_io.h
View file @
cb952932
...
...
@@ -35,9 +35,10 @@
* @param h_grpsph The HDF5 group in which to write
*/
__attribute__
((
always_inline
))
INLINE
static
void
cooling_write_flavour
(
hid_t
h_grp
sph
)
{
hid_t
h_grp
,
const
struct
cooling_function_data
*
cooling
)
{
io_write_attribute_s
(
h_grpsph
,
"Cooling Model"
,
"Constant Lambda"
);
io_write_attribute_s
(
h_grp
,
"Cooling Model"
,
"Constant Lambda"
);
io_write_attribute_d
(
h_grp
,
"Lambda [cgs]"
,
cooling
->
lambda_cgs
);
}
#endif
...
...
src/parallel_io.c
View file @
cb952932
...
...
@@ -1020,7 +1020,7 @@ void prepare_file(struct engine* e, const char* baseName, long long N_total[6],
h_grp
=
H5Gcreate
(
h_file
,
"/SubgridScheme"
,
H5P_DEFAULT
,
H5P_DEFAULT
,
H5P_DEFAULT
);
if
(
h_grp
<
0
)
error
(
"Error while creating subgrid group"
);
cooling_write_flavour
(
h_grp
);
cooling_write_flavour
(
h_grp
,
e
->
cooling_func
);
chemistry_write_flavour
(
h_grp
);
H5Gclose
(
h_grp
);
...
...
src/serial_io.c
View file @
cb952932
...
...
@@ -870,7 +870,7 @@ void write_output_serial(struct engine* e, const char* baseName,
h_grp
=
H5Gcreate
(
h_file
,
"/SubgridScheme"
,
H5P_DEFAULT
,
H5P_DEFAULT
,
H5P_DEFAULT
);
if
(
h_grp
<
0
)
error
(
"Error while creating subgrid group"
);
cooling_write_flavour
(
h_grp
);
cooling_write_flavour
(
h_grp
,
e
->
cooling_func
);
chemistry_write_flavour
(
h_grp
);
H5Gclose
(
h_grp
);
...
...
src/single_io.c
View file @
cb952932
...
...
@@ -720,7 +720,7 @@ void write_output_single(struct engine* e, const char* baseName,
h_grp
=
H5Gcreate
(
h_file
,
"/SubgridScheme"
,
H5P_DEFAULT
,
H5P_DEFAULT
,
H5P_DEFAULT
);
if
(
h_grp
<
0
)
error
(
"Error while creating subgrid group"
);
cooling_write_flavour
(
h_grp
);
cooling_write_flavour
(
h_grp
,
e
->
cooling_func
);
chemistry_write_flavour
(
h_grp
);
H5Gclose
(
h_grp
);
...
...
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new file
.
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment