# Effect of shape anisotropy on transport in a 2-d computational model: Numerical simulations showing experimental features observed in biomembranes

## 1 Introduction

The problem of packing of spheres plays a major role in the modeling of many physical systems and has been studied for more than four decades. Some of the early examples  of the computer simulations of hard sphere liquids suggest the existence of a first order freezing transition. The problem of packing of spheres in two and three dimensions is of great interest. Recent investigations of such systems have focused on the study of the statistical geometry of the dense sphere packing. Such studies are important in the understanding of physical properties of many systems, composed of a large number of particles .

In this context we pose a question, with the motivation of studying the transport across a two dimensional structure of packed circular disks (membrane), how does the packing change when the membrane is doped with objects of various shapes and sizes (e.g. spheres arranged rigidly in the form of rods of different lengths, L, T, X shapes etc. See Fig. 1) ? In particular we investigate the effect of these shapes on the distribution of “voids”. The “anisotropy” in the interaction potential appears to play a key role in the induction of large voids.

As pointed out by Sastri et. al. , no algorithm is available to compute void statistics for the packing of shapes other than spheres. In this paper we propose a simple numerical algorithm to compute void statistics. Unlike a probabilistic algorithm (Monte Carlo), our algorithm is based on digitization and cell counting.

The paper is organized as follows. In Sec. 2, we describe the model system. A definition of “void” and an algorithm to compute void statistics is given in Sec. 3. The results of numerical simulations and their relevance in lipid biomembranes is discussed in Sec. 4 We summarize the paper in Sec. 5.

## 2 The model system

The configuration space of the model system (membrane) is considered as a two dimensional space with periodic (toroidal) boundary conditions. The constituents of the membrane are disks and dopants.

### 2.1  The basic model

We consider a membrane made up of only circular disks interacting pairwise via the Lennard-Jones potential:

$`V_{LJ}(r_{ij})=4ϵ{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \underset{j=i+1}{\overset{N}{}}}\left(({\displaystyle \frac{\sigma }{r_{ij}}})^{12}({\displaystyle \frac{\sigma }{r_{ij}}})^6\right)`$   
where, $`r_{ij}`$ is the distance between the centers of the $`i^{\mathrm{th}}`$ and $`j^{\mathrm{th}}`$ disks, $`\sigma `$ determines the range of hard core part in the potential and $`ϵ`$ signifies the depth of the attractive part. We choose the number of disks such that the area occupied by these disks is around $`70\%`$, which is less than that of the close-packed structure but still large enough to produce some closed voids.

### 2.2  The model with impurities

Further, we consider different shape anisotropic combinations (dopants) consisting of $`\kappa `$ number of circular disks. We treat each of these combinations as a single rigid cluster. Several such dopants (impurities) are considered. Fig. 1 shows some of these impurities. The interaction between impurities and disks or other impurities is obtained by superposing the Lennard-Jones potential corresponding to each of the constituent disk in impurity. We consider a membrane with circular disks and impurities amounting to $`10\%`$ of the total number of circular disks, such that the area occupied is still $`70\%`$.

These membranes are brought to an equilibrium configuration by the Monte Carlo method  at a fixed temperature. Fig. 2 and Fig. 3 show typical equilibrium configurations of membrane without and with impurities respectively (The impurity in Fig. 3 is a rod shaped structure made up of five disks (Rod<sub>5</sub>), in general Rod<sub>κ</sub> for rod made up of $`\kappa `$ number of disks). In the simulation the temperature is so chosen that $`k_BT<4ϵ`$, where $`k_B`$ is the Boltzmann constant. The equilibrium is confirmed by simulated annealing.

## 3 Voids and an algorithm for void statistics

Now, we introduce the notion of an “$`r`$-void” in a membrane which is suitable for the description of transport across membrane and further, propose an algorithm to compute statistical quantities such as the number of voids in the membrane, the void size distribution etc.

We define an $`r`$-void as a closed area in a membrane devoid of disks or impurities, and big enough to accommodate a circular disk of radius $`r`$. Of course an $`r`$-void is also an $`r^{}`$-void if $`r^{}<r`$.

### 3.1 The algorithm to compute void statistics

To compute the void statistics for $`r`$-voids, we increase the radii of the disks forming the membrane (including the disks in the impurities, without altering the positions of the centers) by an amount $`r`$ (See Fig. 4). Then we digitize the entire membrane on a suitably chosen grid. The choice of grid size depends on the required accuracy and the typical sizes of the voids. The digitization of circular disks is carried out by the Bressenham circle drawing algorithm , modified to incorporate periodic boundary conditions. The number of voids in the membrane are computed by flood filling  every closed void with a different color and then counting the number of colors. The sizes of various voids can be obtained by counting the number of grid-cells filled by the corresponding color. The termination of flood fill algorithm is ensured since the voids are closed. In our case this condition is automatically fulfilled in view of periodic boundary conditions.

The geometric algorithms involving Vorenoi polygons  are mathematically satisfying and are expected to be accurate but would take much more computation time. On the other hand, as pointed in , the probabilistic algorithm is time efficient but requires a very large sample size while dealing with small voids.

