diff --git a/configure.ac b/configure.ac
index 6ec1ece5ed2fef20ceab484eff4e09fce808e635..b494a969db2589039e6ab808b7619a553d3e9e60 100644
--- a/configure.ac
+++ b/configure.ac
@@ -270,6 +270,18 @@ elif test "$gravity_force_checks" != "no"; then
AC_DEFINE_UNQUOTED([SWIFT_GRAVITY_FORCE_CHECKS], [$enableval] ,[Enable gravity brute-force checks])
fi
+# Check whether we want to switch on glass making
+AC_ARG_ENABLE([glass-making],
+ [AS_HELP_STRING([--enable-glass-making],
+ [Activate the glass-making procedure by reversing the sign of gravity @<:@yes/no@:>@]
+ )],
+ [gravity_glass_making="$enableval"],
+ [gravity_glass_making="no"]
+)
+if test "$gravity_glass_making" == "yes"; then
+ AC_DEFINE([SWIFT_MAKE_GRAVITY_GLASS], 1, [Make the code run in a way to produce a glass file for gravity/cosmology])
+fi
+
# Check if we want to zero the gravity forces for all particles below some ID.
AC_ARG_ENABLE([no-gravity-below-id],
[AS_HELP_STRING([--enable-no-gravity-below-id],
@@ -1408,6 +1420,7 @@ AC_MSG_RESULT([
Gravity scheme : $with_gravity
Multipole order : $with_multipole_order
No gravity below ID : $no_gravity_below_id
+ Make gravity glass : $gravity_glass_making
External potential : $with_potential
Cooling function : $with_cooling
diff --git a/doc/RTD/source/GettingStarted/index.rst b/doc/RTD/source/GettingStarted/index.rst
index d15a8eee2f3b9089a1c8ee033f9aa3ee7ad92a5f..36de8ea740490c16bc9d6b69d871290e80dc2091 100644
--- a/doc/RTD/source/GettingStarted/index.rst
+++ b/doc/RTD/source/GettingStarted/index.rst
@@ -23,3 +23,4 @@ and keep on your desk.
parameter_file
what_about_mpi
running_on_large_systems
+ special_modes
diff --git a/doc/RTD/source/GettingStarted/parameter_file.rst b/doc/RTD/source/GettingStarted/parameter_file.rst
index 1bf4cbd3940a104c126aa256759edc24ca996338..550040ed25ec307633d6fade81eced58ed65a254 100644
--- a/doc/RTD/source/GettingStarted/parameter_file.rst
+++ b/doc/RTD/source/GettingStarted/parameter_file.rst
@@ -19,38 +19,51 @@ In the sections ``Snapshots`` and ``Statistics``, you can specify the options ``
The ``output_list_on`` enable or not the output list and ``output_list`` is the filename containing the output times.
With the file header, you can choose between writing redshifts, scale factors or times.
-Example of file containing with times (in internal units):
-::
- # Time
- 0.5
- 1.5
- 3.0
- 12.5
-
-Example of file with scale factors:
-::
- # Scale Factor
- 0.1
- 0.2
- 0.3
-
-Example of file with redshift:
-::
- # Redshift
- 20
- 15
- 10
- 5
+Example of file containing with times (in internal units)::
+
+ # Time
+ 0.5
+ 1.5
+ 3.0
+ 12.5
+
+Example of file with scale factors::
+
+ # Scale Factor
+ 0.1
+ 0.2
+ 0.3
+
+Example of file with redshift::
+
+ # Redshift
+ 20
+ 15
+ 10
+ 5
Output Selection
~~~~~~~~~~~~~~~~
-With SWIFT, you can select the particle fields to output in snapshot using the parameter file.
-In section ``SelectOutput``, you can remove a field by adding a parameter formatted in the
-following way ``field_parttype`` where ``field`` is the name of the field that you
-want to remove (e.g. ``Masses``) and ``parttype`` is the type of particles that
-contains this field (e.g. ``Gas``, ``DM`` or ``Star``). For a parameter, the only
-values accepted are 0 (skip this field when writing) or 1 (default, do not skip
-this field when writing).
+With SWIFT, you can select the particle fields to output in snapshot
+using the parameter file. In section ``SelectOutput``, you can remove
+a field by adding a parameter formatted in the following way
+``field_parttype`` where ``field`` is the name of the field that you
+want to remove (e.g. ``Masses``) and ``parttype`` is the type of
+particles that contains this field (e.g. ``Gas``, ``DM`` or ``Star``).
+For a parameter, the only values accepted are 0 (skip this field when
+writing) or 1 (default, do not skip this field when writing). By
+default all fields are written.
+
+This field is mostly used to remove unnecessary output by listing them
+with 0's. A classic use-case for this feature is a DM-only simulation
+(pure n-body) where all particles have the same mass. Outputing the
+mass field in the snapshots results in extra i/o time and unnecessary
+waste of disk space. The corresponding section of the ``yaml``
+parameter file would look like this::
+
+ SelectOutput:
+ Masses_DM: 0
-You can generate a ``yaml`` file containing all the possible fields with ``./swift -o output.yml``. By default, all the fields are written.
+You can generate a ``yaml`` file containing all the possible fields
+available for a given configuration of SWIFT by running ``./swift -o output.yml``.
diff --git a/doc/RTD/source/GettingStarted/special_modes.rst b/doc/RTD/source/GettingStarted/special_modes.rst
new file mode 100644
index 0000000000000000000000000000000000000000..636fc5e5d2237a6eaf4a2f4811cd17fa3e7010f5
--- /dev/null
+++ b/doc/RTD/source/GettingStarted/special_modes.rst
@@ -0,0 +1,43 @@
+.. Special modes
+ Matthieu Schaller, 20/08/2018
+
+Special modes
+=============
+
+SWIFT comes with a few special modes of operating to perform additional tasks.
+
+Static particles
+~~~~~~~~~~~~~~~~
+
+For some test problems it is convenient to have a set of particles that do not
+perceive any gravitational forces and just act as sources for the force
+calculation. This can be achieved by configuring SWIFT with the option
+``--enable-no-gravity-below-id=N``. This will zero the *accelerations* of all
+particles with ``id`` (strictly) lower than ``N`` at every time-step. Note that
+if these particles have an initial velocity they will keep moving at that
+speed.
+
+This will also naturally set their time-step to the maximal value
+(``TimeIntegration:dt_max``) set in the parameter file.
+
+A typical use-case for this feature is to study the evolution of one particle
+orbiting a static halo made of particles. This can be used to assess the
+quality of the gravity tree and time integration. As more particles are added
+to the halo, the orbits will get closer to the analytic solution as the noise
+in the sampling of the halo is reduced.
+
+Note also that this does not affect the hydrodynamic forces. This mode is
+purely designed for gravity-only accuracy tests.
+
+Gravity glasses
+~~~~~~~~~~~~~~~
+
+For many problems in cosmology, it is important to start a simulation with no
+initial noise in the particle distribution. Such a "glass" can be created by
+starting from a random distribution of particles and running with the sign of
+gravity reversed until all the particles reach a steady state. To run SWIFT in
+this mode, configure the code with ``--enable-glass-making``.
+
+Note that this will *not* generate the initial random distribution of
+particles. An initial condition file with particles still has to be provided.
+
diff --git a/examples/Gravity_glass/README b/examples/Gravity_glass/README
new file mode 100644
index 0000000000000000000000000000000000000000..2df2eb1e72cd2979750f9232936a6e183e8636da
--- /dev/null
+++ b/examples/Gravity_glass/README
@@ -0,0 +1,8 @@
+This example can be used to generate a glass file for gravity calculation.
+The makeIC.py script will generate a uniform Poisson distribution of particles
+in a cubic box with zero initial velocities.
+
+By running the code with the SWIFT configuration option --enable-glass-making
+the code will run with gravity as a repulsive force and the particles will
+moe towards a state of minimal energy. These glass files can be used to then
+start simulations with a minimal level of noise.
diff --git a/examples/UniformDMBox/makeIC.py b/examples/Gravity_glass/makeIC.py
similarity index 77%
rename from examples/UniformDMBox/makeIC.py
rename to examples/Gravity_glass/makeIC.py
index 8f3cd943b3cf19c4ae231d125c5ef97d076e0e8e..1a3fde9e2868c8881923fa61d1c308bca0f2f095 100644
--- a/examples/UniformDMBox/makeIC.py
+++ b/examples/Gravity_glass/makeIC.py
@@ -20,7 +20,7 @@
import h5py
import sys
-from numpy import *
+import numpy as np
# Generates a swift IC file containing a cartesian distribution of DM particles
# with a density of 1
@@ -30,7 +30,7 @@ periodic= 1 # 1 For periodic box
boxSize = 1.
rho = 1.
L = int(sys.argv[1]) # Number of particles along one axis
-fileName = "uniformDMBox_%d.hdf5"%L
+fileName = "uniform_DM_box.hdf5"
#---------------------------------------------------
numPart = L**3
@@ -69,27 +69,16 @@ grp.attrs["Unit temperature in cgs (U_T)"] = 1.
#Particle group
grp = file.create_group("/PartType1")
-v = zeros((numPart, 3))
-ds = grp.create_dataset('Velocities', (numPart, 3), 'f')
-ds[()] = v
-v = zeros(1)
+v = np.zeros((numPart, 3))
+ds = grp.create_dataset('Velocities', (numPart, 3), 'f', data=v)
-m = full((numPart, 1), mass)
-ds = grp.create_dataset('Masses', (numPart,1), 'f')
-ds[()] = m
-m = zeros(1)
+m = np.full((numPart, 1), mass)
+ds = grp.create_dataset('Masses', (numPart,1), 'f', data=m)
-ids = linspace(0, numPart, numPart, endpoint=False).reshape((numPart,1))
+ids = np.linspace(0, numPart, numPart, endpoint=False).reshape((numPart,1))
ds = grp.create_dataset('ParticleIDs', (numPart, 1), 'L')
ds[()] = ids + 1
-x = ids % L;
-y = ((ids - x) / L) % L;
-z = (ids - x - L * y) / L**2;
-coords = zeros((numPart, 3))
-coords[:,0] = z[:,0] * boxSize / L + boxSize / (2*L)
-coords[:,1] = y[:,0] * boxSize / L + boxSize / (2*L)
-coords[:,2] = x[:,0] * boxSize / L + boxSize / (2*L)
-ds = grp.create_dataset('Coordinates', (numPart, 3), 'd')
-ds[()] = coords
+coords = np.random.rand(numPart, 3) * boxSize
+ds = grp.create_dataset('Coordinates', (numPart, 3), 'd', data=coords)
file.close()
diff --git a/examples/Gravity_glass/uniform_DM_box.yml b/examples/Gravity_glass/uniform_DM_box.yml
new file mode 100644
index 0000000000000000000000000000000000000000..8f3ef6f025afb1a92320eeb702b50e8bf4befce6
--- /dev/null
+++ b/examples/Gravity_glass/uniform_DM_box.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
+
+# Let's overwrite G to make this more effective
+PhysicalConstants:
+ G: 1.
+
+# Parameters governing the time integration
+TimeIntegration:
+ time_begin: 0.
+ time_end: 100.
+ dt_min: 1e-6
+ dt_max: 1.
+
+Scheduler:
+ max_top_level_cells: 8
+
+# Parameters governing the snapshots
+Snapshots:
+ basename: uniform_DM_box
+ time_first: 0.
+ delta_time: 1.
+ compression: 4
+
+# Parameters for the self-gravity scheme
+Gravity:
+ eta: 0.025
+ theta: 0.3
+ mesh_side_length: 32
+ comoving_softening: 0.001
+ max_physical_softening: 0.001
+
+# Parameters governing the conserved quantities statistics
+Statistics:
+ delta_time: 0.1
+
+# Parameters related to the initial conditions
+InitialConditions:
+ file_name: ./uniform_DM_box.hdf5
diff --git a/examples/UniformDMBox/plot_gravity_checks.py b/examples/UniformDMBox/plot_gravity_checks.py
deleted file mode 100644
index 5efd5847ca9749fffaee48e586c0a1976fbac9d5..0000000000000000000000000000000000000000
--- a/examples/UniformDMBox/plot_gravity_checks.py
+++ /dev/null
@@ -1,219 +0,0 @@
-#!/usr/bin/env python
-
-import sys
-import glob
-import re
-import numpy as np
-import matplotlib.pyplot as plt
-
-params = {'axes.labelsize': 14,
-'axes.titlesize': 18,
-'font.size': 12,
-'legend.fontsize': 12,
-'xtick.labelsize': 14,
-'ytick.labelsize': 14,
-'text.usetex': True,
-'figure.figsize': (10, 10),
-'figure.subplot.left' : 0.06,
-'figure.subplot.right' : 0.99 ,
-'figure.subplot.bottom' : 0.06 ,
-'figure.subplot.top' : 0.985 ,
-'figure.subplot.wspace' : 0.14 ,
-'figure.subplot.hspace' : 0.14 ,
-'lines.markersize' : 6,
-'lines.linewidth' : 3.,
-'text.latex.unicode': True
-}
-plt.rcParams.update(params)
-plt.rc('font',**{'family':'sans-serif','sans-serif':['Times']})
-
-min_error = 1e-6
-max_error = 1e-1
-num_bins = 51
-
-# Construct the bins
-bin_edges = np.linspace(np.log10(min_error), np.log10(max_error), num_bins + 1)
-bin_size = (np.log10(max_error) - np.log10(min_error)) / num_bins
-bins = 0.5*(bin_edges[1:] + bin_edges[:-1])
-bin_edges = 10**bin_edges
-bins = 10**bins
-
-# Colours
-cols = ['b', 'g', 'r', 'm']
-
-# Time-step to plot
-step = int(sys.argv[1])
-
-# Find the files for the different expansion orders
-order_list = glob.glob("gravity_checks_step%d_order*.dat"%step)
-num_order = len(order_list)
-
-# Get the multipole orders
-order = np.zeros(num_order)
-for i in range(num_order):
- order[i] = int(order_list[i][26])
-
-# Start the plot
-plt.figure()
-
-# Get the Gadget-2 data if existing
-gadget2_file_list = glob.glob("forcetest_gadget2.txt")
-if len(gadget2_file_list) != 0:
-
- gadget2_data = np.loadtxt(gadget2_file_list[0])
- gadget2_ids = gadget2_data[:,0]
- gadget2_pos = gadget2_data[:,1:4]
- gadget2_a_exact = gadget2_data[:,4:7]
- gadget2_a_grav = gadget2_data[:, 7:10]
-
- # Sort stuff
- sort_index = np.argsort(gadget2_ids)
- gadget2_ids = gadget2_ids[sort_index]
- gadget2_pos = gadget2_pos[sort_index, :]
- gadget2_a_exact = gadget2_a_exact[sort_index, :]
- gadget2_a_grav = gadget2_a_grav[sort_index, :]
-
- # Compute the error norm
- diff = gadget2_a_exact - gadget2_a_grav
-
- norm_diff = np.sqrt(diff[:,0]**2 + diff[:,1]**2 + diff[:,2]**2)
- norm_a = np.sqrt(gadget2_a_exact[:,0]**2 + gadget2_a_exact[:,1]**2 + gadget2_a_exact[:,2]**2)
-
- norm_error = norm_diff / norm_a
- error_x = abs(diff[:,0]) / norm_a
- error_y = abs(diff[:,1]) / norm_a
- error_z = abs(diff[:,2]) / norm_a
-
- # Bin the error
- norm_error_hist,_ = np.histogram(norm_error, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_x_hist,_ = np.histogram(error_x, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_y_hist,_ = np.histogram(error_y, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_z_hist,_ = np.histogram(error_z, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
-
- norm_median = np.median(norm_error)
- median_x = np.median(error_x)
- median_y = np.median(error_y)
- median_z = np.median(error_z)
-
- norm_per95 = np.percentile(norm_error,95)
- per95_x = np.percentile(error_x,95)
- per95_y = np.percentile(error_y,95)
- per95_z = np.percentile(error_z,95)
-
- plt.subplot(221)
- plt.semilogx(bins, norm_error_hist, 'k--', label="Gadget-2")
- plt.plot([norm_median, norm_median], [2.7, 3], 'k-', lw=1)
- plt.plot([norm_per95, norm_per95], [2.7, 3], 'k:', lw=1)
- plt.subplot(222)
- plt.semilogx(bins, error_x_hist, 'k--', label="Gadget-2")
- plt.plot([median_x, median_x], [1.8, 2], 'k-', lw=1)
- plt.plot([per95_x, per95_x], [1.8, 2], 'k:', lw=1)
- plt.subplot(223)
- plt.semilogx(bins, error_y_hist, 'k--', label="Gadget-2")
- plt.plot([median_y, median_y], [1.8, 2], 'k-', lw=1)
- plt.plot([per95_y, per95_y], [1.8, 2], 'k:', lw=1)
- plt.subplot(224)
- plt.semilogx(bins, error_z_hist, 'k--', label="Gadget-2")
- plt.plot([median_z, median_z], [1.8, 2], 'k-', lw=1)
- plt.plot([per95_z, per95_z], [1.8, 2], 'k:', lw=1)
-
-
-# Plot the different histograms
-for i in range(num_order-1, -1, -1):
- data = np.loadtxt(order_list[i])
- ids = data[:,0]
- pos = data[:,1:4]
- a_exact = data[:,4:7]
- a_grav = data[:, 7:10]
-
- # Sort stuff
- sort_index = np.argsort(ids)
- ids = ids[sort_index]
- pos = pos[sort_index, :]
- a_exact = a_exact[sort_index, :]
- a_grav = a_grav[sort_index, :]
-
- # Cross-checks
- if not np.array_equal(ids, gadget2_ids):
- print "Comparing different IDs !"
-
- if not np.array_equal(pos, gadget2_pos):
- print "Comparing different positions ! max difference:", np.max(pos - gadget2_pos)
-
- if not np.array_equal(a_exact, gadget2_a_exact):
- print "Comparing different exact accelerations ! max difference:", np.max(a_exact - gadget2_a_exact)
-
-
- # Compute the error norm
- diff = a_exact - a_grav
-
- norm_diff = np.sqrt(diff[:,0]**2 + diff[:,1]**2 + diff[:,2]**2)
- norm_a = np.sqrt(a_exact[:,0]**2 + a_exact[:,1]**2 + a_exact[:,2]**2)
-
- norm_error = norm_diff / norm_a
- error_x = abs(diff[:,0]) / norm_a
- error_y = abs(diff[:,1]) / norm_a
- error_z = abs(diff[:,2]) / norm_a
-
- # Bin the error
- norm_error_hist,_ = np.histogram(norm_error, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_x_hist,_ = np.histogram(error_x, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_y_hist,_ = np.histogram(error_y, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
- error_z_hist,_ = np.histogram(error_z, bins=bin_edges, density=False) / (np.size(norm_error) * bin_size)
-
- norm_median = np.median(norm_error)
- median_x = np.median(error_x)
- median_y = np.median(error_y)
- median_z = np.median(error_z)
-
- norm_per95 = np.percentile(norm_error,95)
- per95_x = np.percentile(error_x,95)
- per95_y = np.percentile(error_y,95)
- per95_z = np.percentile(error_z,95)
-
- plt.subplot(221)
- plt.semilogx(bins, norm_error_hist, color=cols[i],label="SWIFT m-poles order %d"%order[i])
- plt.plot([norm_median, norm_median], [2.7, 3],'-', color=cols[i], lw=1)
- plt.plot([norm_per95, norm_per95], [2.7, 3],':', color=cols[i], lw=1)
- plt.subplot(222)
- plt.semilogx(bins, error_x_hist, color=cols[i],label="SWIFT m-poles order %d"%order[i])
- plt.plot([median_x, median_x], [1.8, 2],'-', color=cols[i], lw=1)
- plt.plot([per95_x, per95_x], [1.8, 2],':', color=cols[i], lw=1)
- plt.subplot(223)
- plt.semilogx(bins, error_y_hist, color=cols[i],label="SWIFT m-poles order %d"%order[i])
- plt.plot([median_y, median_y], [1.8, 2],'-', color=cols[i], lw=1)
- plt.plot([per95_y, per95_y], [1.8, 2],':', color=cols[i], lw=1)
- plt.subplot(224)
- plt.semilogx(bins, error_z_hist, color=cols[i],label="SWIFT m-poles order %d"%order[i])
- plt.plot([median_z, median_z], [1.8, 2],'-', color=cols[i], lw=1)
- plt.plot([per95_z, per95_z], [1.8, 2],':', color=cols[i], lw=1)
-
-
-plt.subplot(221)
-plt.xlabel("$|\delta \overrightarrow{a}|/|\overrightarrow{a}_{exact}|$")
-plt.ylabel("Density")
-plt.xlim(min_error, 2*max_error)
-plt.ylim(0,3)
-plt.legend(loc="upper left")
-plt.subplot(222)
-plt.xlabel("$\delta a_x/|\overrightarrow{a}_{exact}|$")
-plt.ylabel("Density")
-plt.xlim(min_error, 2*max_error)
-plt.ylim(0,2)
-#plt.legend(loc="center left")
-plt.subplot(223)
-plt.xlabel("$\delta a_y/|\overrightarrow{a}_{exact}|$")
-plt.ylabel("Density")
-plt.xlim(min_error, 2*max_error)
-plt.ylim(0,2)
-#plt.legend(loc="center left")
-plt.subplot(224)
-plt.xlabel("$\delta a_z/|\overrightarrow{a}_{exact}|$")
-plt.ylabel("Density")
-plt.xlim(min_error, 2*max_error)
-plt.ylim(0,2)
-#plt.legend(loc="center left")
-
-
-
-plt.savefig("gravity_checks_step%d.png"%step)
diff --git a/examples/UniformDMBox/uniformBox.yml b/examples/UniformDMBox/uniformBox.yml
deleted file mode 100644
index 1abb256671f1cc8c87daa711bd63f7ea6abdbbab..0000000000000000000000000000000000000000
--- a/examples/UniformDMBox/uniformBox.yml
+++ /dev/null
@@ -1,38 +0,0 @@
-# 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
-
-# Parameters governing the time integration
-TimeIntegration:
- time_begin: 0. # The starting time of the simulation (in internal units).
- time_end: 100. # The end time of the simulation (in internal units).
- dt_min: 1e-6 # The minimal time-step size of the simulation (in internal units).
- dt_max: 1. # The maximal time-step size of the simulation (in internal units).
-
-Scheduler:
- max_top_level_cells: 8
- cell_split_size: 50
-
-# Parameters governing the snapshots
-Snapshots:
- basename: uniformDMBox # Common part of the name of output files
- time_first: 0. # Time of the first output (in internal units)
- delta_time: 10. # Time difference between consecutive outputs (in internal units)
-
-# Parameters for the self-gravity scheme
-Gravity:
- eta: 0.025 # Constant dimensionless multiplier for time integration.
- theta: 0.8 # Opening angle (Multipole acceptance criterion)
- epsilon: 0.01 # Softening length (in internal units).
-
-# Parameters governing the conserved quantities statistics
-Statistics:
- delta_time: 5. # Time between statistics output
-
-# Parameters related to the initial conditions
-InitialConditions:
- file_name: ./uniformDMBox_16.hdf5 # The file to read
diff --git a/src/runner.c b/src/runner.c
index 6093f1194d24cb6a0af1b2fe31d04be3196930e4..7e67c6e943f7493bf2e472da622291653cfe9001 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -1757,6 +1757,14 @@ void runner_do_end_force(struct runner *r, struct cell *c, int timer) {
/* Finish the force calculation */
gravity_end_force(gp, G_newton, potential_normalisation, periodic);
+#ifdef SWIFT_MAKE_GRAVITY_GLASS
+
+ /* Negate the gravity forces */
+ gp->a_grav[0] *= -1.f;
+ gp->a_grav[1] *= -1.f;
+ gp->a_grav[2] *= -1.f;
+#endif
+
#ifdef SWIFT_NO_GRAVITY_BELOW_ID
/* Get the ID of the gpart */