Access the full text.
Sign up today, get DeepDyve free for 14 days.
Loading next page...
/lp/wiley/analyte-driven-clustering-of-bio-conjugated-magnetic-nanoparticles-WcnwT55P60
- Publisher
- Wiley
- Copyright
- © 2023 Wiley‐VCH GmbH
- eISSN
- 2513-0390
- DOI
- 10.1002/adts.202200796
- Publisher site
- See Article on Publisher Site
Abstract
IntroductionEarly diagnosis and monitoring the progression of a disease are crucial for its successful medical treatment. In many cases the method of choice for the detection are immunoassays, which provide high accuracy and specificity. For example, for screening of cardiac biomarkers an enzyme‐linked immunoabsorbent assay (ELISA) kits are commonly used.[1] The limitation of this method is that it is time consuming[2] as well as it is not sensitive enough to very low biomarker concentrations. It is therefore of paramount importance to design rapid and sensitive immunoassays.Owing to their excellent physico‐chemical properties, magnetic nanoparticles are becoming increasingly popular as transducers in the field of nanomedicine.[3] These so called immunomagnetic assays are promising candidates to overcome the above mentioned limitations of other assay techniques. The basic idea behind immunomagnetic assays is to attach recognition ligands onto magnetic nanoparticles to specifically bind target biomarkers. Then nanoparticles begin to form clusters due to the molecular recognition between biomarkers and conjugated ligands. The formed nanoassemblies have macroscopic magnetic response, which differs from the one related to individual nanoparticles, allowing to display the biomarker detection.Different experimental techniques are used to detect these changes such as AC susceptometry[3] and magnetometry,[4] relaxometry[5,6] or giant magneto‐resistance.[7] Using magnetorelaxometry it was shown[5,6] that relaxation amplitudes for a dispersion of biotinylated magnetic nanoparticle increase as one increases the concentration of avidin (analyte). Furthermore, saturation magnetization was observed to increase as one increases the concentration of c‐reactive proteins.[8] Remanence magnetization measurements[9] have also been performed to quantify the biomarker concentration.In case of immunomagnetic reduction assays (IMR) the AC susceptibility measurements[3,10–13] are performed. Due to the clustering and the consequently increased rotational viscosity one generally observes a decrease in the magnetic susceptibility as one increases the concentration of the analyte. Another related quantity that similarly decreases as one increases the concentration of the analyte is the area of the dynamic hysteresis.While there have been many experimental studies on IMR, numerical simulations, which could provide additional insight into the mechanisms involved, are so far still lacking. In this paper, we couple a mesoscopic particle‐based simulation method (DPD) with the spin dynamics (SD) of the nanoparticle core. SD methods coupled to molecular dynamics (MD), have been used in studies of coupling of spin and lattice dynamics.[14,15] See ref. [16] for a symplectic algorithm based on Trotter–Suzuki decomposition. A related MD–SD hybrid method has been developed in studies of giant magnetoresistance effect of nanoparticles suspended in a liquid gel,[17–19] where in contrast to our work, point‐particle approximation was used. In addition we distinguish between the dynamics of the orientation of the nanoparticle and its magnetic moment, which is crucial in a certain range of amplitudes and frequencies of the external field.[20,21]The presented hybrid DPD–SD method serves as a first step toward realistic IMR simulations and other biomedical applications based on magnetic nanoparticles.[22] We lay the experimental groundwork by characterizing the iron oxide nanoparticles (IONPs) to be used in subsequent studies, that is, we measure their size distribution and magnetic response, see Methods. In this work we establish the computational model of bio‐conjugated nanoparticles suspended in a solvent containing target proteins. We first discuss the behavior of a single nanoparticle under nonequilibrium conditions induced by a pressure difference across the system to gain insight into its rheological properties. In the next step we study the cluster formation as a function of the ratio of the analyte and recognition ligand concentrations. Using the developed simulation technique, we aim to show the connection between the amount of analyte and the change in the structure of the nanoparticle suspension. We characterize this change with the hydrodynamic radius, on which the magnetic response strongly depends.Results and DiscussionWe first discuss the rheological properties of the solvent. Then we consider a single coated nanoparticle with superparamagnetic behavior (i.e., zero magnetization in absence of external magnetic field) suspended in the solvent. We calculate its hydrodynamic radius and study the dynamics of the nanoparticle under non‐equilibrium conditions. Finally, we consider a suspension of multiple nanoparticles with complete magnetic interaction and study the process of nanoparticle clustering.SolventTo put things into a physical context we first examine the physical properties of the DPD solvent, that is, water. We determine the self‐diffusion constant of the water beads as well as its viscosity. The viscosity is an important parameter as it influences the response of the dispersion of magnetic nanoparticles under oscillatory magnetic fields.[20,23]To calculate the self‐diffusion constant we determined the mean squared displacement ri(t)$r_i(t)$ (MSD) of the water beads in a simulation box with periodic boundary conditions and volume/dimensions V=30rc×30rc×30rc$V = 30r_\text{c}\times 30 r_\text{c}\times 30 r_\text{c}$1Dw=limt→∞16Nwt∑i=1Nw⟨[ri(t)−ri(0)]2⟩$$\begin{equation} D_\text{w}=\lim _{t\rightarrow \infty }\frac{1}{6N_\text{w}t}\sum _{i=1}^{N_\text{w}}\langle [r_i(t)-r_i(0)]^2\rangle \end{equation}$$where Nw$N_\text{w}$ is the number of DPD water beads and rc$r_\text{c}$ is the DPD length scale, Equation (20). We obtain Dw=(0.0507±0.0002)rc2/τ$D_\text{w}=(0.0507\pm 0.0002){r_\text{c}^2}/{\tau }$, with τ the DPD time scale Equation (24). This is equal to Dw=(5.96±0.02)×10−9$D_\text{w}=(5.96\pm 0.02)\times 10^{-9}$,m2 s−1 in physical units, which is about 2.6 times larger than experimental self‐diffusion of water at 25°.[24]The dynamic viscosity η is calculated using the open boundary molecular dynamics (OBMD) algorithm, described in Methods. Here we use a simulation box with dimensions 30rc×20rc×20rc$30 r_\text{c}\times 20 r_\text{c}\times 20 r_\text{c}$ and buffer sizes 9rc$9r_\text{c}$, Figure 1. A pair of oppositely equal forces along the y axis is introduced into the left and right buffers. In this way, a shear flow is induced across the sample, Figure 2. As can be seen, the velocity profile is linear in the region of interest and flat in the buffers with relatively small error bars, except near x=0$x=0$ and x=30rc$x=30r_\text{c}$, which is due to the low density of the particles.1FigureA schematic of the regions in the OBMD algorithm. The simulation box is open along the horizontal direction (x axis) and periodic in the other two directions.2FigureVelocity profile vy$v_y$ along the x axis of the simulation box at three different values of the applied shear stress σxy$\sigma _{xy}$. The dashed vertical lines at x=9rc$x=9r_\text{c}$ and x=21rc$x=21r_\text{c}$ represent the edges of the left and the right buffer, respectively.We then measure the shear rate in the middle of the simulation box at different values of the imposed shear stress, Figure 3, and fit the numerical values corresponding to low shear stress σxy<1.0kBT/rc3$\sigma _{xy}<1.0\, {k_{\text{B}}}T/{r_\text{c}^3}$ to the equation2σxy=ηγ̇$$\begin{equation} \sigma _{xy}=\eta \dot{\gamma } \end{equation}$$where γ̇=∂vy∂x$\dot{\gamma }=\frac{\partial v_y}{\partial x}$ is the shear rate. We obtain η=(5.23±0.04)$\eta =(5.23\pm 0.04)$ τkBT/rc3=(0.202±0.002)${\tau k_\text{B}T}/{r_\text{c}^3}=(0.202\pm 0.002)$ mPa s, which is 4.4 times smaller than real viscosity of water at 25°. Such low viscosity is consistent with the standard DPD model, which is based on a rather soft potential.[25] Figure 3 shows a noticeable shear‐thinning above γ̇≈0.2$\dot{\gamma }\approx 0.2$ τ−1$\tau ^{-1}$. Since the clustering of the nanoparticles is limited by diffusion, the typical velocities of the nanoparticles are rather low. One can therefore neglect this phenomenon.3FigureShear stress as a function of shear rate of DPD water. The numerical values are fit to Equation (2).Nanoparticle Hydrodynamic RadiusThe magnetization dynamics of superparamagnetic nanoparticles dispersed in solution is determined by two relaxation processes. The first is the Néel rotation of the magnetic moments inside the nanoparticle[26]3τN=τ0expKVckBT$$\begin{equation} \tau _{\text{N}} = \tau _0 \exp {\left(\frac{KV_\text{c}}{k_\text{B}T}\right)} \end{equation}$$where τ0 is typically between 10−12 and 10−8 s.[27] The second mechanism is the Brownian rotation of the nanoparticle, described by the Brownian relaxation time[28]4τB=4πηRh3kBT$$\begin{equation} \tau _B = \frac{4\pi \eta R_h^3}{k_BT} \end{equation}$$which depends strongly on the hydrodynamic radius Rh$R_\text{h}$ of the nanoparticle. One can change this radius by modifying the nanoparticle surface or by clustering. In this paper, we focus on the latter.The hydrodynamic radius can be estimated using several different approaches. One of these is based on the Kirkwood formula, which relates the inverse of the hydrodynamic radius to the thermal averages of the inverse distances rij$r_{ij}$ between the beads of the particle51Rh=1N2∑i≠j1rij$$\begin{equation} \frac{1}{R_\text{h}} = \frac{1}{N^2}\sum _{i\ne j}{\left\langle \frac{1}{r_{ij}}\right\rangle} \end{equation}$$where N is the number of beads of the conjugated nanoparticle, Figure 4. As expected, the hydrodynamic radius increases with the nanoparticle core radius rp$r_\text{p}$ and the polyethylene glycol (PEG) chain length nl$n_\text{l}$, Figure 5. The influence of the latter is more visible at lower particle radii, due to the lower total number of nanoparticle beads for smaller nanoparticle sizes.4FigureThe nanoparticle (red beads) is coated by a dextran shell (yellow beads) and PEG (blue beads). Recognition parts (pink beads) are located at the end of the recognition ligands.5FigureHydrodynamic radius Rh$R_\text{h}$ as a function of the nanoparticle core radius rp$r_\text{p}$ at three different lengths nl$n_\text{l}$ of the PEG chains.An alternative way of calculating the hydrodynamic radius is through the Einstein‐Stokes equation, which relates the particle self‐diffusion constant Dmd$D_{md}$ and the solvent viscosity η6Dmd=limt→∞16Npt∑i=1Np⟨[ri(t)−ri(0)]2⟩=kBT6πRhη$$\begin{equation} D_{\text{md}}=\lim _{t\rightarrow \infty }\frac{1}{6N_{\text{p}}t}\sum _{i=1}^{N_\text{p}}\langle [r_i(t)-r_i(0)]^2\rangle =\frac{k_\text{B}T}{6\pi R_\text{h} \eta } \end{equation}$$where ri(t)$r_i(t)$ is the distance of the ith nanoparticle from its initial position ri(0)$r_i(0)$ and Np$N_\text{p}$ the number of the nanoparticles. Equation (6) holds for infinitely large systems. For a finite size cubic simulation box with periodic boundary conditions the correction reads[29,30]7D=Dmd+ξkBT6πηL$$\begin{equation} D=D_{\text{md}}+\frac{\xi k_\text{B}T}{6\pi \eta L} \end{equation}$$where D is the corrected self‐diffusion constant, L is the box size and ξ is the analogue of the Madelung constant and is equal to ξ≈2.837298$\xi \approx 2.837298$. For L=60rc$L=60r_\text{c}$ this correction yields 0.0005 rc2/τ${r_\text{c}^2}/{\tau }$.The dependence of MSD is shown in Figure 6. A system of Np=30$N_{\text{p}}=30$ nanoparticles in a box with dimensions 60rc×60rc×60rc$60r_\text{c}\times 60 r_\text{c}\times 60r_\text{c}$ was simulated for t=1500τ$t=1500\tau$. Averaging and error estimation was performed on blocks of length 50τ. Due to the small number of particles, a rather large error of the MSD curve at high lag times is observed, see Figure 6. Using Equations (6) and (7), we get for the hydrodynamic radii Rh=(2.92±0.12)rc$R_\text{h}=(2.92\pm 0.12)\,r_\text{c}$, Rh=(3.96±0.16)rc$R_\text{h}=(3.96\pm 0.16)\,r_\text{c}$, Rh=(5.64±0.24)rc$R_\text{h}=(5.64\pm 0.24)\,r_\text{c}$, which are larger than what Equation (5) gives. This could be due to the hydrodynamic interactions between the nanoparticles. Assuming an average box size L/Np1/3$L/N_\text{p}^{1/3}$ for a single particle, the finite‐size correction becomes larger and equals 0.0015 rc2/τ${r_\text{c}^2}/{\tau }$. The corrected hydrodynamic radii are then a lot closer to those given by Equation (5): Rh=(2.11±0.07)rc$R_\text{h}=(2.11\pm 0.07)\,r_\text{c}$, Rh=(2.61±0.08)rc$R_\text{h}=(2.61\pm 0.08)\,r_\text{c}$, Rh=(3.25±0.10)rc$R_\text{h}=(3.25\pm 0.10)\,r_\text{c}$.6FigureMean squared displacement of the nanoparticles as a function of lag time t at three different nanoparticle core radii rp$r_\text{p}$.Non‐Equilibrium ConditionsTo study the behavior of the nanoparticle under non‐equilibrium conditions, we impose a constant flow along the x axis. We can achieve such a flow by applying, in addition to the ambient pressure force, a uniformly distributed force on the left buffer (Figure 1) pointing toward the region of interest. Throughout this numerical experiment we fix the position of the nanoparticle in the middle of the box to prevent it from leaving the simulation box. Experimentally one could achieve this by applying a sufficiently strong inhomogeneous magnetic field.The flow is approximately constant, except in the center near the nanoparticle, Figure 7. The surfactant chains on the nanoparticle surface are bent toward the flow direction, Figure 8.7FigureVelocity field of water in the xy$xy$ plane averaged in a thin slice around the layer at z=Lz/2$z=L_z/2$. The nanoparticle core of radius rp=2.0rc$r_\text{p}=2.0\,r_\text{c}$ is depicted here as a red circle in the middle, while the dashed circle represents its hydrodynamic size Rh≈2.7rc$R_\text{h}\approx 2.7\,r_\text{c}$.8FigureShape of nanoparticle coating when exposed to flow along the x axis. The PEG chains (blue beads) are visibly bent in the flow direction, in contrast to dextran (yellow beads).From the induced flow around the nanoparticle, one can in principle determine the hydrodynamic radius.[31] More precisely, the tangential flow velocity as a function of distance perpendicular to the flow direction is given by[32]8vt(r)=vs+(v∞−vs)1−3Rh4re−rReRh−Rh34r3,$$\begin{equation} v_\text{t}(r) = v_\text{s} + (v_\infty -v_\text{s}){\left(1-\frac{3R_\text{h}}{4r}e^{-\frac{r\text{Re}}{R_\text{h}}}-\frac{R_\text{h}^3}{4r^3}\right)}, \end{equation}$$where vs$v_\text{s}$ is the slip velocity, v∞$v_\infty$ is the velocity far away from the nanoparticle and Re is the Reynolds number. Equation (8) includes an inertial correction given by Oseen in 1910. Another hallmark of the finite Re is a slightly asymmetric velocity profile as can be seen in Figure 7. All these parameters can be seen as fit parameters for the calculated tangential velocity, Figure 9. We used the previously determined hydrodynamic radius Equation (5) to determine the slip velocity at the surface of the nanoparticle. As can be seen from this figure, there is a rather thick layer around the nanoparticle with zero tangential velocity. This corresponds to a depletion layer with a thickness of about rc$r_\text{c}$ with almost zero density of water beads (not shown) around the nanoparticle. This is caused by the relative high density of the nanoparticle beads as well as the dextran layer. One can also see that for our choice of model parameters the no‐slip condition is satisfied.9FigureTangential velocity vt$v_\text{t}$ of water as a function of the distance from the center of the nanoparticle r in the direction perpendicular to the flow. The velocity profiles for three different nanoparticle radii are fitted to Equation (8).In passing, we comment on the question of slip or no‐slip condition at the boundary of the nanoparticle and its physical relevance. It has been shown that the interaction between the solvent and the nanostructures and the corresponding boundary conditions play an important role on the nanoscale, such as an enhancement of flow in carbon nanotubes[33,34] as well as tunability of the dielectric response of water.[35] The boundary condition at the nanoparticle surface can strongly influence the behavior of magnetic nanoparticle suspensions under certain conditions (e.g., frequency and amplitude of the external magnetic field). For a slip‐corrected solution,[36] the Brownian relaxation time Equation (4) should read9τB=4πηRh4(Rh+3λ)kBT$$\begin{equation} \tau _\text{B} = \frac{4\pi \eta R_\text{h}^4}{(R_\text{h}+3\lambda )k_\text{B}T} \end{equation}$$where λ is the slip length, which can vary from 0 for no‐slip flow and up to infinity for perfect slip. For a large slip length, this expression therefore yields a very small relaxation time. Given that the overall relaxation dynamics of the magnetization is determined by the faster process (either Brownian or Néel),[37] this fact could be experimentally realized by, for example, modifying the nanoparticle surface.Nanoparticle Cluster FormationWe now consider a suspension of multiple nanoparticles with complete magnetic interaction and magnetization dynamics treated by the Landau–Lifshitz–Gilbert equation (details in Methods). This empowers us to study the process of nanoparticle clustering, which is central to the detection of analyte amount through magnetic measurements.The key feature of the model that enables clustering is the dynamic binding of biomarkers to the recognition ligands. Clustering, in turn, alters the rheological properties and the response to the external magnetic field. In this work, we set the valency of the biomarkers to 4. This means that a single biomarker can bind to four different recognition ligands, not necessarily from different nanoparticles. On the other hand, a single recognition ligand can bind to only one biomarker.In the present study, we simulated Np=20$N_\text{p}=20$ coated nanoparticles with nanoparticle radius rp=1rc$r_\text{p}=1r_\text{c}$ in a simulation box of dimensions 20rc×20rc×20rc$20r_\text{c}\times 20r_\text{c} \times 20r_\text{c}$ and periodic boundary conditions. This yields a volume concentration of about 1%.The time dependences of the number of unbound biomarkers (Figure 10) and nanoparticles not bound to other nanoparticles (Figure 11) show a quick decrease at small times and a slow decrease at large times. The binding dynamics and the resulting cluster configurations depend on the ratio cr=cbio/clig$c_\text{r}={c_{\text{bio}}}/{c_{\text{lig}}}$ of the concentrations of biomarkers cbio$c_{\text{bio}}$ and recognition ligands clig$c_{\text{lig}}$. For low cr$c_\text{r}$, one expects mostly singlets (unbound nanoparticles) and doublets (pairs of bound nanoparticles). For larger cr$c_\text{r}$, triplets, quadruplets, … start to form and the number of doublets decreases, Figure 12.10FigureNormalized concentration of free unbound biomarkers as a function of time at three different ratios cr$c_\text{r}$.11FigureNormalized concentration of singlets cs$c_\text{s}$ (unbound nanoparticles) as a function of time at three different ratios cr$c_\text{r}$.12FigureNormalized concentration of singlets cs$c_\text{s}$ (unbound nanoparticles) and doublets cd$c_\text{d}$ (pairs of bound nanoparticles) as a function of time at cr=0.73$c_\text{r}=0.73$.Moreover, we find that the concentration of singlets decreases with increasing cr$c_\text{r}$ until cr≈1$c_\text{r} \approx 1$, where there are no more unbound nanoparticles, Figure 13. If cr$c_\text{r}$ is further increased, the number of singlets increases again. Namely, if there are too many biomarkers, by the time the nanoparticles have diffused to each other, the biomarkers will have occupied all available sites. Nanoparticles with all recognition ligands occupied cannot bind.13FigureNormalized concentration of singlets cs$c_\text{s}$ as a function of the concentration ratio cr$c_\text{r}$ at two different values of magnetic interaction parameter λ.We studied the clustering of nanoparticles in absence of external magnetic fields. The simulations at two different values of the nanoparticle magnetic moment, corresponding to dipolar interaction strengths λ=0.006$\lambda =0.006$ and λ=60$\lambda =60$, Equation (34), do not show any major differences.The characteristic time of the clustering process can be roughly estimated by considering the times required for two nanoparticles to approach each other (ttr$t_{\text{tr}}$) and the time for the nanoparticles to rotate sufficiently (trot$t_{\text{rot}}$) for binding of a biomarker, bound to one of the nanoparticles, to a free recognition ligand on another nanoparticle.The characteristic time ttr$t_{\text{tr}}$ is given by the translational diffusion and the concentration of the nanoparticles[38]10ttr=x26πRhηkBT$$\begin{equation} t_{\text{tr}} = \frac{x^26\pi R_\text{h} \eta }{k_\text{B}T} \end{equation}$$where x=(VNp)1/3$x=(\frac{V}{N_\text{p}})^{1/3}$. As can be seen from Equation (10), the clustering is faster (i.e., tr$t_\text{r}$ is smaller) when using solvent with lower viscosity η. In addition to a different magnetic response, another physical consequence of lower viscosity is a larger translational self‐diffusion coefficient of the nanoparticles as well as that of the biomarkers. We note that the viscosity of our DPD solvent can be matched to the experimental viscosity of water by using an extended DPD thermostat, see ref. [39], which will be considered in future work.The characteristic time trot$t_{\text{rot}}$ is much more difficult to estimate, since it depends on the number of ligands per nanoparticle Nl$N_\text{l}$, the distribution of the number of bound biomarkers per nanoparticle as well as the radius of the nanoparticle. In the limits crNl≲1$c_\text{r} N_\text{l} \lesssim 1$ (multiple biomarker bonds per nanoparticle improbable) and cr≳1$c_\text{r} \gtrsim 1$ (more than one free ligand per nanoparticle improbable) it can be estimated from the rotational diffusion constant alone:11trot=8πηRh3kBT.$$\begin{equation} t_{\text{rot}} = \frac{8\pi \eta R_\text{h}^3}{k_\text{B}T}. \end{equation}$$For Rh=2rc$R_\text{h}=2\,r_\text{c}$ we get ttr≈104τ=68$t_{\text{tr}}\approx 10^4\,\tau = 68$ ns and trot≈1050τ=7$t_{\text{rot}}\approx 1050\,\tau =7$ ns. In this regime of small nanoparticles the clustering is therefore determined by translational diffusion.As expected, the hydrodynamic radius of clusters also depends strongly on the concentration ratio cr$c_\text{r}$, Figure 14. It increases with increasing cr$c_\text{r}$ until about cr≈1$c_\text{r}\approx 1$, where it begins to decrease due to the over‐saturation of biomarkers. This is a well‐known phenomenon in nanoparticle‐based immunoassays and is commonly referred to as the hook effect[40] or the postzone (in the case of antigen excess) and prozone effect (in the case of recognition ligand excess). Such an effect may lead to false‐negative results.14FigureNormalized hydrodynamic radius as a function of the concentration ratio cr$c_\text{r}$ at two different values of magnetic interaction parameter λ.Increasing the magnetic moment does not seem to change the simulation results significantly, which could be due to the small number and size of the nanoparticles. The hydrodynamic radius at cr≈1$c_\text{r}\approx 1$ might indicate that dipolar interactions do not favor a large cluster. More specifically, for λ=0.006$\lambda =0.006$ we observe that all 20 nanoparticles in the system cluster together to form a complete graph of bonds, Figure 15. This is known as gelation, a phenomenon often observed in coagulation dynamics and described by the Smoluchowski equations.[41,42] A more thorough study of the gelation and the effects of the nanoparticle size and magnetic fields on the hydrodynamic radius will be discussed in our future work. We mention that the nanoparticle radius r=1.0rc$r=1.0\,r_\text{c}$ used for simulating the clustering dynamics of the nanoparticles is rather small for considering Brownian dominating nanoparticles, which are typically used in immunomagnetic assay. We used a small radius simply to increase the total number of particles involved in the clustering dynamics. The fact that there exists an optimal value of the concentration ratio of the biomarkers and ligands for a large hydrodynamic radius, should not be affected significantly by the nanoparticle size. On the other hand, the nanoparticle size affects the timescales of the clustering dynamics and the behavior of the cluster in oscillatory magnetic fields. This will be covered in our subsequent study.15FigureNanoparticle system at cr≈1$c_\text{r}\approx 1$ and λ=0.006$\lambda =0.006$ forms a complete graph of connections (gelation). Due to the periodic boundary conditions the nanoparticle on the left hand side appears to be isolated, but is actually a part of the same cluster.ConclusionsIn this work we have combined the dissipative particle dynamics with the Landau–Lifshitz–Gilbert equation to simulate a dispersion of conjugated magnetic nanoparticles. To the best of our knowledge such a hybrid DPD‐SD solver that takes into account also the orientation (easy axis) and the rigid body dynamics of nanoparticles as well as their binding in the presence of analytes has not been considered before.With such a numerical set‐up we were able to simulate the nanoparticle dispersion in presence and absence of biomarkers. We have demonstrated that the presence of biomarkers promotes the clustering of magnetic particles. For low biomarker concentrations, increasing the concentration increases the hydrodynamic radius. Above a certain critical biomarker concentration the hydrodynamic size starts to decrease again, which is also known as the hook effect in immunoassay studies.As a future prospect we have developed a coupling between the lattice‐Boltzmann method[43] (LBM) to describe the solvent, the molecular dynamics of the nanoparticle cores and conjugated surface and magnetic dynamics of nanoparticle cores. The nanoparticles are coupled to the LBM fluid through the so called immersed boundary method (IBM).[44,45] This allows one to simulate larger spatio‐temporal scales, which is essential to simulate dispersions of realistic 50 nm sized conjugated nanoparticles. The simulation results will be compared to our experimental measurements of the magnetic response at different analyte concentrations.MethodsModel of Conjugated Magnetic NanoparticlesThe simulated system consisted of coated spherical magnetic nanoparticles, suspended in water. The particles were surrounded by biomarkers, which may attach to the nanoparticle coating.The construction of magnetic nanoparticles is first discussed. To properly capture the rotational behavior, more than one coupling point per nanoparticle was needed.[46] The nanoparticle cores were constructed by placing a large number of DPD beads at a fixed radius rp$r_\text{p}$ from the center. This is the so called raspberry model,[47,48] which was convenient for the attachment of the nanoparticle coating and ligands. These beads were approximately homogeneously distributed on the surface of the nanoparticle, see Figure 16. In addition, one bead was put at the center of the particle, which carried the magnetic moment m$\mathbf {m}$ of the entire nanoparticle core. This is the so called macrospin approximation, which was relevant for single‐domain nanoparticles smaller than 40 nm.[49,50]16FigureThe nanoparticle core is constructed by placing many beads at a fixed radius from the center. The core is magnetic and is within the macrospin approximation represented by a magnetic moment m$\mathbf {m}$ (blue arrow). The preffered orientation of the magnetic moment of an isolated magnetic particle in absence of external fields points along the easy axis n$\mathbf {n}$ (orange arrow).The nanoparticles were decorated with a surfactant, in the case a shell of dextran, Figure 4. The physical advantages of such a shell were biocompatibility with the aqueous environment, a decreased toxicity and increased stability of the nanoparticle dispersion.[51,52] The nanoparticles were in addition conjugated with multiple recognition ligands. Each of these ligands consisted of a longer chain of a certain polymer, in the case polyethylene glycol (PEG), and a recognition part at the end.The surfactants and the ligands were modeled by constructing chains of beads, connected by harmonic springs12Vb(r)=Kb(r−r0)2$$\begin{equation} V_\text{b}(r) = K_\text{b}(r-r_0)^2 \end{equation}$$where Kb$K_\text{b}$ is the elastic constant, r0 is the equilibrium bond distance, and r is the distance between neighboring beads within a chain. A harmonic potential was also added between a chosen nanoparticle bead and the chains. The nanoparticle core beads to which the dextran chains are attached were chosen randomly at the start of the simulation and were covering 75% of the nanoparticle core surface. The same bond constants were used for all types of beads. Note that the usual factor 12$\frac{1}{2}$ is already included in the bond constant Kb$K_\text{b}$. The number of beads used for each chain depended on their length. In the case the dextran shell consisted of chains with two beads, while the ligands consisted of seven beads of PEG and one bead corresponding to the recognition part. In all of the simulations six such recognition ligands were used per nanoparticle, which were distributed on the sphere according to the Thomson problem of point charges[53] to prevent any effects of inhomogeneous distribution of ligands. Furthermore, a bending potential Vθ$V_{\theta }$ is added between neighboring bonds within the chain13Vθ(θ)=Kθ(θ−θ0)2$$\begin{equation} V_{\theta }(\theta ) = K_\theta (\theta -\theta _0)^2 \end{equation}$$where θ is the bond angle and θ0 is the equilibrium bond angle. The equilibrium shape of the chain was determined by the complicated interplay of the interactions between the rest of the particles in the systems, especially the interaction between the ligands and the solvent. For PEG, the bond parameters were Kb=8500$K_\text{b}=8500$ kJ mol−1 nm−2, Kθ=42.5$K_\theta =42.5$ kJ mol−1, r0=3.30$r_0=3.30$ Å, and θ0=130∘$\theta _0=130^\circ$, which are taken from the MARTINI coarse‐grained model.[54] In the DPD model, three monomer units of PEG are represented as one bead (HO–[CH2CH2O]3–H). This was to ensure that different types of beads (water, dextran, and PEG) have a similar bead volume. A similar model for PEG as a surfactant has also been used in DPD simulations of endocytosis of PEGylated nanoparticles.[55] For dextran, the so called M3B coarse grained model was often used, where each monomer unit of dextran was represented by three beads.[56] Since the dextran chains were rather small compared to the ligands and were mainly there to prevent overlap of nanoparticles and to some extent to prevent water leaking inside nanoparticles, the same bond parameters Kb,Kθ,r0$K_\text{b}, K_\theta , r_0$ and θ0 were taken as for PEG.Dissipative Particle DynamicsThe simulations were performed using DPD method,[57,58] which was a mesoscopic method often used for simulating polymers,[25] biological systems such as cell membranes[59,60] and colloids.[61] It was recently also used in modeling of ultrasound propagation in water.[62]The system was represented by DPD beads, which interacted with each other via a soft repulsive potential, also called a conservative force FijC$F_{ij}^\text{C}$, where i and j are the indices of the ith and jth bead, respectively. The force, which acted on a given bead consisted also of a dissipative force FijD$F_{ij}^\text{D}$ and a random force FijR$F_{ij}^\text{R}$, which acted together effectively as a thermostat for the system. DPD thermostat conserved linear and angular momentum, which was suitable for modeling hydrodynamic systems. The total non‐bonded force on the ith particle, Fi$\mathbf {F}_i$, is thus a sum of these three forces14Fi=∑i≠jFijC+FijD+FijR$$\begin{equation} \mathbf {F}_i = \sum _{i\ne j} \mathbf {F}_{ij}^\text{C} + \mathbf {F}_{ij}^\text{D} + \mathbf {F}_{ij}^\text{R} \end{equation}$$where15FijC=aijω(rij)r̂ij$$\begin{align} \mathbf {F}_{ij}^\text{C} &= a_{ij}\omega (r_{ij})\mathbf {\hat{r}}_{ij} \end{align}$$16FijD=−γωD(rij)(r̂ij·vij)r̂ij$$\begin{align} \mathbf {F}_{ij}^\text{D} &=-\gamma \omega _\text{D}(r_{ij})(\mathbf {\hat{r}}_{ij}\cdot \mathbf {v}_{ij})\mathbf {\hat{r}}_{ij} \end{align}$$17FijR=σωR(rij)Θijr̂ij$$\begin{align} \mathbf {F}_{ij}^\text{R} & = \sigma \omega _\text{R}(r_{ij})\Theta _{ij}\mathbf {\hat{r}}_{ij} \end{align}$$where aij$a_{ij}$ is the interaction parameter between beads i and j, γ is the friction parameter, σ is the noise strength, ωR$\omega _\text{R}$ and ωD$\omega _\text{D}$ are weight functions, r̂ij$\mathbf {\hat{r}}_{ij}$ denotes the normalized vector of the inter‐particle axis rij=ri−rj$\mathbf {r}_{ij} = \mathbf {r}_i-\mathbf {r}_j$, vij=vi−vj$\mathbf {v}_{ij}=\mathbf {v}_i-\mathbf {v}_j$ is the relative velocity of the beads and Θij=Θji$\Theta _{ij}=\Theta _{ji}$ is a Gaussian random variable with zero mean, ⟨Θij(t)⟩=0$\langle \Theta _{ij}(t) \rangle = 0$, and the following second moment18⟨Θij(t)Θkl(t′)⟩=(δikδjl+δilδjk)δ(t−t′)$$\begin{equation} \langle \Theta _{ij}(t) \Theta _{kl}(t^{\prime })\rangle = (\delta _{ik}\delta _{jl}+\delta _{il}\delta _{jk})\delta (t-t^{\prime }) \end{equation}$$From the fluctuation–dissipation theorem one obtained the relation ωD(r)=ωR2(r)$\omega _\text{D}(r) = \omega _\text{R}^2(r)$ and σ2=2kBTγ$\sigma ^2 = 2k_\text{B}T\gamma$. Typically, ωR(rij)$\omega _\text{R}(r_{ij})$ is chosen to be equal to ω(rij)$\omega (r_{ij})$19ω2(rij)=ωD(rij)=ωR2(rij)=(1−rijrc)2,rij<rc0,rij≥rc$$\begin{equation} \omega ^2(r_{ij}) =\omega _\text{D}(r_{ij}) = \omega _\text{R}^2(r_{ij}) = {\left\lbrace \def\eqcellsep{&}\begin{array}{ll}(1-\frac{r_{ij}}{r_\text{c}})^2, & \quad r_{ij} < r_\text{c} \\[3pt] 0, & \quad r_{ij} \ge r_\text{c} \end{array} \right.} \end{equation}$$DPD is a mesoscopic method, which involved representing multiple atoms or molecules as a single bead. To find the length scale of the system one matched the mass densities of a DPD coarse grained description and the real density of water ρw$\rho _\text{w}$[59]20rc=(vmρDPDNm)1/3≈3.107(ρDPDNm)1/3[Å]$$\begin{equation} r_\text{c} = (v_\text{m} \rho _{\text{DPD}} N_\text{m})^{1/3} \approx 3.107\,(\rho _{\text{DPD}}N_\text{m})^{1/3} [\mathrm{\mbox{{\AA}}}] \end{equation}$$where vm=MwρwNA≈30$v_\text{m}=\frac{M_\text{w}}{\rho _\text{w} N_\text{A}}\approx 30$ Å3 is the molecular volume of water with Mw≈18$M_\text{w}\approx 18$ g mol−1 the molar mass of water, NA$N_\text{A}$ is the Avogadro constant, ρw=997$\rho _\text{w}=997$ kg m−3 is the mass density of water at 25∘$25\text{ }^\circ$C and 1 atm. In Equation (20) the mapping number Nm$N_\text{m}$ represents the number of molecules per DPD bead while ρDPD$\rho _{\text{DPD}}$ is the number density of the DPD beads. It was noted that the mass of a single bead is m=NmmH2O$m = N_\text{m} m_{\text{H}_2\text{O}}$, where mH2O≈18$m_{\text{H}_2\text{O}} \approx 18$ Da. The mass of a single water bead was taken as unit mass in this work. Beads representing PEG and dextran had masses m=0.92$m=0.92$ and m=1.25$m=1.25$, respectively. Equation (20) can also be derived by matching the volume occupied by one DPD bead21VDPD=rc3ρDPD=Nmvm$$\begin{equation} V_{\text{DPD}} = \frac{r_\text{c}^3}{\rho _{\text{DPD}}} = N_\text{m} v_\text{m} \end{equation}$$The DPD number density of beads ρDPD$\rho _{\text{DPD}}$ as well as the mapping Nm$N_\text{m}$ were free parameters that determine the length scale rc$r_\text{c}$. It should be mentioned that one can start the other way around by first setting the scale rc$r_\text{c}$ and then deriving the mapping number Nm$N_\text{m}$. The density was set to ρDPD=3.0rc−3$\rho _{\text{DPD}} = 3.0\,r_\text{c}^{-3}$ throughout the work. For densities above ρDPD=2rc−3$\rho _{\text{DPD}} = 2\,r_\text{c}^{-3}$, the following equation of state is valid[25]22p=ρDPDkBT+αaWWρDPD2$$\begin{equation} p = \rho _{\text{DPD}}k_\text{B}T + \alpha \, a_{\text{WW}}\rho _{\text{DPD}}^2 \end{equation}$$where α=(0.101±0.001)rc4$\alpha = (0.101 \pm 0.001)\,r_\text{c}^4$ is a numerical factor and aWW$a_{\text{WW}}$ is the interaction parameter for the beads of the solvent (e.g., water).To determine the interaction parameter aWW$a_{\text{WW}}$ one typically matched the dimensionless isothermal compressibility κ to the experimental one for water (κ−1≈16$\kappa ^{-1}\approx 16$). Using Equation (22) and κ−1=1kBT(∂p∂n)exp=1NmkBT(∂p∂ρDPD)sim$\kappa ^{-1} = \frac{1}{k_\text{B}T}(\frac{\partial p}{\partial n})_{\text{exp}}=\frac{1}{N_\text{m}k_\text{B}T}(\frac{\partial p}{\partial \rho _{\text{DPD}}})_{\text{sim}}$, with n the number density of water molecules, one can solve for aii$a_{ii}$23aWW=Nmκ−1−12αρDPDkBT$$\begin{equation} a_{\text{WW}} = \frac{N_\text{m} \kappa ^{-1}-1}{2\alpha \rho _{\text{DPD}}}k_\text{B}T \end{equation}$$The time scale is then24τ=rcmkBT$$\begin{equation} \tau = r_\text{c} \sqrt {\frac{m}{k_\text{B}T}} \end{equation}$$For the dissipative parameter γ, γ=4.5mikBT/rc=4.5τkBT/rc2=4.5$\gamma = 4.5\,\sqrt {m_ik_\text{B}T}/r_\text{c}=4.5\, \tau k_\text{B}T/r_\text{c}^2=4.5$ was chosen, which gave σ2=9τ(kBT)2/rc2=9$\sigma ^2 = 9\,\tau (k_\text{B}T)^2/r_\text{c}^2=9$. In this work the unit of energy was measured in kBT$k_\text{B}T$, the unit of mass is NmmH2O$N_\text{m} m_{\text{H}_2\text{O}}$, the length is measured in rc$r_\text{c}$ and the time in τ. For the case of Nm=8$N_\text{m}=8$ one gets: rc≈0.89$r_\text{c} \approx 0.89$ nm, τ≈6.8$\tau \approx 6.8$ ps and aWW≈210kBTrc$a_{\text{WW}} \approx 210\,\frac{k_\text{B}T}{r_\text{c}}$. The integration time step δt$\delta t$ was chosen to be 0.005τ$0.005\,\tau$.The choice of the interaction parameter for water aWW$a_{\text{WW}}$ was therefore determined by the experimental value of the isothermal compressibility κ and the level of coarse‐graining, given by the mapping number Nm$N_\text{m}$. It should be mentioned that this mapping was limited to about Nm=10$N_\text{m}=10$, above which one typically observed solidification effects.[63] A different approach was needed to estimate other interaction parameters, see Table 1. Following the model of Groot and Warren, interaction parameters between the beads of the same type were set to the same value, in this case aii=aWW≈210kBTrc$a_{ii} = a_{\text{WW}} \approx 210 \frac{k_\text{B}T}{r_\text{c}}$ was taken. An alternative method estimated the interaction parameters for the like‐beads directly from solubility parameters.[64]1TableTable of interaction parameters aij$a_{ij}$ (units kBTrc$\frac{k_\text{B}T}{r_\text{c}}$). Acronyms: W = water, Mag = nanoparticle beads, Dex = dextran, Lig = ligand, Ab = recognition parts, Bio = biomarkerWMagDexLigAbBioW210105240320320320Mag210105105105210Dex210230230230Lig210210210Ab210210Bio210For the nanoparticle surface one must ensure that it was impermeable to the solvent and the surfactants. Since the surface number density σs$\sigma _\text{s}$ was rather large σs=4$\sigma _\text{s}=4$ rc−2$r_\text{c}^{-2}$ it was found that the value 0.5aWW$0.5\,a_{\text{WW}}$ was sufficient. Using an interaction parameter that was too large leads to large density oscillations near the rigid walls.[65,66]The interaction between the ligand beads and the water beads was estimated using the Flory–Huggins parameter χij$\chi _{ij}$25aij=aii+3.497χijkBTrc$$\begin{equation} a_{ij}=a_{ii}+3.497\chi _{ij}\frac{k_\text{B}T}{r_\text{c}} \end{equation}$$which was valid for the value of the bead number density used in this work, ρDPD=3$\rho _{\text{DPD}}=3$.[25,59] The Flory–Huggins parameter can be estimated from Hildebrand solubility parameters δi$\delta _i$,[67,68] see Table 226χij=(δi−δj)2(Vi+Vj)2kBT$$\begin{equation} \chi _{ij}=\frac{(\delta _i-\delta _j)^2(V_i+V_j)}{2k_\text{B}T} \end{equation}$$where Vi$V_i$ is the volume of bead i.2TableSolubility parameters and bead volumesBead typeδ [(Jcm−3)1/2$(\mathrm{J}\,\mathrm{cm}^{-3})^{1/2}$]Bead volume [Å3]Water[69]47.5 ± 0.4240.0PEG[69]23.7 ± 0.3221.3Dextran[70]34.8200.0Open Boundary Molecular DynamicsTo study the behavior of the conjugated nanoparticle under nonequilibrium conditions the OBMD method was employed.[71–73] This was an algorithm that allowed one to simulate grand‐canonical ensembles, where the system exchanges energy and mass with the environment. The number of particles in the simulation box was therefore not constant. In this paper open boundaries were employed along the x axis, while the simulation box had periodic boundary conditions in the other two directions. The simulation box was split into three regions: a region of interest in the middle surrounded by two buffers (left and right), where the insertion and the deletion of particles was performed, see Figure 1. At each time step the chemical potential was maintained by adding the particles to the buffers according to the following dynamic equation27ΔN=−δtτB(N−αN0)$$\begin{equation} \Delta N = -\frac{\delta t}{\tau _\text{B}}(N-\alpha N_0) \end{equation}$$where τB$\tau _\text{B}$ is the buffer relaxation time, which is typically τB≈100δt$\tau _\text{B} \approx \ 100\,\delta t$,[71] ΔN$\Delta N$ is the number of particles to be inserted or deleted, N0 is the equilibrium number of particles in either of the buffers and α is a user‐defined parameter, which was set for stability reasons to 0.8 in this work. In short, particles were inserted if ΔN>0$\Delta N > 0$. If ΔN<0$\Delta N < 0$ the algorithm was not performed, but particles were allowed to diffuse over the open boundaries and were deleted as soon as they leave the simulation box.The insertion was performed using USHER algorithm.[74,75] In USHER, the particle was in the first step randomly inserted in the chosen buffer region. Then an energy minimization was performed until the total energy of the inserted particle was below a certain threshold. Since DPD potential was rather soft an implementation of adaptive resolution scheme (AdResS)[76,77] was not necessary.One imposed boundary conditions by adding external forces fi$\mathbf {f}_i$ to the particles in the buffer regions and were determined by the momentum balance28JA·n=∑ifi+∑i′Δ(mi′vi′)δt$$\begin{equation} \textsf {J}A\cdot \mathbf {n} = \sum _i \mathbf {f}_i +\sum _{i^{\prime }} \frac{\Delta (m_{i^{\prime }}\mathbf {v}_{i^{\prime }})}{\delta t} \end{equation}$$where J$\textsf {J}$ is the momentum flux tensor, A is the surface area of the interface between a buffer and the region of interest, n$\mathbf {n}$ is normal to this interface and points toward the region of interest, while Δ(mi′vi′)$\Delta (m_{i^{\prime }}\mathbf {v}_{i^{\prime }})$ represents the momentum change when the particle was inserted or deleted from the buffer. The momentum flux tensor was composed of the pressure contribution and the tangential shear force29JA·n=pxxAn+σxyAt$$\begin{equation} \textsf {J}A\cdot \mathbf {n} = p_{xx}A \mathbf {n} + \sigma _{xy}A\mathbf {t} \end{equation}$$where pxx$p_{xx}$ is the pressure, t$\mathbf {t}$ points along the desired shear force (êy$\hat{\mathbf {e}}_y$) and σxy$\sigma _{xy}$ is the imposed shear stress.Hybrid SolverTo model the dynamics of the magnetic moments of the nanoparticles the Landau–Lifshitz–Gilbert equation was employed. In short, this equation described the relaxation of the magnetization toward the local magnetic field and the precession around it. Since the magnetic interactions were of long range, the calculation of this local magnetic fields involved summation of the fields of all the surrounding magnetic particles and their periodic images. The simulations were carried out using LAMMPS simulation package.[78]A solver was implemented, which incorporated the translation of the particles, their rigid body dynamics as well as the dynamics of the magnetization. The equations for the positions and the velocities of all the particles in the system were updated using the well known velocity‐Verlet algorithm[79]30v(k)(t+12δt)=v(k)(t)+12δtf(k)(t)m(k)$$\begin{eqnarray} \mathbf {v}^{(k)}(t+{\textstyle {\frac{1}{2}}}\delta t) &=& \mathbf {v}^{(k)}(t)+ \frac{1}{2}\delta t\frac{\mathbf {f}^{(k)}(t)}{m^{(k)}} \end{eqnarray}$$31r(k)(t+δt)=r(k)(t)+δtv(k)(t+12δt)$$\begin{eqnarray} \mathbf {r}^{(k)}(t+ \delta t)&=&\mathbf {r}^{(k)}(t)+\delta t\mathbf {v}^{(k)}(t+\textstyle {\frac{1}{2}}\delta t) \end{eqnarray}$$32v(k)(t+δt)=v(k)(t+12δt)+12δtf(k)(t+δt)m(k)$$\begin{eqnarray} \mathbf {v}^{(k)}(t+\delta t) &=& \mathbf {v}^{(k)}(t+{\textstyle {\frac{1}{2}}}\delta t)+ \frac{1}{2}\delta t\frac{\mathbf {f}^{(k)}(t+\delta t)}{m^{(k)}} \end{eqnarray}$$where m(k),r(k),v(k),f(k)$m^{(k)}, \mathbf {r}^{(k)}, \mathbf {v}^{(k)}, \mathbf {f}^{(k)}$ are the mass, the position, the velocity, and the total force acting on the kth bead in the system. In the last step, Equation (32), the force is calculated using positions at t+δt$t+\delta t$ and velocity at t+12δt$t+\frac{1}{2}\delta t$, that is f(k)(t+δt)≡f(k)(t+δt,r(k)(t+δt),v(k)(t+12δt))$\mathbf {f}^{(k)}(t+\delta t)\equiv \mathbf {f}^{(k)}(t+\delta t,\mathbf {r}^{(k)}(t+\delta t),\mathbf {v}^{(k)}(t+\frac{1}{2}\delta t))$. It should be emphasized that the magnetic force acting between the nanoparticles was already included in this step in the force f(k)$\mathbf {f}^{(k)}$33fm(l)=∑i≠l3μ04πr5(m(i)·r)(m(l)·r)r2(m(i)·r)m(l)+(m(l)·r)m(i)+(m(i)·m(l))r−5(m(i)·r)(m(l)·r)r2r$$\begin{eqnarray} \mathbf {f}_m^{(l)} &=& \sum _{i\ne l}\frac{3\mu _0}{4\pi r^5}{\left[\vphantom{\frac{(\mathbf {m}^{(i)}\cdot \mathbf {r})(\mathbf {m}^{(l)}\cdot \mathbf {r})}{r^2}}(\mathbf {m}^{(i)}\cdot \mathbf {r})\mathbf {m}^{(l)}+(\mathbf {m}^{(l)}\cdot \mathbf {r})\mathbf {m}^{(i)}\right.}\nonumber \\ &&+\,(\mathbf {m}^{(i)}\cdot \mathbf {m}^{(l)} )\mathbf {r} -5{\left. \frac{(\mathbf {m}^{(i)}\cdot \mathbf {r})(\mathbf {m}^{(l)}\cdot \mathbf {r})}{r^2}\mathbf {r} \right]} \end{eqnarray}$$where f(l)$\mathbf {f}^{(l)}$ stands for the magnetic force acting on the center bead of the lth nanoparticle with m(l)$\mathbf {m}^{(l)}$ the corresponding magnetic moment, μ0 is the magnetic constant μ0=4π×10−7$\mu _0 = 4\pi \times 10^{-7}$ N A−2 and r=|r|$r=|\mathbf {r}|$ is the distance between the centers of the ith and lth nanoparticle, with r=r(l)−r(i)$\mathbf {r} = \mathbf {r}^{(l)} - \mathbf {r}^{(i)}$. The magnetic moments were measured in units of m∗=4πrc3kBT/μ0≈5.4×10−21$m^*=\sqrt {4\pi r_\text{c}^3k_\text{B}T/\mu _0}\approx 5.4\times 10^{-21}$ A m2. The sum in the expression for the dipolar force, Equation (33), extended to periodic images as well. For this the particle–particle particle–mesh method was used,[80] which split the dipolar force contributions into short range, and long‐range contributions, which were efficiently calculated in the Fourier space using fast Fourier transform.It was noted that the strength of the dipolar interactions compared to the thermal energy is given by34λ=μ0|m|24πkBT(2r)3$$\begin{equation} \lambda = \frac{\mu _0|\mathbf {m}|^2}{4\pi k_\text{B}T(2r)^3} \end{equation}$$where 2r is the inter‐particle distance of two touching spheres. λ can be estimated from the typical saturation magnetizations Ms$M_\text{s}$ of the nanoparticle core materials, for example, for magnetite nanoparticles Ms≈400$M_\text{s}\approx 400$ kA m−1.[81] For a nanoparticle with radius r in DPD units, this yields |m|=0.22r3$|\mathbf {m}|=0.22\,r^3$ in units of m∗$m^*$ and λ=0.006r3$\lambda = 0.006\,r^3$. Dipolar effects become important when λ≫1$\lambda \gg 1$, which is in the case of monodomain magnetite particles equivalent to r≫5.5rc≈5$r \gg 5.5\,r_\text{c} \approx 5$ nm. In the analysis of nanoparticle clustering rather small nanoparticles with radius r=1$r=1$ rc$r_c$ was used, which means λ=0.006$\lambda =0.006$. To study also the effects of larger magnetic parameters the magnetization of the nanoparticle core was artificially increased by a factor of 100. This is equivalent to λ=60$\lambda =60$. The study of more realistic larger particle sizes will be covered in a future work.The nanoparticle cores were treated as rigid bodies with three translational and three rotational degrees of freedom. The dynamics of a given nanoparticle was then determined by the total force and the total torque on the center of mass, which yielded the updated coordinates, velocities, and angular velocities. The numerical algorithm for the angular velocity and the orientation is similar to the velocity‐Verlet algorithm35L(t+12δt)=L(t)+12δtτ(t)$$\begin{equation} \mathbf {L}(t+{\textstyle {\frac{1}{2}}}\delta t) = \mathbf {L}(t)+ \frac{1}{2}\delta t{\mbox{\boldmath $\tau $}}(t)\end{equation}$$36q(t)⟶q(t+δt)seeEquation(38)$$\begin{equation} \mathbf {q}(t) \longrightarrow \mathbf {q}(t+\delta t) \hspace{8.5359pt} \mathrm{see\nobreakspace Equation\nobreakspace (38)} \end{equation}$$37L(t+δt)=L(t+12δt)+12δtτ(t+δt)$$\begin{equation} \mathbf {L}(t+\delta t) = \mathbf {L}(t+{\textstyle {\frac{1}{2}}}\delta t)+ \frac{1}{2}\delta t{\mbox{\boldmath $\tau $}}(t+\delta t) \end{equation}$$where L$\mathbf {L}$ is the angular momentum, τ${\mbox{\boldmath $\tau $}}$ is the total torque on a nanoparticle and q$\mathbf {q}$ is a quaternion which conveniently describes the orientation and evolves according to the angular velocity Ω${\mbox{\boldmath $\Omega $}}$[82]38q̇=12Ωq$$\begin{equation} \dot{\bf q} = \frac{1}{2}{\mbox{\boldmath $\Omega $}}{\bf q} \end{equation}$$where on the right‐hand side Ω${\mbox{\boldmath $\Omega $}}$ is represented by a quaternion Ω≡(0,Ωx,Ωy,Ωz)${\mbox{\boldmath $\Omega $}}\equiv (0,\Omega _x,\Omega _y,\Omega _z)$.The central particle carried in addition the magnetic moment and an orientation, defined by the easy axis n$\mathbf {n}$, which determined the direction of the magnetic moment of an isolated nanoparticle in equilibrium. In general, due to the dipolar forces, external magnetic fields and the clustering m∦n$\mathbf {m} \nparallel \mathbf {n}$.At each step the additional orientational degree of freedom n(l)$\mathbf {n}^{(l)}$ was taken into account, by updating it according to the angular velocities Ω(l)${\mbox{\boldmath $\Omega $}}^{(l)}$ of the corresponding nanoparticle[20]39∂n(l)∂t+n(l)×Ω(l)=0$$\begin{eqnarray} \frac{\partial \mathbf {n}^{(l)}}{\partial t} +\mathbf {n}^{(l)}\times {\mbox{\boldmath $\Omega $}}^{(l)} &=& 0 \end{eqnarray}$$The rotation of the magnetic moment was determined by the magnetic relaxation and precession processes, described by the Landau–Lifshitz–Gilbert equation[83–85]40∂α(l)∂t+α(l)×Ω(l)+X=0$$\begin{eqnarray} \frac{\partial \mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}}{\partial t} + \mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times {\mbox{\boldmath $\Omega $}}^{(l)} + \mathbf {X} &=& 0 \end{eqnarray}$$where α(l)=m(l)|m(l)|${\mbox{\boldmath $\alpha $}}^{(l)} = \frac{\mathbf {m}^{(l)}}{|\mathbf {m}^{(l)}|}$ is the unit vector representing the direction of the magnetic moment of the lth nanoparticle, X$\mathbf {X}$ is the magnetization current41X=γR1+η2α(l)×H+ηγR1+η2α(l)×(α(l)×H)$$\begin{equation} \mathbf {X} = \frac{\gamma ^R}{1+\eta ^2}\mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times \mathbf {H} + \frac{\eta \gamma ^R}{1+\eta ^2} \,\mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times (\mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times \mathbf {H}) \end{equation}$$with Hi$H_i$ being the local magnetic field at the position of the central bead of the lth nanoparticle, composed of the external field Hiext$H_i^{\text{ext}}$, anisotropic contribution Hiani$H_i^{\text{ani}}$, dipole field Hidip$H_i^{\text{dip}}$, demagnetization field Hidem$H_i^{\text{dem}}$ and stochastic thermal contribution Hith$H_i^{\text{th}}$42Hi=Hiext+Hiani+Hidip+Hidem+Hith−Ωi(l)γR$$\begin{equation} H_i = H_i^{\text{ext}}+H_i^{\text{ani}}+H_i^{\text{dip}}+H_i^{\text{dem}}+H_i^{\text{th}}-\frac{\Omega _i^{(l)}}{\gamma ^\text{R}} \end{equation}$$In Equation (41) γR$\gamma ^\text{R}$ represents the Gilbert gyromagnetic ratio, γR=|γe|μ0≈2.21×105$\gamma ^\text{R}=|\gamma _\text{e}| \mu _0 \approx 2.21\times 10^5$ m A−1 s−1 with γe$\gamma _\text{e}$ the electron gyromagnetic ratio γe≈−1.76×1011$\gamma _\text{e} \approx -1.76\times 10^{11}$ T−1 s−1. Expressed in DPD units, the gyromagnetic factor equals γR=11.42$\gamma ^\text{R}=11.42$. The parameter η was the dimensionless Gilbert damping constant, which is set to η=0.01$\eta =0.01$ in this work.The anisotropic contribution Hiani$H_i^{\text{ani}}$ can be derived from a potential Hiani,(l)=−1μ0MsVc∂Eani∂αi(l)$H_i^{\text{ani},(l)}=-\frac{1}{\mu _0M_\text{s}V_\text{c}}\frac{\partial E^{\text{ani}}}{\partial \alpha _i^{(l)}}$43Eani=−∑lKaVc(n(l)·α(l))2$$\begin{eqnarray} E^{\text{ani}} &=& -\sum _l K^aV_\text{c}(\mathbf {n}^{(l)}\cdot {\mbox{\boldmath $\alpha $}}^{(l)})^2\end{eqnarray}$$44Hiani,(l)=2Kaμ0Ms(n(l)·α(l))n(l)$$\begin{eqnarray} H_i^{\text{ani},(l)} &=& \frac{2K^a}{\mu _0M_\text{s}}(\mathbf {n}^{(l)}\cdot {\mbox{\boldmath $\alpha $}}^{(l)})\mathbf {n}^{(l)} \end{eqnarray}$$where Ka$K^a$ is the anisotropy constant and Vc$V_\text{c}$ is the volume of the nanoparticle magnetic core. In this work Ka=223$K^a=223$ kJ m−3 was used, which corresponded to cobalt ferrite nanoparticles.[86]The dipole contribution was a summation of the magnetic fields of all the nanoparticles and their periodic images (including the periodic images of a given nanoparticle for which the force was calculated)45Hidip(r)=∑j≠l3(m(j)·r)r−r2m(j)4πr5$$\begin{equation} H_i^{\text{dip}}(\mathbf {r}) = \sum _{j\ne l}\frac{3(\mathbf {m}^{(j)}\cdot \mathbf {r})\mathbf {r}-r^2\mathbf {m}^{(j)}}{4\pi r^5} \end{equation}$$where r=r(l)−r(j)$\mathbf {r} = \mathbf {r}^{(l)} - \mathbf {r}^{(j)}$Related to the dipolar contribution is the demagnetization field of a given nanoparticle Hidem$H_i^{\text{dem}}$. In general it can be written as Hidem=DijMj$H_i^{\text{dem}}=D_{ij}M_j$ where Dij$D_{ij}$ is the demagnetization tensor. For the case of a sphere this tensor was isotropic, Dij=−13δij$D_{ij}=-\frac{1}{3}\delta _{ij}$. Since the dynamic equation conserved the magnitude of the magnetic moment, this contribution did not have an effect on the dynamics of spherical nanoparticles.The field Hith$H_i^{\text{th}}$ was a stochastic contribution[85,87,88] with zero mean ⟨Hith⟩=0$\langle H_i^{\text{th}}\rangle = 0$ and the following second moment46⟨Hith(t)Hjth(t′)⟩=2ηkBTμ0γRMsVcδijδ(t−t′)$$\begin{equation} \langle H_i^{\text{th}}(t)H_j^{\text{th}}(t^{\prime })\rangle = \frac{2\eta k_\text{B}T}{\mu _0\gamma ^\text{R}M_\text{s}V_\text{c}}\delta _{ij}\delta (t-t^{\prime }) \end{equation}$$It has been shown in ref. [89] that both Itô and Stratonovich interpretation of the stochastic term lead to the same result when applied to Equation (40). Therefore, forth order Runge–Kutta method was used to solve this equation.The last term in Equation (42) was needed to correctly account for the conservation of the total angular momentum,[21,90] J=L+S$\mathbf {J}=\mathbf {L}+\mathbf {S}$, where L=IΩ$\mathbf {L}=\textsf {I}{\mbox{\boldmath $\Omega $}}$ is the angular momentum of the nanoparticle with I$\textsf {I}$ the moment of inertia tensor and S=−μ0MsVcγRα$\mathbf {S}=-\frac{\mu _0M_\text{s}V_\text{c}}{\gamma ^\text{R}} {\mbox{\boldmath $\alpha $}}$ is the spin momentum. This contribution is also known as the Barnett field,[91,92] and is related to the Barnett effect[93] where magnetization is induced by rotation. The moment of inertia tensor was in the case isotropic and equal to Iij=35ρmVcr2δij$\textsf {I}_{ij}=\frac{3}{5}\rho _\text{m}V_\text{c}r^2 \delta _{ij}$ for a hollow sphere with radius r and nanoparticle core mass density ρm$\rho _\text{m}$.Since the particles are small and single domain, the exchange interaction term Eex≈A|∇M|2$E^{\text{ex}} \approx \ A|\nabla \mathbf {M}|^2$ with A the spin stiffness constant was discarded.The conservative part of the torque τmc${\mbox{\boldmath $\tau $}}_\text{m}^\text{c}$ acting on the nanoparticle can be calculated by inducing an infinitesimal rotation (δn=δφ×n$\delta \mathbf {n} = \delta {\mbox{\boldmath $\varphi $}}\times \mathbf {n}$). A nonzero torque appeared due to the misalignment of the easy axis and the magnetic moment47τmc=εijkδEaniδnjnk=−2KaVc(n(l)·α(l))α(l)×n(l)$$\begin{equation} {\mbox{\boldmath $\tau $}}_\text{m}^\text{c} = \epsilon _{ijk}\frac{\delta E^{\text{ani}}}{\delta n_j}n_k = -2K^aV_\text{c}(\mathbf {n}^{(l)}\cdot {\mbox{\boldmath $\alpha $}}^{(l)}) {\mbox{\boldmath $\alpha $}}^{(l)}\times \mathbf {n}^{(l)} \end{equation}$$To guarantee the conservation of the total angular momentum J$\mathbf {J}$ there is also a dissipative torque, which is oppositely equal to the dissipative torque, second term of Equation (41)48τmd=−ηMsVc1+η2α(l)×(α(l)×H)$$\begin{equation} {\mbox{\boldmath $\tau $}}_\text{m}^\text{d} = -\frac{\eta M_\text{s}V_\text{c}}{1+\eta ^2} \,\mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times (\mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times \mathbf {H}) \end{equation}$$The overall rigid body dynamics of the nanoparticle core is then given by49ẋ(l)=v(l)$$\begin{align} &\dot{\mathbf {x}}^{(l)}=\mathbf {v}^{(l)}\end{align}$$50ρmVcẍ(l)=Ftot(l)$$\begin{align} &\rho _\text{m}V_\text{c}\ddot{\mathbf {x}}^{(l)}=\mathbf {F}_{\text{tot}}^{(l)}\end{align}$$51IΩ̇(l)=τw+τmc+τmd=τw+μ0MsVcγR∂α(l)∂t+μ0MsVcα(l)×(Hext+Hdip+Hth)$$\begin{align} &\textsf {I}\dot{{\mbox{\boldmath $\Omega $}}}^{(l)} = {\mbox{\boldmath $\tau $}}_\text{w} + {\mbox{\boldmath $\tau $}}_\text{m}^c +{\mbox{\boldmath $\tau $}}_\text{m}^d \nonumber \\ &= {\mbox{\boldmath $\tau $}}_\text{w} + \frac{\mu _0M_\text{s}V_\text{c}}{\gamma ^\text{R}}\frac{\partial \mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}}{\partial t} + \mu _0M_\text{s}V_\text{c} \mathbf {{\mbox{\boldmath $\alpha $}}}^{(l)}\times (\mathbf {H}^{\text{ext}}+\mathbf {H}^{\text{dip}}+\mathbf {H}^{\text{th}}) \end{align}$$where x(l)${\mathbf {x}}^{(l)}$ is the position of the lth nanoparticle, v(l)$\mathbf {v}^{(l)}$ is its velocity, Ftot$\mathbf {F}_{\text{tot}}$ is the total force acting on the nanoparticle core, τw${\mbox{\boldmath $\tau $}}_\text{w}$ is the total torque of the solvent and other nanoparticles, analytes, and the surfactants acting on the nanoparticle surface. Since the solvent was treated explicitly and was a source of random uncorrelated kicks to the nanoparticle surface beads, adding a stochastic torque to τw${\mbox{\boldmath $\tau $}}_\text{w}$ in the spirit of rotational Brownian motion ref. [94, 95] was not needed.As a side note, it was mentioned that the rotational dynamics of the nanoparticles was characterized by low Reynolds numbers, which meant that the angular velocity of the nanoparticle could in principle be estimated from the balance of hydrodynamic and magnetic anisotropy torques, yielding52Ω(l)=−2KaVcγn(n(l)·α(l))α(l)×n(l)$$\begin{equation} {\mbox{\boldmath $\Omega $}}^{(l)} = -\frac{2K^aV_\text{c}}{\gamma ^n}(\mathbf {n}^{(l)}\cdot {\mbox{\boldmath $\alpha $}}^{(l)}) {\mbox{\boldmath $\alpha $}}^{(l)}\times \mathbf {n}^{(l)} \end{equation}$$with γn$\gamma ^n$ the rotational viscosity, which is proportional to the viscosity of the solvent. In Equation (52) it was assumed that the magnetization has already relaxed toward the direction of the magnetic field. According to the Einstein–Stokes equations, γn=8πηRh3$\gamma ^n=8\pi \eta R_\text{h}^3$. Since γn$\gamma ^n$ depended on the hydrodynamic radius, which was in principle not known beforehandthis approach will not be used. Therefore, Equation (47) was used to calculate the torque and corresponding forces on the nanoparticle, which did not require knowledge of any additional mesoscopic parameters, such as γn$\gamma ^n$.Dynamic BindingTo model the attachment of the biomarkers to the recognition parts, which were located at the end of the recognition ligands, the so called dynamic binding was employed. This is a model of the binding process, with which it was already able to be described the cluster formation of nanoparticles. In short, if the biomarker came sufficiently close (e.g., less than a critical distance dc$d_\text{c}$) to the recognition part of the recognition ligand, a harmonic bond Equation (12) was formed between the biomarker and the recognition part. dc=rc$d_c = r_c$ was set in this work. For the newly formed bond, the same parameters was used as for PEG. Due to the known high affinity between the analyte and the ligand, this bond was not allowed to break once it was formed. A single biomarker (analyte) can bind to several recognition ligands (e.g., 1, 2, 4, …) thereby modeling the analytes as multivalent. This last aspect was crucial for the cluster formation. In ref. [96] multivalent binding of DNA oligonucleotide probes attached to a surface was shown to increase the sensitivity of bacterial genome detection.Magnetic Nanoparticles and Colloidal CharacterizationThe employed magnetic nanoparticle were iron oxide nanoparticles (IONPs) coated with dextran and carboxylic PEG (product code is 104‐56‐701) synthesized by micromod Partikeltechnologie GmbH (Germany). The IONP crystal size was (30 ± 4) nm, the diffusion coefficient was (8.2×106±5.2×105)$(8.2\times 10^6 \pm 5.2\times 10^5)$ nm2 s−1, the hydrodynamic diameter was 64 nm whose PDI = 0.08. IONP size and shape of the studied IONPs were evaluated by transmission electron microscopy (TEM) (see Figure 17) in a JEM1400 Flash (Tokio, Japan) TEM operating at a 100 kV located at Servicio de Microscopía Electrónica del Centro de Biología Molecular “Severo Ochoa” CSIC‐UAM. TEM images were examined through manual analysis of more than 150 particles randomly selected in different areas of TEM micrographs using Image‐J software to obtain the mean size and size distribution. The hydrodynamic size (Dh$D_\text{h}$) of IONPs was determined by dynamic light scattering (DLS) with a Zetasizer Nano ZS equipment (Malvern Instruments, USA). IONP suspensions were diluted in double distilled water (DDW) to a final concentration of 0.05 gFe L−1 in a commercial cuvette. The energy source was a laser emitting at 633 nm, and the angle between sample and detector was 173°. For the determination of the diffusion coefficient of IONPs per unit volume in magnetic suspensions with iron contents 1 gFe L−1, Nanosight NS300 (Malvern Instruments, USA) was employed. Samples were diluted 1:5000 in DDW and injected into the instrument chamber using a 1 mL syringe. Camera settings were adjusted in order to focus the objective. The video data was collected for 60 s and repeated three times per sample.17Figurea) Representative TEM image of morphology IONPs. b) Intensity weighted hydrodynamic size distribution of IONPs dispersed in DDW.Magnetic CharacterizationSQUID measurements under quasi‐static conditions were performed in a range of temperature between 4 and 300 K in IONP suspensions with a magnetic mass concentration of 1 gFe L−1, in a commercial SQUID model MPMS‐XL from Quantum Design at Servicio de Tecnicas Físicas en la Facultad de Físicas de la Universidad Complutense de Madrid, see Figure 18.18Figurea) SQUID measurements of magnetization loops for IONPs suspension at 1 gFe L−1 at 4 K (red line) and 300 K (black line). b) Zoom in low magnetic field regions.AcknowledgementsFinancial support of the Ministry of Education, Science and Sport (MIZŠ, Slovenia, grant No. C3330‐19‐252009), and Ministerio de Ciencia e Innovación (Spain) through the M‐ERA Net BioMag project and of the Slovenian and Spanish Research Agencies (grants No. P1‐0002 and PCI2019‐103600, respectively) is gratefully acknowledged.Conflict of InterestThe authors declare no conflict of interest.Data Availability StatementThe data that support the findings of this study are available from the corresponding author upon reasonable request.H. A. Katus, A. Remppis, S. Looser, K. Hallermeier, T. Scheffold, W. Kübler, J. Mol. Cell. Cardiol. 1989, 21, 1349.S. Hosseini, P. Vázquez‐Villegas, M. Rito‐Palomares, S. O. Martinez‐Chapa, Advantages, Disadvantages and Modifications of Conventional ELISA, Springer Singapore, Singapore 2018, pp. 67–115.A. Fornara, P. Johansson, K. Petersson, S. Gustafsson, J. Qin, E. Olsson, D. Ilver, A. Krozer, M. Muhammed, C. Johansson, Nano Lett. 2008, 8, 3423.D. Cabrera, A. Aires, E. J. Artés, J. Camarero, A. L. Cortajarena, F. J. Teran, Method for detection of an analyte, WO2019092131, 2018.R. Kötitz, W. Weitschies, L. Trahms, W. Brewer, W. Semmler, J. Magn. Magn. Mater. 1999, 194, 62.H. C. Yang, S. Y. Yang, G. L. Fang, W. H. Huang, C. H. Liu, S. H. Liao, H. E. Horng, C.‐Y. Hong, J. Appl. Phys. 2006, 99, 124701.K. Wu, D. Tonini, S. Liang, R. Saha, V. K. Chugh, J.‐P. Wang, ACS Appl. Mater. Interfaces 2022, 14, 9945.H. E. Horng, S. Y. Yang, C.‐Y. Hong, C. M. Liu, P. S. Tsai, H. C. Yang, C. C. Wu, Appl. Phys. Lett. 2006, 88, 252506.K. Enpuku, D. Kuroda, A. Ohba, T. Yang, K. Yoshinaga, T. Nakahara, H. Kuma, N. Hamasaki, Jpn. J. Appl. Phys. Part 2: Lett. 2003, 42, L1436.C.‐Y. Hong, C. C. Wu, Y. C. Chiu, S. Y. Yang, H. E. Horng, H. C. Yang, Appl. Phys. Lett. 2006, 88, 212512.D. Eberbeck, F. Wiekhorst, U. Steinhoff, L. Trahms, Appl. Phys. Lett. 2009, 95, 213701.C. C. Yang, S. Y. Yang, H. H. Chen, W. L. Weng, H. E. Horng, J. J. Chieh, C. Y. Hong, H. C. Yang, J. Appl. Phys. 2012, 112, 024704.Y. Wu, Y. Yeh, C. Wu, C. Huang, Y. Chang, C. Wu, Acta Cardiol. Sin. 2020, 36, 649.D. Perera, D. P. Landau, D. M. Nicholson, G. Malcolm Stocks, M. Eisenbach, J. Yin, G. Brown, J. Appl. Phys. 2014, 115, 17D124.D. Perera, D. M. Nicholson, M. Eisenbach, G. M. Stocks, D. P. Landau, Phys. Rev. B 2017, 95, 014431.J. Tranchida, S. Plimpton, P. Thibaudeau, A. Thompson, J. Comput. Phys. 2018, 372, 406.L. Teich, D. Kappe, T. Rempel, J. Meyer, C. Schröder, A. Hütten, Sensors 2015, 15, 9251.L. Teich, C. Schröder, IEEE Trans. Magn. 2015, 51, 7209204.L. Teich, C. Schröder, Sensors 2015, 15, 28826.N. A. Usov, B. Y. Liubimov, J. Appl. Phys. 2012, 112, 023901.K. D. Usadel, C. Usadel, J. Appl. Phys. 2015, 118, 234303.K. Wu, D. Su, J. Liu, R. Saha, J.‐P. Wang, Nanotechnology 2019, 30, 502003.D. Cabrera, A. Lak, T. Yoshida, M. E. Materia, D. Ortega, F. Ludwig, P. Guardia, A. Sathya, T. Pellegrino, F. J. Teran, Nanoscale 2017, 9, 5094.M. Holz, S. R. Heil, A. Sacco, Phys. Chem. Chem. Phys. 2000, 2, 4740.R. D. Groot, P. B. Warren, J. Chem. Phys. 1997, 107, 4423.L. Néel, C. R. Hebd. Seances Acad. Sci. 1949, 228, 664.L. Bessais, L. Ben Jaffel, J. L. Dormann, Phys. Rev. B 1992, 45, 7805.W. T. Coffey, P. J. Cregg, Y. U. P. Kalmykov, 1992, 263.H. Hasimoto, J. Fluid Mech. 1959, 5, 317.B. Dünweg, K. Kremer, J. Chem. Phys. 1993, 99, 6983.J. H. Walther, M. Praprotnik, E. M. Kotsalis, P. Koumoutsakos, J. Comput. Phys. 2012, 231, 2677.P. Kundu, I. Cohen, Fluid Mechanics, Elsevier Science, 2010.A. Popadić, J. H. Walther, P. Koumoutsakos, M. Praprotnik, New J. Phys. 2014, 16, 082001.E. R. Cruz‐Chú, E. Papadopoulou, J. H. Walther, A. Popadić, G. Gengyun Li, M. Praprotnik, P. Koumoutsakos, Nat. Nanotechnol. 2017, 12, 1106.E. Papadopoulou, J. Zavadlav, R. Podgornik, M. Praprotnik, P. Koumoutsakos, ACS Nano 2021, 15, 20311.B. Padmavathi, T. Amaranath, S. Nigam, Fluid Dyn. Res. 1993, 11, 229.B. K. P. Scaife, J. Phys. D: Appl. Phys. 1986, 19, L195.N. R. Anderson, J. Davidson, D. R. Louie, D. Serantes, K. L. Livesey, Nanomaterials 2021, 11, 11.C. Junghans, M. Praprotnik, K. Kremer, Soft Matter 2008, 4, 156.J. Tate, G. Ward, Clin. Biochem. Rev. 2004, 25, 105.J. A. Wattis, Phys. D: Nonlinear Phenom. 2006, 222, 1.M. Zidar, D. Kuzman, M. Ravnik, Soft Matter 2018, 14, 6001.S. Succi, The Lattice Boltzmann Equation: For Complex States of Flowing Matter, Oxford University Press, Oxford 2018.C. S. Peskin, Acta Numer. 2002, 11, 479.P. Ahlrichs, B. Dünweg, J. Chem. Phys. 1999, 111, 8225.P. Kreissl, C. Holm, R. Weeber, Soft Matter 2021, 17, 174.L. P. Fischer, T. Peter, C. Holm, J. de Graaf, J. Chem. Phys. 2015, 143, 084107.J. de Graaf, T. Peter, L. P. Fischer, C. Holm, J. Chem. Phys. 2015, 143, 084108.A. Hubert, R. Schäfer, Magnetic Domains: The Analysis of Magnetic Microstructures, Springer, Berlin, Heidelberg 1998.D. Berkov, J. Magn. Magn. Mater. 2006, 300, 159.Z. Shaterabadi, G. Nabiyouni, M. Soleymani, Mater. Sci. Eng., C 2017, 75, 947.R. Hong, B. Feng, L. Chen, G. Liu, H. Li, Y. Zheng, D. Wei, Biochem. Eng. J. 2008, 42, 290.J. J. Thomson, Lond., Edinburgh, Dublin Philos. Mag. J. Sci. 1904, 7, 237.H. Lee, A. H. de Vries, S.‐J. Marrink, R. W. Pastor, J. Phys. Chem. B 2009, 113, 13186.Y. Li, M. Kröger, W. K. Liu, Biomaterials 2014, 35, 8467.V. Molinero, W. A. Goddard, J. Phys. Chem. B 2004, 108, 1414.P. J. Hoogerbrugge, J. M. V. A. Koelman, Europhys. Lett. 1992, 19, 155.P. Español, Phys. Rev. E 1995, 52, 1734.R. Groot, K. Rabone, Biophys. J. 2001, 81, 725.M. Venturoli, B. Smit, PhysChemComm 1999, 2, 45.W. Pan, A. Tartakovsky, Adv. Water Resour. 2013, 58, 41.P. Papež, M. Praprotnik, J. Chem. Theory Comput. 2022, 18, 1227.I. V. Pivkin, G. E. Karniadakis, J. Chem. Phys. 2006, 124, 184101.K. P. Travis, M. Bankhead, K. Good, S. L. Owens, J. Chem. Phys. 2007, 127, 014109.I. V. Pivkin, G. E. Karniadakis, Phys. Rev. Lett. 2006, 96, 206001.E. I. Barcelos, S. Khani, A. Boromand, L. F. Vieira, J. A. Lee, J. Peet, M. F. Naccache, J. Maia, Comput. Phys. Commun. 2021, 258, 107618.J. H. Hildebrand, J. Am. Chem. Soc. 1935, 57, 866.J. H. Hildebrand, Proc. Natl. Acad. Sci. U. S. A. 1979, 76, 6040.A. Jarray, V. Gerbaud, M. Hemati, Powder Technol. 2016, 291, 244.J. Gupta, C. Nunes, S. Vyas, S. Jonnalagadda, J. Phys. Chem. B 2011, 115, 2014.R. Delgado‐Buscalioni, J. Sablić, M. Praprotnik, Eur. Phys. J.: Spec. Top. 2015, 224, 2331.J. Sablić, M. Praprotnik, R. Delgado‐Buscalioni, Soft Matter 2016, 12, 2416.J. Sablić, M. Praprotnik, R. Delgado‐Buscalioni, Soft Matter 2017, 13, 4971.R. Delgado‐Buscalioni, P. V. Coveney, J. Chem. Phys. 2003, 119, 978.G. De Fabritiis, R. Delgado‐Buscalioni, P. V. Coveney, J. Chem. Phys. 2004, 121, 12139.M. Praprotnik, L. Delle Site, K. Kremer, J. Chem. Phys. 2005, 123, 224106.M. Praprotnik, L. Delle Site, K. Kremer, Phys. Rev. E 2006, 73, 066701.A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, P. J. in 't Veld, A. Kohlmeyer, S. G. Moore, T. D. Nguyen, R. Shan, M. J. Stevens, J. Tranchida, C. Trott, S. J. Plimpton, Comp. Phys. Commun. 2022, 271, 108171.M. E. Tuckerman, Statistical Mechanics: Theory and Molecular Simulation, Oxford University Press, New York 2010.J. J. Cerdó, V. Ballenegger, O. Lenz, C. Holm, J. Chem. Phys. 2008, 129, 234104.S. J. Kemp, R. M. Ferguson, A. P. Khandhar, K. M. Krishnan, RSC Adv. 2016, 6, 77452.D. Fincham, Mol. Simul. 1992, 8, 165.L. D. Landau, E. Lifshitz, Phys. Z. Sowjet. 1935, 8, 153.T. L. Gilbert, IEEE Trans. Magn. 2004, 40, 3443.W. F. Brown, Phys. Rev. 1963, 130, 1677.A. J. Rondinone, A. C. S. Samia, Z. J. Zhang, Appl. Phys. Lett. 2000, 76, 3624.J. L. García‐Palacios, F. J. Lázaro, Phys. Rev. B 1998, 58, 14937.W. T. Coffey, Y. P. Kalmykov, J. Appl. Phys. 2012, 112, 121301.D. V. Berkov, N. L. Gorn, J. Phys.: Condens. Matter 2002, 14, L281.N. Usov, B. Y. Liubimov, J. Magn. Magn. Mater. 2015, 385, 339.H. Keshtgar, S. Streib, A. Kamra, Y. M. Blanter, G. E. W. Bauer, Phys. Rev. B 2017, 95, 4447.T. V. Lyutyy, S. I. Denisov, P. Hänggi, Phys. Rev. B 2019, 100, 134403.S. J. Barnett, Phys. Rev. 1915, 6, 239.D. B. Reeves, J. B. Weaver, Appl. Phys. Lett. 2015, 107, 223106.P. Ilg, Phys. Rev. B 2017, 95, 214427.T. Curk, C. A. Brackley, J. D. Farrell, Z. Xing, D. Joshi, S. Direito, U. Bren, S. Angioletti‐Uberti, J. Dobnikar, E. Eiser, D. Frenkel, R. J. Allen, Proc. Natl. Acad. Sci. U. S. A. 2020, 117, 8719.
Journal
Advanced Theory and Simulations – Wiley
Published: May 1, 2023
Keywords: biosensors; clustering; magnetic dynamics; mesoscopic simulation; nanomedicine; nanoparticles