Our algorithm is quite efficient as well as suitable even when there are small voids in the membrane. We further note that the algorithm can be easily generalized to higher dimensions. We expect that the efficiency of this algorithm can be further enhanced by the use of a multi-resolution adaptive grid.

## 4 Results and Discussions

The simulations were carried out for membranes of different compositions. Fig. 5 shows the graphs of the number of $`r`$-voids as a function of $`r`$ measured in units of the radius of the constituent disks. Curve (a) shows void distribution in absence of impurities. Curve (b) represents the void distribution in a membrane with rod shaped impurities made up of two disks (Rod<sub>2</sub>). Curves (c) and (d) show the void distribution with L shaped impurities made up of four disks (L<sub>4</sub>) and rod like impurities made up of four disks (Rod<sub>4</sub>) respectively. It is clear from the graph that the number of large voids increases with an increase in the anisotropy of the impurity. Even though L<sub>4</sub> and Rod<sub>4</sub> occupy the same area, Rod<sub>4</sub> being more anisotropic induces a larger number of big voids than L<sub>4</sub>. This fact can be clearly seen in Fig. 5, curves (c) and (d). Moreover, the Fig. 2 and Fig. 3 demonstrate the fact that the voids are mostly found in the neighborhood of the centers of anisotropy. Further, to strengthen our claim that the shape anisotropy induces voids, we compared two membranes. In one case we added rod impurities made up of two disks (Rod<sub>2</sub>) in the assembly of circular disks, and in the other case we added circular impurities of larger size, which occupied the same area as that of Rod<sub>2</sub>. We found that the former, being more anisotropic, induced larger and more numerous voids as compared to the later, though they occupied the same area.

Thus, reduced to the bare essentials, the anisotropy in the interaction potential of the constituents, is seen to be responsible for the induction of large voids. If studied from the perspective of energy minimization, as the potential becomes direction dependent, some positions of the constituents are preferred over the other positions. This induces large voids.

These features show a remarkable similarity with the observations reported in certain biological experiments . These experiments deal with the size-dependent permeation of non-electrolytes across biological membranes. The effect of doping on the permeation of large molecules was studied in these experiments. The liposome-membrane used in these experiments was made up of mixture of two types of lipids (cardiolipins and phosphatidylcholine) in a proportion 1:10. The understanding of the enhancement of transport in doped membranes needed an algorithmic statement. The ingredients at the algorithmic level involved:

1. consideration of the structure as a strictly 2–dimensional assembly
2. the cross sections of molecules being considered as constituents
3. interactions of the constituents via the Lennard Jones potential
4. permeating particles being considered as hard disks.

The features reported in  bear a similarity with the simulation carried out with Rod<sub>2</sub> as dopants. We have already seen in numerical simulations (See Fig. 5, curves (a) and (b)) that the Rod<sub>2</sub> type of impurities induced large voids in the membrane. The appearance of larger voids naturally enhances the transport of large particles. Thus an enhancement in the transport of large non-electrolytes like glucose, which was observed in the lipid mixture  can possibly be understood using our simple approach.

Further, apart from the biological implications, the model discussed is general enough to incorporate the studies of transport in various weakly bound granular media.

## 5 Summary

We have presented a numerical algorithm to compute the entire void statistics in a two dimensional membrane consisting of circular disks and dopants. We found that our simple two dimensional model has shown results consistent with features observed in a complex biological system. The biological justification of the model and implications are discussed elsewhere . Nevertheless, our model and the proposed numerical algorithm which finds out the void statistics in the model system are quite general and use no specific features of any particular system. Therefore it is possible to use this method effectively in various systems from diverse disciplines. The result that the shape anisotropy induces large voids in mixtures may be used as a tool for achieving controlled selective permeability across such a system by merely changing the shape of the constituents of the mixture.

Acknowledgments We thank N.V. Joshi, Deepak Dhar, H. E. Stanley and S.S. Manna for fruitful discussions.

Figure captions:

Some examples of the impurities.

Rod type impurity made up of three circles (Rod<sub>3</sub>).

L type impurity made up of four circles (L<sub>4</sub>).

X type impurity made up of five circles (X<sub>5</sub>).

T type impurity made up of five circles (T<sub>5</sub>).

Typical equilibrium configuration of a membrane without impurity. There are 556 circular disks used to form this membrane. The number is so chosen that the area occupied is $`70\%`$. The $`\sigma `$ in Lennard-Jones potential is chosen as two times the radius of a circular disk.

Typical equilibrium configuration of a membrane with impurity of type Rod<sub>5</sub>. The amount of impurity is 1:10 proportion. All the other parameters are same as Fig. 2.

Figure describes the algorithm to compute void statistics. The radius of a circular disks (black disks) is $`R`$. These disks are expanded by amount $`r`$, so that the region $`V`$ is the void for particle of size $`r`$.

The graphs of number of $`r`$-voids as a function of $`r`$ measured in units of the radius of the constituents.

The void distribution without impurities.

The void distribution with impurity of type Rod<sub>2</sub>.

The void distribution with impurity of type L<sub>4</sub>.

The void distribution with impurity of type Rod<sub>4</sub>.

Typically 10000 Monte Carlo steps are thrown away as thermalisation, and it is ensured that total energy is minimized. The curves are averaged over 100 Monte Carlo steps.