diff --git a/solvers/DPStokes/mobility.h b/solvers/DPStokes/mobility.h index 3384b23..7da91fc 100644 --- a/solvers/DPStokes/mobility.h +++ b/solvers/DPStokes/mobility.h @@ -148,14 +148,17 @@ class DPStokes : public libmobility::Mobility { this->dppar.zmax += 1.5 * this->dppar.w * h / 2; } real Lz = this->dppar.zmax - this->dppar.zmin; + if (Lz <= 2 * this->dppar.hydrodynamicRadius) { + throw std::runtime_error("[DPStokes] The box size in z is too small to " + "fit the particles. Try increasing zmax."); + } real H = Lz / 2; // sets chebyshev node spacing at its coarsest (in the middle) to be h real nz_actual = M_PI / (asin(h / H)) + 1; - // pick nearby N such that 2(Nz-1) has two factors of 2 and is FFT - // friendly - this->dppar.nz = floor(nz_actual); - this->dppar.nz += (int)ceil(nz_actual) % 2; + // pick nearby N such that 2(Nz-1) has two factors of 2 and is FFT friendly + this->dppar.nz = (int)floor(nz_actual); + this->dppar.nz += (this->dppar.nz - 1) % 2; dpstokes->initialize(dppar); Mobility::initialize(ipar); diff --git a/tests/test_initialization.py b/tests/test_initialization.py index a001ed6..93a7b12 100644 --- a/tests/test_initialization.py +++ b/tests/test_initialization.py @@ -97,3 +97,10 @@ def test_nbody_wall_height_parameter(): ) # single_wall periodicity requires wall height param solver.setParameters(wallHeight=0.5) + + +def test_dpstokes_narrow_channel(): + solver = DPStokes("periodic", "periodic", "two_walls") + solver.setParameters(Lx=3.0, Ly=3.0, zmin=0.0, zmax=1.0) + with pytest.raises(RuntimeError, match="DPStokes"): + solver.initialize(hydrodynamicRadius=1.0, viscosity=1.0)