{'fraction_non_alphanumeric': 0.03321178120617111, 'fraction_numerical': 0.012566619915848527, 'mean_word_length': 4.637571157495256, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Modelling Meso-Scale Diffusion Processes in Stochastic Fluid Bio-Membranes

The space-time dynamics of rigid inhomogeneities (inclusions) free to move in a randomly fluctuating fluid bio-membrane is derived and numerically simulated as a function of the membrane shape changes. Both vertically placed (embedded) inclusions and horizontally placed (surface) inclusions are considered. The energetics of the membrane, as a two-dimensional (2D) meso-scale continuum sheet, is described by the Canham-Helfrich Hamiltonian, with the membrane height function treated as a stochastic process. The diffusion parameter of this process acts as the link coupling the membrane shape fluctuations to the kinematics of the inclusions. The latter is described via Ito stochastic differential equation. In addition to stochastic forces, the inclusions also experience membrane-induced deterministic forces. Our aim is to simulate the diffusion-driven aggregation of inclusions and show how the external inclusions arrive at the sites of the embedded inclusions. The model has potential use in such emerging fields as designing a targeted drug delivery system.

Amphiphilic molecules, such as lipids and proteins, can self-assemble themselves into a variety of exotic structures in an aqueous environment [1]. Computer-based simulation of the dynamics of these structures forms an interesting research area in the computational statistical mechanics of bio-material systems. One such structure is the phospholipid bilayer which represents the generic structure of all bio-membrane systems, both natural and artificial. These membranes can have thicknesses of only few nano-metres but linear sizes of up to tens of micro-metres, and can therefore be regarded as highly flexible, fluid-like, 2D continuum sheets embedded in a three-dimensional space.
Thermal fluctuations can induce shape fluctuations and shape transformations in membranes. For example, in the so-called budding transition [2] a spherical vesicle transforms into a prolate ellipsoid as the temperature increases. There is also the possibility that the spherical geometry becomes oblate, producing a shape similar to the biconcave rest shape of a red blood cell.
Bio-membranes regulate the recognition and transport processes to and from the interior of the cells, as well as between the cells, forming a barrier which all external particles arriving at the cell must cross. They contain a variety of integral (embedded) inhomogeneities, such as proteins and other macromolecules [3], that penetrate the thickness of the membrane and act as transport channels. We shall refer to these as the M-type inclusions. These inclusions are mobile and can freely diffuse across the membrane. Their presence, however, force the bilayer to adjust its thickness locally so as to match the thickness of the hydrophobic region of these inclusions [4,5], causing local deformations in the membrane geometry. The perturbations produced in the membrane shape due to these local deformations give rise to both short-range and long-range membrane-induced indirect forces between the inclusions. These forces act in addition to direct Van der Waals and electrostatic forces between the inclusions. The Long-range forces originate from the perturbations associated with the long wavelength shape fluctuations [6], whereas the short-range forces are associated with the local deformations in the immediate vicinity of the inclusions [4]. These membrane-induced forces between the inclusions play a far more significant role, viz the direct molecular interactions, when the length scales involved are comparable to the size of the membrane. In addition to this mode of deformation, a membrane can also deform due to the tension at the amphiphilic molecules-water interfaces.
This tension results in a change in the overall surface area of the membrane. A third mode of deformation also exists and this is associated with the bending elastic-curvature property of the membrane which distinguishes it from a sheet of simple fluid dominated by surface tension [2]. Accordingly, two models to study the inclusion-induced local deformations have been developed. In the first model, the membrane energy is taken to consist of a contribution from the molecular expansion/compression due to the change in the thickness at the inclusion boundary, and also a contribution from the overall change in the surface area. Using this model, it is shown that [7,8,9,10,11] the inclusion-induced deformations cause exponential decays in the thickness of the membrane, extending from the inclusion-imposed value to the equilibrium thickness value, as shown schematically in Fig.1 for two rod-like inclusions. In the second model [4,12], the contribution of the membrane bending property is taken into account in the energy term, and it is found that this significantly affects the perturbation profile at the inclusion boundary as well as modifying the membrane-induced interactions. Evidently, an object supported by surface tension would have a different dynamics than one supported by the bending elasticity of the surface [13].
In addition to the M-type inclusions, membranes can also carry inclusions that lie on their surfaces [14] as shown schematically in Fig.2. These surface inclusions can represent objects that have arrived at the membrane from the outside and are therefore referred to as external inclusions. We shall refer to these as the S-type inclusions.
At the meso-scale, i.e. when the detailed molecular architecture of the membrane can be subsumed into a background 2D sheet, the free elastic energy of a symmetric membrane is described by the Canham-Helfrich Hamiltonian [15,16] where are respectively the Gaussian and mean curvatures of the sheet, R 1 and R 2 are the two principle radii of curvature of the sheet, σ 0 is the surface tension, κ is the bending rigidity,κ is the Gaussian rigidity, g is determinant of the metric tensor and σ = (σ 1 , σ 2 ) is the 2D local coordinate on the sheet as opposed to the coordinates on the embedding space. The last term in (1) is, by Gauss-Bonett theorem, an invariant for closed surfaces implying that the dynamics of a membrane is not influenced by this term if its topology remains fixed. In what follows, we concentrate on membranes with fixed topology and drop this term. We then have The study of a membrane whose free energy is described by (3) is facilitated by considering it to be nearly flat, i.e. its thickness to be much smaller than its linear size L. This is indeed what we mean by a meso-scale model of a membrane. We therefore take the membrane to be almost parallel to the (x 1 , x 2 ) plane, regarded as the reference plane. The position of a point on the membrane can then be described by a single-valued function h(x 1 , x 2 ) representing the height of that point. This simplification is achieved by writing the Hamiltonian (3) in the Monge representation [17] which gives for the mean curvature where ∂ i ≡ ∂ ∂x i , i = 1, 2. We assume that the area of the membrane can fluctuate without constraint by setting σ 0 = 0 in (3). Consequently, using (4), the Hamiltonian (3) to leading order in derivatives of h becomes This is the Canham-Helfrich Hamiltonian in Monge representation, expressed in terms of the height function of the membrane. It is the expression that we employ to describe the energetics of the membrane.
Employing a statistical mechanics based on (5) only, and ignoring the contributions from the expansion/compression and interfacial energies, the potential energy function was constructed [6] to describe the membrane-induced temperature-dependent long-range forces between a pair of disk shape M-type inclusions that can freely tilt with respect to each other. Another function was also constructed for long-range interaction between two S-type where A = πr 2 0 is the area of an M-type inclusion of radius r 0 , k B is the Boltzmann constant, R ij is distance between the centres of mass of two inclusions i and j, L i and L j are the lengths of two S-type inclusions making the angles θ i and θ j respectively with the line joining their centres of mass (see Fig.2) and T is the membrane temperature. It is evident that both of these membrane-induced potentials are attractive and fall off as R −4 with the distance. These expressions are derived for rod-like inclusions that are assumed to be much more rigid than the ambient membrane so that these inclusions can not move coherently with the membrane.
The only degrees of freedom for the rods are rigid translations and rotations while they remain attached to the membrane.
So far, the modelling of bio-membrane dynamics decorated with inclusions has been mainly concerned with constructing potential energy functions such as those given in (6) and (7).
An interesting problem, however, would be to use this information to simulate the spacetime behaviour of inclusions in a membrane described by (5) and undergoing stochastic shape fluctuations. This is the problem that we address in this paper. This type of simulation can establish a direct link between the randomly changing membrane shape on the one hand and the inclusion dynamics on the other. In such a simulation, the thermodynamic phase behaviour of inclusions, such as their temperature-dependent aggregation, can be directly computed.
This phase behaviour plays a crucial role in the functional specialisation of a membrane [4]. To proceed, let us consider a 2D bio-membrane described by (5) containing both the Mtype and the S-type inclusions. To make the membrane a stochastically fluctuating medium, we treat the height function in (5) as a stochastic Wiener process with a Gaussian distribution, whose mean and variance can be written as [19] h( where D is the diffusion constant associated with the height fluctuations at the local position (x 1 , x 2 ) and represents the measure with which random fluctuations propagates in the local geometry. Such random height changes would cause a roughening of the membrane surface on molecular scales, and this has been observed in NMR experiments [20].
Assuming that this is the only stochastic process present in the membrane, it is then reasonable to suppose that this stochastic dynamics is communicated to the inclusions as well, and that their ensuing random motions are contingent only on these fluctuations. This implies that the mathematical point representing the centre of mass of an inclusion coinciding with the membrane point (x 1 , x 2 ) would also experience the same fluctuations and would diffuse with the same diffusion constant. To derive an expression for D, based on (5), we start with the static height-height correlation function obtained from (5). This is given by [2,21] h(q; where · · · is the thermal averaging with respect to the Boltzmann weight, exp(−H 0 /k B T ), and q is the wave vector of magnitude q. The corresponding dynamic correlation function can be obtained [2] by writing giving where the damping rate, γ 0 (q), reflecting the long-range character of the hydrodynamic damping, is defined as and η denotes the coefficient of viscosity of the fluid membrane. In real space, (12) transforms [17,21] to where L is the length of the membrane. This is the equal-time correlation function for membrane fluctuations. A similar model of an active fluctuating membrane in which the vertical displacements of the membrane satisfy a Langevin equation in the q space has also been proposed [22], and is its shown that a term similar to the static version of (12) contributes to the correlation function which also contains a contribution from non-equilibrium fluctuations.
The latter is in the form of a q −5 term which dominates at long distances. Comparison of (9) and (14) yields the desired result It should be emphasised that the association of a diffusive process with the membrane height function, and the resulting diffusion constant, is not analogous to the usual model of a diffusion process in which, for example, a particle diffuses through a medium, such as a fluid. Rather, what is suggested here is that the magnitude of a mathematical function representing the height of a mathematical point in the membrane is subject to random stochastic variations, and the diffusion constant is a measure of this variation.
When the M-type inclusions are present they produce exponentially decaying local deformations in the membrane geometry (see Fig.1). Correspondingly, the correlation function can be modified by a multiplicative exponential factor to where r 0 is the radius of an M-type inclusion, R + r 0 is the radius of the circular region around the inclusion with its centre coinciding with that of the inclusion, and NI stands for near inclusion. It is evident that outside this region the exponential decay of the profile is negligible. Accordingly, within this circular region of radius R + r 0 , the diffusion constant is also modified to This equation implies that when the centre of mass an S-type inclusion enters a circular region of radius R + r 0 its diffusion coefficient goes over to D m and progressively approaches zero as the boundary of an M-type inclusion is approached. We can ascertain, as a first approximation, that this is how an M-type inclusion interacts with an S-type inclusion.
In our simulations, the equations motion of both the S-type and the M-type inclusions are represented by the differential equation of the the Ito stochastic calculus [19] dr(t) = A[r(t), t] dt + D 1/2 dW(t).
This equation describes the stochastic trajectory, r(t), of the centres of mass of the inclusions in terms of a dynamical variable of the inclusions, A[r(t), t], which is referred to as the drift velocity, and a term, dW(t), which is a given Gaussian Wiener process with the mean and variance given by Equation (18)  we can re-write (18) as where F(t) is the instantaneous systematic force experienced by the i-th inclusion and is obtained from the inter-inclusion potentials, given in (6) and (7), according to We implemented (20) for our 2D simulations according to the iterative scheme [24] where R G X and R G Y are standard random Gaussian variables chosen separately and independently for each inclusion according to the procedure given in [23], and F X , F Y are the X and Y components of the force F. For the S-type inclusions, we treated the angles in (7) as independent stochastic variables described by where τ is the torque experienced by an S-type inclusion and is given by and θ G is the angular counterpart of R G X and R G Y .
In the numerical simulations, recently reported in their broad outlines [25], we use a square membrane with L = 40µm on its side. The other parameters used were set at κ = 10 −19 J and η = 10 −3 J sec m −3 [2]. These values correspond to the condition in which the bending mode of the membrane is important. From these data the damping coefficient, γ 0 , in the real space of the membrane, can be obtained from (13). The simulation temperature was set at T = 300 • K, and the correlation (delay) time, t, over which the diffusion coefficient in (15) was calculated, was set at t = 10 −4 sec. These data gave D = 2.6 × 10 −9 m 2 sec −1 , which is in close agreement with the value of D ≈ 4.4 × 10 −9 m 2 sec −1 obtained at the molecular level via an MD simulation of a fully hydrated phospholipid dipalmitoylphosphatidylcholine (DPPC) bilayer diffusing in the z-direction [26]. To justify our choice of the correlation time, t = 10 −4 sec, we recall that the time scale of a stochastic particle, t D , of mass m is usually determined from the relation [27] Since t D is normally of the order of 10 −9 sec, then for the criterion of long-time dynamics, employed in our model (cf (20)), to be justified the correlation (diffusive) time scale, t, in (15) has to satisfy the condition [27] t ≫ t D .
For our calculated value of D and our choice of the inclusion mass m = 10 −12 µg, corresponding to an inclusion of length L i = 0.1µm, we obtained a value of t D = 0.6 × 10 −9 sec, showing that our choice of the correlation time was appropriate to satisfy the condition in (26). The radius of an M-type inclusion was set at r 0 = 0.01µm, and the inclusions were all equal in length.
The stochastic trajectories of the inclusions were obtained in a set of five simulations. The simulation time step, dt, in (22) was set at dt = 10 −9 sec, and each simulation was performed for 4 × 10 6 time steps, i.e for a mesoscopic interval of 4000µsec. The total number of inclusions considered was 36, consisting of 13 S-type and 23 M-type.
In the first simulation, we computed the random motions of the S-type inclusions in a membrane devoid of the M-type inclusions. This was done in order to observe the details of the drift-diffusion motion over mesoscopic time scales. Figure 3 shows To examine the membrane response to temperature changes, two of the simulations were performed at different temperatures. Figure 5 shows the snapshots of these simulations at T = 100 • K (a to c) and at T = 350 • K (d to f). Figures 5d to 5f clearly show that both the number of captured inclusions and the aggregation of the M-type inclusions were affected by these temperature differences, as can be seen by comparing figures 5c and 5f.
To sum up, although many dynamical aspects of membrane-like surfaces have been addressed in the past [28], it is only relatively recently that attention has focused on the dynamics of membranes with inclusions. To our knowledge no computer-based simulation of this dynamics has been reported so far. In this paper we constructed a meso-scale model of a generic bio-membrane based on the Canham-Helfrich curvature-energy formalism. We treated the height function of the membrane as a stochastic Wiener process whose correlation function