diff --git a/INSTALL.swift b/INSTALL.swift index 8e1635b0715426512503fd9dcde32f59a7ad1b62..90096bc1b88d34ab86c04f623f1d3f04ca5a2997 100644 --- a/INSTALL.swift +++ b/INSTALL.swift @@ -143,6 +143,8 @@ before you can build it. - DOXYGEN: the doxygen library is required to create the SWIFT API documentation. + - python: Examples and solution script use python and rely on the + numpy library version 1.8.2 or higher. SWIFT Coding style diff --git a/examples/BigCosmoVolume/makeIC.py b/examples/BigCosmoVolume/makeIC.py index c141337c06fb28aa4049e2823fcc7cd3e9d5513c..8f83b564a3f747d164fd03b7dddf3cd9c609d9eb 100644 --- a/examples/BigCosmoVolume/makeIC.py +++ b/examples/BigCosmoVolume/makeIC.py @@ -107,7 +107,7 @@ if n_copy > 1: for i in range(n_copy): for j in range(n_copy): for k in range(n_copy): - coords = np.append(coords, coords_tile + np.array([ i * boxSize, j * boxSize, k * boxSize ]), axis=0) + coords = np.append(coords, coords_tile + np.array([ i * boxSize[0], j * boxSize[1], k * boxSize[2] ]), axis=0) v = np.append(v, v_tile, axis=0) m = np.append(m, m_tile) h = np.append(h, h_tile) diff --git a/examples/CoolingHalo/makeIC.py b/examples/CoolingHalo/makeIC.py index 0b542e200da709e2cc7f668ab8b62b94e0bf95ee..3ec1be6f7b5e568ebe8e0fefe508ef8287edb29c 100644 --- a/examples/CoolingHalo/makeIC.py +++ b/examples/CoolingHalo/makeIC.py @@ -227,7 +227,7 @@ ds[()] = u u = np.zeros(1) # Particle IDs -ids = 1 + np.linspace(0, N, N, endpoint=False, dtype='L') +ids = 1 + np.linspace(0, N, N, endpoint=False) ds = grp.create_dataset('ParticleIDs', (N, ), 'L') ds[()] = ids diff --git a/examples/CoolingHaloWithSpin/makeIC.py b/examples/CoolingHaloWithSpin/makeIC.py index a6d57868ad7542498b27007a5c3ef9234b9feb84..2cf3127c743f61756b3ff6c4a7738c83d185f9cd 100644 --- a/examples/CoolingHaloWithSpin/makeIC.py +++ b/examples/CoolingHaloWithSpin/makeIC.py @@ -233,7 +233,7 @@ ds[()] = u u = np.zeros(1) # Particle IDs -ids = 1 + np.linspace(0, N, N, endpoint=False, dtype='L') +ids = 1 + np.linspace(0, N, N, endpoint=False) ds = grp.create_dataset('ParticleIDs', (N, ), 'L') ds[()] = ids diff --git a/examples/DiscPatch/GravityOnly/makeIC.py b/examples/DiscPatch/GravityOnly/makeIC.py index 42cd26e235deb17a899a65050ef5caa9c832c59c..6e4e148392eb7ca2fbf8c29c3f737d029916c59b 100644 --- a/examples/DiscPatch/GravityOnly/makeIC.py +++ b/examples/DiscPatch/GravityOnly/makeIC.py @@ -150,7 +150,7 @@ ds[()] = m m = numpy.zeros(1) -ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False, dtype='L') +ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False) ds = grp1.create_dataset('ParticleIDs', (numPart, ), 'L') ds[()] = ids diff --git a/examples/DiscPatch/HydroStatic/makeIC.py b/examples/DiscPatch/HydroStatic/makeIC.py index e2846d82a8cfa8bf08de83632b19ae2e7818f3c1..6ba1ccd06fed84ca728aadaa5922dbba536b6881 100644 --- a/examples/DiscPatch/HydroStatic/makeIC.py +++ b/examples/DiscPatch/HydroStatic/makeIC.py @@ -205,7 +205,7 @@ if (entropy_flag == 1): else: ds[()] = u -ids = 1 + numpy.linspace(0, numGas, numGas, endpoint=False, dtype='L') +ids = 1 + numpy.linspace(0, numGas, numGas, endpoint=False) ds = grp0.create_dataset('ParticleIDs', (numGas, ), 'L') ds[()] = ids diff --git a/examples/HydrostaticHalo/density_profile.py b/examples/HydrostaticHalo/density_profile.py index 52bebb9ffefa77dae66af155fb31fed539dcde13..d0afd399f951cf3b727e869ca8571a3a802c2e8d 100644 --- a/examples/HydrostaticHalo/density_profile.py +++ b/examples/HydrostaticHalo/density_profile.py @@ -28,7 +28,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 #lambda_cgs = float(params.attrs["LambdaCooling:lambda_cgs"]) #X_H = float(params.attrs["LambdaCooling:hydrogen_mass_abundance"]) diff --git a/examples/HydrostaticHalo/internal_energy_profile.py b/examples/HydrostaticHalo/internal_energy_profile.py index a1b2bda314a66eb965974d34519f66c544ee8aed..ea52cf8fc5fd098a46f05eaa58494529a868000c 100644 --- a/examples/HydrostaticHalo/internal_energy_profile.py +++ b/examples/HydrostaticHalo/internal_energy_profile.py @@ -46,7 +46,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 #lambda_cgs = float(params.attrs["LambdaCooling:lambda_cgs"]) #X_H = float(params.attrs["LambdaCooling:hydrogen_mass_abundance"]) diff --git a/examples/HydrostaticHalo/makeIC.py b/examples/HydrostaticHalo/makeIC.py index f33387e18dd0ab523684227f5c745b5c8b807b7f..d5081ac84473edc87857c6872278b4d0ca6389b1 100644 --- a/examples/HydrostaticHalo/makeIC.py +++ b/examples/HydrostaticHalo/makeIC.py @@ -227,7 +227,7 @@ ds[()] = u u = np.zeros(1) # Particle IDs -ids = 1 + np.linspace(0, N, N, endpoint=False, dtype='L') +ids = 1 + np.linspace(0, N, N, endpoint=False) ds = grp.create_dataset('ParticleIDs', (N, ), 'L') ds[()] = ids diff --git a/examples/HydrostaticHalo/velocity_profile.py b/examples/HydrostaticHalo/velocity_profile.py index f6a7350b9731d660b2092266d4d6ad3730bab48c..9133195d942233514148aa419003ee0ab7923494 100644 --- a/examples/HydrostaticHalo/velocity_profile.py +++ b/examples/HydrostaticHalo/velocity_profile.py @@ -46,7 +46,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] diff --git a/examples/IsothermalPotential/makeIC.py b/examples/IsothermalPotential/makeIC.py index 7d1c5361f9a255365517226e49c55a8a50c4d6ce..27ddf15fe63d69d4172b927cfca8b8662d11ea4e 100644 --- a/examples/IsothermalPotential/makeIC.py +++ b/examples/IsothermalPotential/makeIC.py @@ -138,7 +138,7 @@ ds = grp1.create_dataset('Masses', (numPart,), 'f') ds[()] = m m = numpy.zeros(1) -ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False, dtype='L') +ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False) ds = grp1.create_dataset('ParticleIDs', (numPart, ), 'L') ds[()] = ids diff --git a/examples/Noh_3D/makeIC.py b/examples/Noh_3D/makeIC.py index ec8d46639eecdefd55b917a955f13efc8ffe4126..0c25a5c8b3e967185cf16bae4b1f21c215266def 100644 --- a/examples/Noh_3D/makeIC.py +++ b/examples/Noh_3D/makeIC.py @@ -35,8 +35,8 @@ glass = h5py.File("glassCube_64.hdf5", "r") vol = 8. -pos = glass["/PartType0/Coordinates"][:,:] * cbrt(vol) -h = glass["/PartType0/SmoothingLength"][:] * cbrt(vol) +pos = glass["/PartType0/Coordinates"][:,:] * vol**(1./3.) +h = glass["/PartType0/SmoothingLength"][:] * vol**(1./3.) numPart = size(h) # Generate extra arrays @@ -65,7 +65,7 @@ file = h5py.File(fileName, 'w') # Header grp = file.create_group("/Header") -grp.attrs["BoxSize"] = [cbrt(vol), cbrt(vol), cbrt(vol)] +grp.attrs["BoxSize"] = [vol**(1./3.), vol**(1./3.), vol**(1./3.)] grp.attrs["NumPart_Total"] = [numPart, 0, 0, 0, 0, 0] grp.attrs["NumPart_Total_HighWord"] = [0, 0, 0, 0, 0, 0] grp.attrs["NumPart_ThisFile"] = [numPart, 0, 0, 0, 0, 0] diff --git a/examples/SineWavePotential_3D/plotSolution.py b/examples/SineWavePotential_3D/plotSolution.py index 7ae5dcd2a52235812c3c11cfcf4d01f4fde647af..13cae037b64eff4ad4fec0003bf0f5ad3ce94896 100644 --- a/examples/SineWavePotential_3D/plotSolution.py +++ b/examples/SineWavePotential_3D/plotSolution.py @@ -55,7 +55,7 @@ rho_x = 1000.*np.exp(-0.5*A/np.pi/cs2*np.cos(2.*np.pi*x)) P = cs2*rho -n1D = int(np.cbrt(len(P))) +n1D = np.ceil(len(P)**(1./3.)) gradP = np.zeros(P.shape) for i in range(len(P)): iself = int(ids[i]/n1D/n1D) diff --git a/src/space.c b/src/space.c index b34571368bc7b5fe6644f3af931dfab3812bbf90..fec47f28b26e0bf22265cb765525b913acd7a638 100644 --- a/src/space.c +++ b/src/space.c @@ -2132,47 +2132,52 @@ void space_split_recursive(struct space *s, struct cell *c, c->split = 0; maxdepth = c->depth; - /* Get dt_min/dt_max. */ + timebin_t time_bin_min = num_time_bins, time_bin_max = 0; + + /* parts: Get dt_min/dt_max and h_max. */ + for (int k = 0; k < count; k++) { +#ifdef SWIFT_DEBUG_CHECKS + if (parts[k].time_bin == time_bin_inhibited) + error("Inhibited particle present in space_split()"); +#endif + time_bin_min = min(time_bin_min, parts[k].time_bin); + time_bin_max = max(time_bin_max, parts[k].time_bin); + h_max = max(h_max, parts[k].h); + } + /* parts: Reset x_diff */ for (int k = 0; k < count; k++) { - struct part *p = &parts[k]; - struct xpart *xp = &xparts[k]; - const float h = p->h; - const integertime_t ti_end = - get_integer_time_end(e->ti_current, p->time_bin); - const integertime_t ti_beg = - get_integer_time_begin(e->ti_current + 1, p->time_bin); - xp->x_diff[0] = 0.f; - xp->x_diff[1] = 0.f; - xp->x_diff[2] = 0.f; - if (h > h_max) h_max = h; - if (ti_end < ti_end_min) ti_end_min = ti_end; - if (ti_end > ti_end_max) ti_end_max = ti_end; - if (ti_beg > ti_beg_max) ti_beg_max = ti_beg; + xparts[k].x_diff[0] = 0.f; + xparts[k].x_diff[1] = 0.f; + xparts[k].x_diff[2] = 0.f; } + /* gparts: Get dt_min/dt_max, reset x_diff. */ for (int k = 0; k < gcount; k++) { - struct gpart *gp = &gparts[k]; - const integertime_t ti_end = - get_integer_time_end(e->ti_current, gp->time_bin); - const integertime_t ti_beg = - get_integer_time_begin(e->ti_current + 1, gp->time_bin); - gp->x_diff[0] = 0.f; - gp->x_diff[1] = 0.f; - gp->x_diff[2] = 0.f; - if (ti_end < ti_end_min) ti_end_min = ti_end; - if (ti_end > ti_end_max) ti_end_max = ti_end; - if (ti_beg > ti_beg_max) ti_beg_max = ti_beg; +#ifdef SWIFT_DEBUG_CHECKS + if (gparts[k].time_bin == time_bin_inhibited) + error("Inhibited s-particle present in space_split()"); +#endif + time_bin_min = min(time_bin_min, gparts[k].time_bin); + time_bin_max = max(time_bin_max, gparts[k].time_bin); + + gparts[k].x_diff[0] = 0.f; + gparts[k].x_diff[1] = 0.f; + gparts[k].x_diff[2] = 0.f; } + /* sparts: Get dt_min/dt_max */ for (int k = 0; k < scount; k++) { - struct spart *sp = &sparts[k]; - const integertime_t ti_end = - get_integer_time_end(e->ti_current, sp->time_bin); - const integertime_t ti_beg = - get_integer_time_begin(e->ti_current + 1, sp->time_bin); - if (ti_end < ti_end_min) ti_end_min = ti_end; - if (ti_end > ti_end_max) ti_end_max = ti_end; - if (ti_beg > ti_beg_max) ti_beg_max = ti_beg; +#ifdef SWIFT_DEBUG_CHECKS + if (sparts[k].time_bin == time_bin_inhibited) + error("Inhibited g-particle present in space_split()"); +#endif + time_bin_min = min(time_bin_min, sparts[k].time_bin); + time_bin_max = max(time_bin_max, sparts[k].time_bin); } + /* Convert into integer times */ + ti_end_min = get_integer_time_end(e->ti_current, time_bin_min); + ti_end_max = get_integer_time_end(e->ti_current, time_bin_max); + ti_beg_max = get_integer_time_begin(e->ti_current + 1, time_bin_max); + /* Construct the multipole and the centre of mass*/ if (s->gravity) gravity_P2M(c->multipole, c->gparts, c->gcount); } diff --git a/src/task.c b/src/task.c index e3c217ebe06f2aba47f4ceaf51af40b2e4557e24..0ac7bc0c59eb678383adba9587a8026cf872ee6d 100644 --- a/src/task.c +++ b/src/task.c @@ -283,6 +283,10 @@ void task_unlock(struct task *t) { /* Act based on task type. */ switch (type) { + case task_type_init: + case task_type_kick1: + case task_type_kick2: + case task_type_timestep: case task_type_drift: cell_unlocktree(ci); cell_gunlocktree(ci); @@ -366,6 +370,10 @@ int task_lock(struct task *t) { #endif break; + case task_type_init: + case task_type_kick1: + case task_type_kick2: + case task_type_timestep: case task_type_drift: if (ci->hold || ci->ghold) return 0; if (cell_locktree(ci) != 0) return 0; diff --git a/tests/test125cells.c b/tests/test125cells.c index aedf0e0b77898cadf97383dfc956a3d6b5ebf39f..b3895ffa9fc0e4c147bfbf58dc1b5a6301b02763 100644 --- a/tests/test125cells.c +++ b/tests/test125cells.c @@ -532,6 +532,7 @@ int main(int argc, char *argv[]) { engine.s = &space; engine.time = 0.1f; engine.ti_current = 8; + engine.max_active_bin = num_time_bins; struct runner runner; runner.e = &engine; diff --git a/tests/test27cells.c b/tests/test27cells.c index 929a148d1f5730b63de79e9a1ab7e25f1ca7311e..599c8713a4b2077444b10f85ed37ee76f5219c5a 100644 --- a/tests/test27cells.c +++ b/tests/test27cells.c @@ -396,6 +396,7 @@ int main(int argc, char *argv[]) { engine.s = &space; engine.time = 0.1f; engine.ti_current = 8; + engine.max_active_bin = num_time_bins; struct runner runner; runner.e = &engine; diff --git a/tests/testPair.c b/tests/testPair.c index 8b23cc419a661f4d50ea53948302729784a129f9..c734424bca58b7a8ce05bba2c3392ca5c600fa9e 100644 --- a/tests/testPair.c +++ b/tests/testPair.c @@ -252,6 +252,7 @@ int main(int argc, char *argv[]) { engine.s = &space; engine.time = 0.1f; engine.ti_current = 8; + engine.max_active_bin = num_time_bins; runner.e = &engine; volume = particles * particles * particles;