Separation of (sub)micro sized particles is of paramount importance due to its vast applications including in particular medical diagnostics [1] . Anomalies in a bioparticle size often indicate various illnesses. This is apparent, for instance, for Alzheimer's, Parkinson's [2] and Huntington's disease [3] , needless to say that often cancer cells noticeably deviates in size from healthy ones [4] . A reliable and effective approach to separating bioparticles is therefore much in demand. They span a large size range covering several orders of magnitude from nanometers to micrometers [5, 6] what crucially complicates the development of such highly anticipated strategies. For instance, deadly viruses like HIV or COVID-19 are approximately 0.1 micrometer in diameter [7] . On the other hand, the soma of a neuron can vary from 4 to 100 micrometers in diameter [8] . For such a broad (sub)micrometer scale, efficient isolation techniques are required to allow tunability of the particle size intended for separation. Unluckily, the latter are rather scarce [9, 10] , however, recently some progress has been made in this direction [11] [12] [13] In this work we aspire to partially fill this significant know-how gap by demonstrating a nonintuitive, yet efficient separation strategy taking advantage of a paradoxical mechanism of negative mobility [14] [15] [16] [17] . We show that under the action of a static bias only particles of a given linear size move in the direction opposite to this net force whereas the others migrate concurrently towards it. This effect creates a possibility of steering different species of particles in opposite directions under identical experimental conditions thus facilitating their separation. A proof-of-principle experiment of a similar isolation scheme was performed with insulator dielectrophoresis in a nonlinear, symmetric microfluidic structure with electrokinetically activated transport [18, 19] . Recently, such a setup allowed to induce the negative mobility not only for a colloidal particle but also for a biological compound in the form of mouse-liver mitochondria [20] .
Motivated by the large size range encountered in biochemical applications, a crucial result of this work is a demonstration of a tunable separation strategy, in which the size of the particle undergoing separation may be effectively controlled by variation of the parameters characterizing the external force applied to the particle, e.g. the magnitude of the static bias. The same setup can be applied to segregate particles with respect to their mass in a similar tunable manner [21] . This approach may provide selectivity required for individual isolation of nano and micro particles, proteins, organelles and cells and thus constitutes an important step towards the development of robust lab-on-a-chip devices exploited in both research and industrial applications, in particular point-of-care medical diagnostics.
The paper is organized as follows. In Sec. II we outline the model of a Brownian particle dwelling in a spatially periodic potential under the action of both an external harmonic driving as well as a constant bias. In the next section we exemplify the negative mobility phenomenon. Sec. IV provides crucial results of the paper, namely, a tunable particle separation strategy. In Sec. V we discuss the possibility of tailoring particle isolation in the considered system. Finally, the last section is devoted to summary and conclusion.
The system considered in this study is a classical inertial Brownian particle of mass M which moves in a spatially periodic one-dimensional potential U (x) = U (x + L) of the period L, additionally subjected to an unbiased timeperiodic force A cos (Ωt) of the amplitude A and the angular frequency Ω, as well as an external static force F . Dynamics of such a particle is described by the Langevin equation [22] Mẍ
where the dot and the prime denote differentiation with respect to time t and the particle coordinate x, respectively. Coupling of the particle with thermal bath of temperature T is modeled by Gaussian white noise of zero mean and unity intensity, namely
The noise intensity factor 2Γk B T (where k B is the Boltzmann constant) follows from the fluctuation-dissipation theorem [23] and ensures that the system reaches the equilibrium state when A = 0 and F = 0. The potential U (x) is spatially periodic with the period L and the barrier height 2∆U ,
At first glance the studied system looks simply, however, it exhibits peculiar transport behaviour including noiseenhanced transport efficiency [16] , anomalous diffusion [24, 25] , amplification of normal diffusion [26, 27] and a non-monotonic temperature dependence of normal diffusion [28] .
As the first step of the analysis we transform the Eq. (1) into its dimensionless form. This aim can be achieved in several ways. It often allows to simplify the setup description as after the rescaling procedure some parameters appearing in the corresponding dimensional version may be eliminated thus reducing the complexity of the problem. Moreover, recasting into the dimensionless variables ensures that the results obtained later are independent of the setup which is essential to facilitate the choice in realizing the best setup for testing theoretical predictions in experiments. Here, we propose the use of the following scales as the characteristic units of length and timê
The Langevin equation (1) transformed according to the above variables readŝ
The rescaled dimensionless friction coefficient γ is the ratio of two characteristic time scales
where τ 1 = M/Γ is characteristic time for the velocity relaxation of the free Brownian particle. The dimensionless mass is set to unity, m = 1. Other rescaled parameters are as follows: a = (L/∆U )A, ω = τ 0 Ω, f = (L/∆U )F . The dimensionless potentialÛ (x) = U (Lx)/∆U = sin(2πx) has the periodL = 1. The dimensionless thermal noise ζ(t) assumes the same statistical properties as ξ(t), i.e. ζ(t) = 0 and ζ(t)ζ(ŝ) = δ(t −ŝ). The dimensionless noise intensity D = k B T /∆U is the ratio of thermal energy and half of the non-rescaled potential barrier. From now on, only the dimensionless variables will be used in this study and therefore, in order to simplify the notation, the ∧-symbol will be omitted in the equation (5) . The dimensionless friction coefficient γ is the most important parameter for later mentioned process of the particle separation with respect to its size. It is due to the fact that even for the simplest model of hydrodynamic interactions occurring in this setup it depends on the linear size R of the particle. For instance, the spherical particle diffusing in the surrounding medium is subjected to Stokes drag −Γẋ where Γ = 6πηR and η is the viscosity of the environment [29] . We stress that (sub)micro sized particles typically possess small mass and therefore for them the dimensionless friction coefficient γ given by Eq. (6) is expected to be either of the order or larger than the dimensionless mass m = 1. A. The directed velocity
The particle mobility describes its ability to move through the medium as a response to the biased force acting on it. Hence, the observable of foremost interest in this study is a directed velocity v of the particle which may be written as
where · indicates averaging over all realizations of thermal noise as well as over initial conditions for the particle position x(0) and its velocityẋ(0). The latter is obligatory for the deterministic limit D ∝ T → 0 when the dynamics may be non-ergodic and results can be affected by the specific choice of initial conditions [30] . The Fokker-Planck equation corresponding to the Langevin Eq. (5) cannot be solved analytically in a closed form. Therefore the system may be analyzed only by numerical simulations. Equation (5) is characterized by a 5-dimensional parameter space {γ, a, ω, f, D} the detailed exploration of which is a very challenging task. All numerical calculations were performed using an innovative computational method which is based on employing graphics processing unit supercomputers. This procedure allowed us to speed up computations by about 10 3 times compared to traditional methods. For technical details we refer the reader to Ref. [31] .
In general, the directed velocity v is an increasing function of the static force f and the resultant particle transport follows the direction of the bias, i.e. v = µ(f )f with a positive nonlinear mobility µ(f ) > 0. However, in the parameter space there are also regimes for which the particle moves on average in the opposite way, namely v < 0 for f > 0. Such anomalous transport behaviour is characterized by a negative mobility µ(f ) < 0 [15, 32, 33] . The key ingredient for the occurrence of the latter effect is that the system is driven away from thermal equilibrium into a time-dependent nonequilibrium state [15] . This fact provides a negation of the Le Chatelier-Braun equilibrium principle [34] stating that the response of the perturbed system occurs into the direction of the applied bias towards a new equilibrium state. In our case this condition is guaranteed by the presence of the external harmonic driving a cos (ωt).
It is known that there exist two fundamentally different mechanisms responsible for the emergence of negative mobility in the considered system. Firstly, it may be generated solely by the deterministic dynamics given by Eq. (5) with D ∝ T = 0 [15, 32] . Secondly, it can be induced by an appropriate dose of thermal fluctuations [15] . Among the first mentioned class we distinguish two completely distinct scenarios. For the deterministic counterpart of the system the negative mobility may be induced either by chaos-assisted dynamics [32] or regular, non-chaotic attractors transporting the particle in the direction opposite to the applied bias [35] . Our numerical research reveals that the most common reason for the occurrence of the negative mobility is rooted in the complexity of the deterministic and chaotic dynamics [35, 36] . This observation is of great importance for uncovering the parameter regimes allowing the particle separation.
In Fig. 1 we illustrate the negative mobility effect. For γ = 1.29 the directed velocity v > 0 assumes the same sign as the force f > 0 leading to the normal particle transport regime with µ > 0. However, for γ = 1.23, the directed velocity v < 0 is opposite to the bias f > 0 and in consequence the negative mobility effect emerges µ < 0. If the value of f is positive and large enough then the sign of v coincides with the force f again. Therefore this anomalous transport behaviour is observed only in the vicinity of the zero bias f = 0 and therefore often is termed as the absolute negative mobility [37] . As it is illustrated in the the panel (b) the presented parameter regime belongs to the class of deterministically induced negative mobility as for the limiting case of vanishing thermal noise intensity D → 0 the directed velocity is negative v < 0. We want to stress that this limit should be considered with utmost care as in such a case attractors transporting the particle in opposite directions may coexist and the dynamics may be non-ergodic. It means that depending on the initial conditions, a particle would either move in the direction of the bias force or opposite to it [32] . Nevertheless, at any finite D > 0, possibly coexisting deterministic attractors turn metastable and due to thermally activated transitions between them the ergodicity of dynamics is restored. Consequently, the directed velocity v is independent of the initial conditions. Moreover, as it is shown in the panel (b) the deterministically induced negative mobility effect is generally quite robust with temperature change and usually survives up to moderate thermal noise intensities.
We now want to harvest the negative mobility phenomenon to separate (sub)micrometer sized particles. This task may be achieved as follows. Imagine that there is a mixture of several species of spherical particles each differing by its linear size R. The friction coefficient γ, via e.g. Stokes formula, depends on the particle radius γ = γ(R). We can extract the given species of particles R * characterized by the friction coefficient γ * ≡ γ(R * ) if only for this particular group of particles the negative mobility effect arises, i.e. v < 0 for γ * and v > 0 for the rest. Therefore this task translates to discovery in the complex four-dimensional parameter subspace {a, ω, f, D} regimes where in the characteristic v (γ) there exists only one interval δγ of the friction coefficient around the desired particle size γ * for which the negative mobility v < 0 emerges. In such a way only particles with the radius R * would be extracted from the mixture.
Motivated by the large size range typically encountered in biochemical applications, we aim to develop a tunable scheme that allows to control the particle size targeted for isolation by changing only one parameter of the system. In [21] it has been shown that the negative mobility effect can be harvested to separate the (sub)micrometer sized particles with respect to their mass. The particle mass targeted for isolation might be effectively controlled over a regime of nearly two orders of magnitude upon changing solely the frequency ω of the external harmonic driving. Moreover, in [36] an efficient separation mechanism based on thermal fluctuation induced negative mobility phenomenon has been proposed. By tuning solely temperature of the system D ∝ T , one can extract from the mixture of particle species differing by size only those of a strictly defined radius. This scheme opens an opportunity to separate particles that carry no charge or dipole, however, it may be inconvenient since temperature variation typically takes too much time to offer a robust experimental implementation. Therefore, in contrast, in this work we harvest the negative mobility effect to develop the particle separation strategy in which the particle size intended for isolation will be controlled by changing only the parameters characterizing the externally applied force, namely, the static bias f or the amplitude a or the frequency ω of the harmonic driving. Unfortunately, there is no clear relationship between the presence of the negative mobility and the model parameter values. A tiny displacement in the parameter space may either cause a sudden emergence of the negative mobility or its rapid decay. Therefore, extensive numerical simulations of Eq. (5) were performed in order to systematically investigate the established parameter space. As the deterministically induced negative mobility is the most populated mechanism in the parameter space we set D = 0. Then Eq. Overall, we considered nearly 10 9 different parameter sets. This exceptional precision was possible only because of our innovative simulation method [31] .
The so collected data was transformed into two-dimensional maps presenting the directed velocity v versus two chosen model parameters to facilitate the further analysis. The results of foremost interest are those with γ dependence since the friction coefficient may be used as an indicator to differentiate particles by their size. We explored the data to discover any correlations between the presence of the negative mobility, the friction coefficient γ as well as the magnitude of the parameters a, ω and f . We exemplify such a situation in Fig. 2 where we depict the directed velocity v as a function of the static bias f and the friction coefficient γ. The color bar in the plot represents the magnitude of the directed velocity v . The occurrence of the negative mobility is marked by blue areas for which v < 0. The reader can observe a linear trend between the friction coefficient γ, the static bias f and the emergence of the anomalous transport. The negative mobility effect occurs for progressively smaller values of γ as the force f increases. Among many so discovered negative mobility regimes we distilled only those where the latter phenomenon is present solely for one indivisible interval δγ of the friction coefficient γ thus permitting the particle separation.
We exemplify this procedure in Fig. 3 where we depict the directed velocity v versus the friction coefficient γ which can be identified with the particle size. As it is illustrated in the panel (a), among many particle sizes corresponding to the friction coefficient γ ∈ [0. 2, 4] only those with the friction coefficient γ * ≈ 1.05 will move in the direction opposite v < 0 to the applied bias f > 0. Other particles will travel concurrently towards it. As a result, only the particles with γ * ≈ 1.05 will be extracted from the heterogeneous mixture. Panel (b) presents a magnification of the interval δγ where the effect of negative mobility emerges. It can be interpreted as a resolution capacity of this method. In this case it reads δγ ≈ 0.0678. The selectivity of the proposed separation scheme is impressive as δγ/γ * ≈ 0.06 . Moreover, as it is illustrated, typically the velocity v is noticeably peaked in the interval δγ where the negative mobility occurs. The mechanism of the latter anomalous transport effect is rooted in the deterministic dynamics and as such thermal fluctuations generally have destructive impact on it [35, 36] . Therefore the above outlined particle separation strategy is viable in low to moderate temperature regimes in which the size targeted for isolation γ * , i.e. the middle of the interval δγ, coincides well with the the friction coefficient γ min for which the directed velocity attains its minimal value v min ≡ v (γ min ). It means that not only the selectivity of this method is impressive but also the particle separation process is quick.
Such approach allowed us to distill parameter regimes which reveal a specific functional dependence between the particle size tailored for separation γ * and the parameters of the external force applied to the system thus allowing for tunable particle isolation. In Fig. 4 we present the friction coefficient γ * versus the static bias f , the amplitude a and the frequency ω all depicted for different temperature of the system D ∝ T . The data were obtained from the characteristics v (γ) computed for many values of the control parameter, c.f. Fig. 3 . Each dot represents the friction coefficient γ * undergoing the separation process at the fixed f , a or ω, see panels (a)-(d). The bars indicate the friction coefficient interval δγ where the negative mobility emerges. Here we note that the curves on the corresponding panels overlap with the one representing the deterministic solution D = 0. Therefore to visualize the impact of temperature on the separation process in the corresponding insets we schematically show the solutions depicting subsequent D values. The reader may observe that using those tailored parameters read off from Fig. 4 one is able to tune the negative mobility to the particle of a given size γ * by changing solely the static bias f or the amplitude a or the frequency ω. In this way it will be separated from the others possessing positive mobility and thus moving concurrently towards the applied bias f .
From the experimental point of view the most convenient way to manipulate the particle separation is presumably by altering the static bias f . It is because in many realistic setups it is implemented via the constant external field, e.g. in the microfludic experiments in the form of a spatially uniform electric field which induces the particle electrophoresis [10] . In most cases its intensity can be changed relatively easily, as opposed to the frequency ω of the external harmonic driving which often requires complete rebuilding of the experimental setup. We note that the parameter sets reported in Fig. 4 allow for the tunable separation of the particles in the regime of moderate-large friction coefficient which is characteristic for low Reynolds numbers, being typical for (sub)micro sized particles immersed in a solution [38] . For example in the panels (a) and (b) corresponding to the isolation driven by the constant force f , the friction coefficient γ * ∈ [0.6, 2.5], whereas in (d) when the separation is controlled by the frequency ω the size γ * ∈ [0.25, 1.75]. The reader can observe that the friction coefficient γ * is a decreasing function of the force f and the frequency ω (panels (a), (b) and (d)) while for the amplitude a it depicts an increasing dependence (panel (c)). Finally, we note that the dimensionless friction coefficient γ in Eq. (6) is influenced not only by the actual friction Γ but also by the parameters ∆U and L of the potential. Therefore, experimentalists may exploit these characteristics of the periodic substrate to further adapt the particle size targeted for separation.
Last but not least, Fig. 4 reveals the impact of temperature D ∝ T on the tunability of the separation process. As it was stated before, since the negative mobility effect derives from the deterministic dynamics of the system, thermal fluctuations have a destructive influence on it. When temperature increases the regions of negative mobility allowing for the controllable particle separation progressively shrink and eventually vanish completely. Therefore the reported tunability is possible in the widest range of the friction coefficient γ * for low to moderate temperature regimes D ∝ T . However, for instance, the panel (c) illustrates an interesting effect of thermal fluctuations. An increase of temperature leads not only to shrinking of the range of the friction coefficient γ * for which the negative mobility is observed but also to significant decrease of the intervals where this phenomenon occurs, thereby optimizing the width δγ. It means that then the tunability of the method is limited but selectivity of the separation process increases. We remark that the particle isolation upon harvesting the negative mobility phenomenon is also present for different parameter regimes, however, the range of its tunability proves somewhat smaller.
Now we consider another, complementary issue. Let us assume that we deal with particles of a given size γ * which we want to extract from the heterogeneous mixture. We address the following question: for how many different sizes γ * taken from the extended interval γ * ∈ [0.1, 10] it is possible to find a parameter set {a, ω, f, D} for which the negative mobility effect emerges in the small interval [γ * − δγ/2, γ * + δγ/2] around the targeted value γ * , therefore allowing its separation in a unique manner, c.f. Fig. 3 ? We found that in most cases the magnitude δγ/γ * is equal to a few percent but frequently is even smaller. We note that for a single γ * there might be several parameter regimes fulfilling this condition thus facilitating the choice in realizing the best parameter set. We present the answer of this experimentally and practically relevant question in Fig. 5 . The distribution of the size γ * targeted for separation is shown there in the parameter plane of the amplitude a and the frequency ω for different values of the static bias f . The color coded scale displays the friction coefficient γ * value. We observe that small particles can be isolated when the static bias f is likewise small and moderate to large values of the frequency ω. On the other hand, medium and large particles are separated for small frequencies ω. We note that the distribution of the particle size targeted for isolation γ * undergoes a stretching when the static bias f decreases. Moreover, in such a case the range of the particles which might be segregated by harvesting the negative mobility effect is extended as well. Finally, even though the considered panel depicts the deterministic D = 0 dynamics, we found that the distribution of γ * depicted there is quite robust with temperature change and survives up to moderate thermal noise intensity, see Fig. 1 (b) .
In this work we provide an efficient method for tunable separation of (sub)micro sized particles via the negative mobility phenomenon. The approach presented here requires only a spatially periodic nonlinear structure in combination with an unbiased external time-periodic driving. In this scheme the particle size intended for isolation can be effectively controlled by changing solely the parameters of the external force applied to the system, namely, the static bias or the amplitude or the frequency of the external harmonic driving. The approach can be further adapted to the needs by proper fabrication of the nonlinear potential landscape determined by its barrier height and period. It allows the possibility to not deflect the separated particles along the different angles but to steer them in the opposite direction making the isolation process robust.
Our theoretical predictions should be used as a guide towards physical reality indicating the direction for future theoretical and experimental research. In particular, one needs to carefully consider higher dimensional systems as well as geometrical constraints together with hydrodynamic interaction which in real experiments may play essential role. We expect that such research would potentially lead to implementation of the proposed scheme in a lab-on-chip device, as it has been recently demonstrated for a similar system [20] . We envision that current lithographic techniques with advantageous fabrication costs may be used to develop high throughput separation applications concerning in particular biophysical and biochemical problems. Taking into account recent progress in 3D printing technologies allowing its scaling down to the nanometer range the proposed scheme may have even significant commercial potential in future [39] .
We employed a weak 2nd order predictor-corrector method [40] to simulate stochastic dynamics given by Eq. (5) . We integrated it with the time step scaled by the fundamental period T = 2π/ω of the external harmonic driving, namely h = 10 −2 ×T, with the exception of smallest ω < 1 values for which the step was chosen to be h = 10 −3 ×T. The initial positions x(0) and velocities v(0) were uniformly distributed over the intervals [0, 1] and [−2, 2], respectively. The directed velocity v was averaged over the ensemble of 2 10 = 1024 trajectories, each starting with a different initial condition according to the above distributions. The number of realisations of stochastic dynamics is not accidental and was chosen carefully to maximise the performance of the numerical simulation, see Ref. [31] for more details. The time span of the simulations was set to [0, 10 4 T] to guarantee that the directed velocity v relaxed to its asymptotic long time stationary value. 
A. S. carried out numerical calculations. All authors contributed to the discussion and analysis of the results. J. S. wrote the manuscript.