Commit 7b68d177 authored by James Willis's avatar James Willis
Browse files

Merge branch 'doself2-vectorisation' into dopair2-vectorisation

parents 9a0ff76d 8a07c59a
......@@ -47,6 +47,8 @@ tests/brute_force_27_perturbed.dat
tests/swift_dopair_27_perturbed.dat
tests/brute_force_125_standard.dat
tests/swift_dopair_125_standard.dat
tests/brute_force_125_perturbed.dat
tests/swift_dopair_125_perturbed.dat
tests/testGreetings
tests/testReading
tests/input.hdf5
......@@ -65,6 +67,7 @@ tests/parser_output.yml
tests/test27cells.sh
tests/test27cellsPerturbed.sh
tests/test125cells.sh
tests/test125cellsPerturbed.sh
tests/testPair.sh
tests/testPairPerturbed.sh
tests/testParser.sh
......
......@@ -853,6 +853,7 @@ AC_CONFIG_FILES([tests/testPairPerturbed.sh], [chmod +x tests/testPairPerturbed.
AC_CONFIG_FILES([tests/test27cells.sh], [chmod +x tests/test27cells.sh])
AC_CONFIG_FILES([tests/test27cellsPerturbed.sh], [chmod +x tests/test27cellsPerturbed.sh])
AC_CONFIG_FILES([tests/test125cells.sh], [chmod +x tests/test125cells.sh])
AC_CONFIG_FILES([tests/test125cellsPerturbed.sh], [chmod +x tests/test125cellsPerturbed.sh])
AC_CONFIG_FILES([tests/testParser.sh], [chmod +x tests/testParser.sh])
# Save the compilation options
......
......@@ -31,8 +31,6 @@ Gravity:
eta: 0.025 # Constant dimensionless multiplier for time integration.
theta: 0.7 # Opening angle (Multipole acceptance criterion)
epsilon: 0.0001 # Softening length (in internal units).
a_smooth: 1000.
r_cut: 4.
# Parameters for the hydrodynamics scheme
SPH:
......
###############################################################################
# This file is part of SWIFT.
# Copyright (c) 2016 Stefan Arridge (stefan.arridge@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 numpy as np
import h5py as h5
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use("Agg")
from pylab import *
import sys
#for the plotting
......@@ -46,7 +67,8 @@ for i in range(n_snaps):
f = h5.File(filename,'r')
coords_dset = f["PartType0/Coordinates"]
coords = np.array(coords_dset)
#translate coords by centre of box
#translate coords by centre of box
header = f["Header"]
snap_time = header.attrs["Time"]
snap_time_cgs = snap_time * unit_time_cgs
......@@ -63,58 +85,46 @@ for i in range(n_snaps):
bin_width = bin_edges[1] - bin_edges[0]
hist = np.histogram(r,bins = bin_edges)[0] # number of particles in each bin
#find the mass in each radial bin
#find the mass in each radial bin
mass_dset = f["PartType0/Masses"]
#mass of each particles should be equal
#mass of each particles should be equal
part_mass = np.array(mass_dset)[0]
part_mass_cgs = part_mass * unit_mass_cgs
part_mass_over_virial_mass = part_mass_cgs / M_vir_cgs
mass_hist = hist * part_mass_over_virial_mass
radial_bin_mids = np.linspace(bin_width/2.,max_r - bin_width/2.,n_radial_bins)
#volume in each radial bin
#volume in each radial bin
volume = 4.*np.pi * radial_bin_mids**2 * bin_width
#now divide hist by the volume so we have a density in each bin
#now divide hist by the volume so we have a density in each bin
density = mass_hist / volume
##read the densities
# density_dset = f["PartType0/Density"]
# density = np.array(density_dset)
# density_cgs = density * unit_mass_cgs / unit_length_cgs**3
# rho = density_cgs * r_vir_cgs**3 / M_vir_cgs
t = np.linspace(10./n_radial_bins,10.0,1000)
rho_analytic = t**(-2)/(4.*np.pi)
#calculate cooling radius
#r_cool_over_r_vir = np.sqrt((2.*(gamma - 1.)*lambda_cgs*M_vir_cgs*X_H**2)/(4.*np.pi*CONST_m_H_CGS**2*v_c_cgs**2*r_vir_cgs**3))*np.sqrt(snap_time_cgs)
#initial analytic density profile
#initial analytic density profile
if (i == 0):
r_0 = radial_bin_mids[0]
rho_0 = density[0]
rho_analytic_init = rho_0 * (radial_bin_mids/r_0)**(-2)
plt.plot(radial_bin_mids,density/rho_analytic_init,'ko',label = "Average density of shell")
#plt.plot(t,rho_analytic,label = "Initial analytic density profile"
plt.xlabel(r"$r / r_{vir}$")
plt.ylabel(r"$\rho / \rho_{init})$")
plt.title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
#plt.ylim((1.e-2,1.e1))
#plt.plot((r_cool_over_r_vir,r_cool_over_r_vir),(0,20),'r',label = "Cooling radius")
plt.xlim((radial_bin_mids[0],max_r))
plt.ylim((0,20))
plt.plot((0,max_r),(1,1))
#plt.xscale('log')
#plt.yscale('log')
plt.legend(loc = "upper right")
figure()
plot(radial_bin_mids,density/rho_analytic_init,'ko',label = "Average density of shell")
#plot(t,rho_analytic,label = "Initial analytic density profile")
xlabel(r"$r / r_{vir}$")
ylabel(r"$\rho / \rho_{init}$")
title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
xlim((radial_bin_mids[0],max_r))
ylim((0,2))
plot((0,max_r),(1,1))
legend(loc = "upper right")
plot_filename = "./plots/density_profile/density_profile_%03d.png" %i
plt.savefig(plot_filename,format = "png")
plt.close()
savefig(plot_filename,format = "png")
close()
###############################################################################
# This file is part of SWIFT.
# Copyright (c) 2016 Stefan Arridge (stefan.arridge@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 numpy as np
import h5py as h5
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use("Agg")
from pylab import *
import sys
def do_binning(x,y,x_bin_edges):
......@@ -48,8 +69,6 @@ unit_velocity_cgs = float(params.attrs["InternalUnitSystem:UnitVelocity_in_cgs"]
unit_time_cgs = unit_length_cgs / unit_velocity_cgs
v_c = float(params.attrs["IsothermalPotential:vrot"])
v_c_cgs = v_c * unit_velocity_cgs
#lambda_cgs = float(params.attrs["LambdaCooling:lambda_cgs"])
#X_H = float(params.attrs["LambdaCooling:hydrogen_mass_abundance"])
header = f["Header"]
N = header.attrs["NumPart_Total"][0]
box_centre = np.array(header.attrs["BoxSize"])
......@@ -64,7 +83,8 @@ for i in range(n_snaps):
f = h5.File(filename,'r')
coords_dset = f["PartType0/Coordinates"]
coords = np.array(coords_dset)
#translate coords by centre of box
#translate coords by centre of box
header = f["Header"]
snap_time = header.attrs["Time"]
snap_time_cgs = snap_time * unit_time_cgs
......@@ -75,11 +95,11 @@ for i in range(n_snaps):
radius_cgs = radius*unit_length_cgs
radius_over_virial_radius = radius_cgs / r_vir_cgs
#get the internal energies
#get the internal energies
u_dset = f["PartType0/InternalEnergy"]
u = np.array(u_dset)
#make dimensionless
#make dimensionless
u /= v_c**2/(2. * (gamma - 1.))
r = radius_over_virial_radius
......@@ -90,21 +110,16 @@ for i in range(n_snaps):
radial_bin_mids = np.linspace(bin_widths / 2. , max_r - bin_widths / 2. , n_radial_bins)
binned_u = u_totals / hist
#calculate cooling radius
#r_cool_over_r_vir = np.sqrt((2.*(gamma - 1.)*lambda_cgs*M_vir_cgs*X_H**2)/(4.*np.pi*CONST_m_H_CGS**2*v_c_cgs**2*r_vir_cgs**3))*np.sqrt(snap_time_cgs)
plt.plot(radial_bin_mids,binned_u,'ko',label = "Numerical solution")
#plt.plot((0,1),(1,1),label = "Analytic Solution")
#plt.plot((r_cool_over_r_vir,r_cool_over_r_vir),(0,2),'r',label = "Cooling radius")
plt.legend(loc = "lower right")
plt.xlabel(r"$r / r_{vir}$")
plt.ylabel(r"$u / (v_c^2 / (2(\gamma - 1)) $")
plt.title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
plt.ylim((0,2))
figure()
plot(radial_bin_mids,binned_u,'ko',label = "Numerical solution")
legend(loc = "lower right")
xlabel(r"$r / r_{vir}$")
ylabel(r"$u / (v_c^2 / (2(\gamma - 1)) $")
title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
ylim((0,2))
plot_filename = "./plots/internal_energy/internal_energy_profile_%03d.png" %i
plt.savefig(plot_filename,format = "png")
plt.close()
savefig(plot_filename,format = "png")
close()
......
#!/bin/bash
# Generate the initial conditions if they are not present.
echo "Generating initial conditions for the isothermal potential box example..."
python makeIC.py 100000
if [ ! -e Hydrostatic.hdf5 ]
then
echo "Generating initial conditions for the isothermal potential box example..."
python makeIC.py 100000
fi
# Run for 10 dynamical times
../swift -g -s -t 2 hydrostatic.yml 2>&1 | tee output.log
../swift -g -s -t 1 hydrostatic.yml 2>&1 | tee output.log
echo "Plotting density profiles"
mkdir plots
......
###############################################################################
# This file is part of SWIFT.
# Copyright (c) 2016 Stefan Arridge (stefan.arridge@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 numpy as np
import h5py as h5
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use("Agg")
from pylab import *
import sys
n_snaps = int(sys.argv[1])
......@@ -24,7 +45,7 @@ unit_mass_cgs = float(params.attrs["InternalUnitSystem:UnitMass_in_cgs"])
unit_length_cgs = float(params.attrs["InternalUnitSystem:UnitLength_in_cgs"])
unit_velocity_cgs = float(params.attrs["InternalUnitSystem:UnitVelocity_in_cgs"])
unit_time_cgs = unit_length_cgs / unit_velocity_cgs
v_c = float(params.attrs["SoftenedIsothermalPotential:vrot"])
v_c = float(params.attrs["IsothermalPotential:vrot"])
v_c_cgs = v_c * unit_velocity_cgs
header = f["Header"]
N = header.attrs["NumPart_Total"][0]
......@@ -45,7 +66,8 @@ for i in range(n_snaps):
f = h5.File(filename,'r')
coords_dset = f["PartType0/Coordinates"]
coords = np.array(coords_dset)
#translate coords by centre of box
#translate coords by centre of box
header = f["Header"]
snap_time = header.attrs["Time"]
snap_time_cgs = snap_time * unit_time_cgs
......@@ -73,7 +95,6 @@ for i in range(n_snaps):
internal_energy_array = np.append(internal_energy_array,total_internal_energy)
#put energies in units of v_c^2 and rescale by number of particles
pe = potential_energy_array / (N*v_c**2)
ke = kinetic_energy_array / (N*v_c**2)
ie = internal_energy_array / (N*v_c**2)
......@@ -82,14 +103,15 @@ te = pe + ke + ie
dyn_time_cgs = r_vir_cgs / v_c_cgs
time_array = time_array_cgs / dyn_time_cgs
plt.plot(time_array,ke,label = "Kinetic Energy")
plt.plot(time_array,pe,label = "Potential Energy")
plt.plot(time_array,ie,label = "Internal Energy")
plt.plot(time_array,te,label = "Total Energy")
plt.legend(loc = "lower right")
plt.xlabel(r"$t / t_{dyn}$")
plt.ylabel(r"$E / v_c^2$")
plt.title(r"$%d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(N,v_c))
plt.ylim((-2,2))
plt.savefig("energy_conservation.png",format = 'png')
figure()
plot(time_array,ke,label = "Kinetic Energy")
plot(time_array,pe,label = "Potential Energy")
plot(time_array,ie,label = "Internal Energy")
plot(time_array,te,label = "Total Energy")
legend(loc = "lower right")
xlabel(r"$t / t_{dyn}$")
ylabel(r"$E / v_c^2$")
title(r"$%d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(N,v_c))
ylim((-2,2))
savefig("energy_conservation.png",format = 'png')
###############################################################################
# This file is part of SWIFT.
# Copyright (c) 2016 Stefan Arridge (stefan.arridge@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 numpy as np
import h5py as h5
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use("Agg")
from pylab import *
import sys
def do_binning(x,y,x_bin_edges):
......@@ -62,7 +83,8 @@ for i in range(n_snaps):
f = h5.File(filename,'r')
coords_dset = f["PartType0/Coordinates"]
coords = np.array(coords_dset)
#translate coords by centre of box
#translate coords by centre of box
header = f["Header"]
snap_time = header.attrs["Time"]
snap_time_cgs = snap_time * unit_time_cgs
......@@ -73,16 +95,15 @@ for i in range(n_snaps):
radius_cgs = radius*unit_length_cgs
radius_over_virial_radius = radius_cgs / r_vir_cgs
#get the internal energies
#get the internal energies
vel_dset = f["PartType0/Velocities"]
vel = np.array(vel_dset)
#make dimensionless
#make dimensionless
vel /= v_c
r = radius_over_virial_radius
#find radial component of velocity
v_r = np.zeros(r.size)
for j in range(r.size):
v_r[j] = -np.dot(coords[j,:],vel[j,:])/radius[j]
......@@ -94,18 +115,13 @@ for i in range(n_snaps):
radial_bin_mids = np.linspace(bin_widths / 2. , max_r - bin_widths / 2. , n_radial_bins)
binned_v_r = v_r_totals / hist
#calculate cooling radius
#r_cool_over_r_vir = np.sqrt((2.*(gamma - 1.)*lambda_cgs*M_vir_cgs*X_H**2)/(4.*np.pi*CONST_m_H_CGS**2*v_c_cgs**2*r_vir_cgs**3))*np.sqrt(snap_time_cgs)
plt.plot(radial_bin_mids,binned_v_r,'ko',label = "Average radial velocity in shell")
#plt.plot((0,1),(1,1),label = "Analytic Solution")
#plt.plot((r_cool_over_r_vir,r_cool_over_r_vir),(0,2),'r',label = "Cooling radius")
plt.legend(loc = "upper right")
plt.xlabel(r"$r / r_{vir}$")
plt.ylabel(r"$v_r / v_c$")
plt.title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
plt.ylim((0,2))
figure()
plot(radial_bin_mids,binned_v_r,'ko',label = "Average radial velocity in shell")
legend(loc = "upper right")
xlabel(r"$r / r_{vir}$")
ylabel(r"$v_r / v_c$")
title(r"$\mathrm{Time}= %.3g \, s \, , \, %d \, \, \mathrm{particles} \,,\, v_c = %.1f \, \mathrm{km / s}$" %(snap_time_cgs,N,v_c))
ylim((-1,1))
plot_filename = "./plots/radial_velocity_profile/velocity_profile_%03d.png" %i
plt.savefig(plot_filename,format = "png")
plt.close()
savefig(plot_filename,format = "png")
close()
......@@ -26,7 +26,7 @@ from numpy import *
# with a density of 1
# Parameters
periodic= 0 # 1 For periodic box
periodic= 1 # 1 For periodic box
boxSize = 1.
rho = 1.
L = int(sys.argv[1]) # Number of particles along one axis
......
......@@ -35,4 +35,4 @@ Statistics:
# Parameters related to the initial conditions
InitialConditions:
file_name: ./uniformDMBox_100.hdf5 # The file to read
file_name: ./uniformDMBox_50.hdf5 # The file to read
......@@ -107,7 +107,7 @@ case $host_cpu in
*2?6[[ad]]?:*:*:*) ax_gcc_arch="sandybridge corei7-avx corei7 core2 pentium-m pentium3 pentiumpro" ;;
*3?6[[ae]]?:*:*:*) ax_gcc_arch="ivybridge core-avx-i corei7-avx corei7 core2 pentium-m pentium3 pentiumpro" ;;
*3?6[[cf]]?:*:*:*|*4?6[[56]]?:*:*:*) ax_gcc_arch="haswell core-avx2 core-avx-i corei7-avx corei7 core2 pentium-m pentium3 pentiumpro" ;;
*3?6d?:*:*:*) ax_gcc_arch="broadwell core-avx2 core-avx-i corei7-avx corei7 core2 pentium-m pentium3 pentiumpro" ;;
*3?6d?:*:*:*|*4?6f?:*:*:*) ax_gcc_arch="broadwell core-avx2 core-avx-i corei7-avx corei7 core2 pentium-m pentium3 pentiumpro" ;;
*1?6c?:*:*:*|*2?6[[67]]?:*:*:*|*3?6[[56]]?:*:*:*) ax_gcc_arch="bonnell atom core2 pentium-m pentium3 pentiumpro" ;;
*3?67?:*:*:*|*[[45]]?6[[ad]]?:*:*:*) ax_gcc_arch="silvermont atom core2 pentium-m pentium3 pentiumpro" ;;
*000?f[[012]]?:*:*:*|?f[[012]]?:*:*:*|f[[012]]?:*:*:*) ax_gcc_arch="pentium4 pentiumpro" ;;
......
......@@ -63,7 +63,7 @@ AM_SOURCES = space.c runner.c queue.c task.c cell.c engine.c \
nobase_noinst_HEADERS = align.h approx_math.h atomic.h cycle.h error.h inline.h kernel_hydro.h kernel_gravity.h \
kernel_long_gravity.h vector.h cache.h runner_doiact.h runner_doiact_vec.h runner_doiact_grav.h runner_doiact_fft.h \
runner_doiact_nosort.h units.h intrinsics.h minmax.h kick.h timestep.h drift.h adiabatic_index.h io_properties.h \
dimension.h equation_of_state.h part_type.h \
dimension.h equation_of_state.h part_type.h periodic.h \
gravity.h gravity_io.h \
gravity/Default/gravity.h gravity/Default/gravity_iact.h gravity/Default/gravity_io.h \
gravity/Default/gravity_debug.h gravity/Default/gravity_part.h \
......
......@@ -29,25 +29,48 @@
#include "timeline.h"
/**
* @brief Check that a cell been drifted to the current time.
* @brief Check that the #part in a #cell have been drifted to the current time.
*
* @param c The #cell.
* @param e The #engine containing information about the current time.
* @return 1 if the #cell has been drifted to the current time, 0 otherwise.
*/
__attribute__((always_inline)) INLINE static int cell_is_drifted(
__attribute__((always_inline)) INLINE static int cell_are_part_drifted(
const struct cell *c, const struct engine *e) {
#ifdef SWIFT_DEBUG_CHECKS
if (c->ti_old > e->ti_current)
if (c->ti_old_part > e->ti_current)
error(
"Cell has been drifted too far forward in time! c->ti_old=%lld (t=%e) "
"and e->ti_current=%lld (t=%e)",
c->ti_old, c->ti_old * e->timeBase, e->ti_current,
c->ti_old_part, c->ti_old_part * e->timeBase, e->ti_current,
e->ti_current * e->timeBase);
#endif
return (c->ti_old == e->ti_current);
return (c->ti_old_part == e->ti_current);
}
/**
* @brief Check that the #gpart in a #cell have been drifted to the current
* time.
*
* @param c The #cell.
* @param e The #engine containing information about the current time.
* @return 1 if the #cell has been drifted to the current time, 0 otherwise.
*/
__attribute__((always_inline)) INLINE static int cell_are_gpart_drifted(
const struct cell *c, const struct engine *e) {
#ifdef SWIFT_DEBUG_CHECKS
if (c->ti_old_gpart > e->ti_current)
error(
"Cell has been drifted too far forward in time! c->ti_old=%lld (t=%e) "
"and e->ti_current=%lld (t=%e)",
c->ti_old_gpart, c->ti_old_gpart * e->timeBase, e->ti_current,
e->ti_current * e->timeBase);
#endif
return (c->ti_old_gpart == e->ti_current);
}
/* Are cells / particles active for regular tasks ? */
......
......@@ -99,7 +99,8 @@ int cell_unpack(struct pcell *pc, struct cell *c, struct space *s) {
c->h_max = pc->h_max;
c->ti_end_min = pc->ti_end_min;
c->ti_end_max = pc->ti_end_max;
c->ti_old = pc->ti_old;
c->ti_old_part = pc->ti_old_part;
c->ti_old_gpart = pc->ti_old_gpart;
c->count = pc->count;
c->gcount = pc->gcount;
c->scount = pc->scount;
......@@ -128,7 +129,8 @@ int cell_unpack(struct pcell *pc, struct cell *c, struct space *s) {
if (k & 1) temp->loc[2] += temp->width[2];
temp->depth = c->depth + 1;
temp->split = 0;
temp->dx_max = 0.f;
temp->dx_max_part = 0.f;
temp->dx_max_gpart = 0.f;
temp->dx_max_sort = 0.f;
temp->nodeID = c->nodeID;
temp->parent = c;
......@@ -239,7 +241,8 @@ int cell_pack(struct cell *c, struct pcell *pc) {
pc->h_max = c->h_max;
pc->ti_end_min = c->ti_end_min;
pc->ti_end_max = c->ti_end_max;
pc->ti_old = c->ti_old;
pc->ti_old_part = c->ti_old_part;
pc->ti_old_gpart = c->ti_old_gpart;
pc->count = c->count;
pc->gcount = c->gcount;
pc->scount = c->scount;
......@@ -1018,7 +1021,7 @@ void cell_clean_links(struct cell *c, void *data) {
}
/**
* @brief Checks that the particles in a cell are at the
* @brief Checks that the #part in a cell are at the
* current point in time
*
* Calls error() if the cell is not at the current time.
......@@ -1026,7 +1029,7 @@ void cell_clean_links(struct cell *c, void *data) {
* @param c Cell to act upon
* @param data The current time on the integer time-line
*/
void cell_check_particle_drift_point(struct cell *c, void *data) {
void cell_check_part_drift_point(struct cell *c, void *data) {
#ifdef SWIFT_DEBUG_CHECKS
......@@ -1035,14 +1038,40 @@ void cell_check_particle_drift_point(struct cell *c, void *data) {
/* Only check local cells */
if (c->nodeID != engine_rank) return;
if (c->ti_old != ti_drift)
error("Cell in an incorrect time-zone! c->ti_old=%lld ti_drift=%lld",
c->ti_old, ti_drift);
if (c->ti_old_part != ti_drift)
error("Cell in an incorrect time-zone! c->ti_old_part=%lld ti_drift=%lld",
c->ti_old_part, ti_drift);
for (int i = 0; i < c->count; ++i)
if (c->parts[i].ti_drift != ti_drift)
error("part in an incorrect time-zone! p->ti_drift=%lld ti_drift=%lld",
c->parts[i].ti_drift, ti_drift);
#else
error("Calling debugging code without debugging flag activated.");
#endif
}
/**
* @brief Checks that the #gpart and #spart in a cell are at the
* current point in time
*
* Calls error() if the cell is not at the current time.
*
* @param c Cell to act upon
* @param data The current time on the integer time-line
*/
void cell_check_gpart_drift_point(struct cell *c, void *data) {
#ifdef SWIFT_DEBUG_CHECKS
const integertime_t ti_drift = *(integertime_t *)data;
/* Only check local cells */
if (c->nodeID != engine_rank) return;
if (c->ti_old_gpart != ti_drift)
error("Cell in an incorrect time-zone! c->ti_old_gpart=%lld ti_drift=%lld",
c->ti_old_gpart, ti_drift);
for (int i = 0; i < c->gcount; ++i)
if (c->gparts[i].ti_drift != ti_drift)
......@@ -1327,7 +1356,7 @@ int cell_unskip_tasks(struct cell *c, struct scheduler *s) {
error("bad flags in sort task.");
#endif