Monte Carlo Simulation of Static and Dynamic Properties of Linear Polymer in a Crowded Environment

In this paper, we investigate the static and dynamic properties of linear polymer in the presence of obstacles. A Monte Carlo (MC) simulation method in two dimensions with a bond ﬂ uctuation model (BFM) was used to achieve this goal. To overcome the entropic barrier, we put the middle monomer of the polymer in the middle of the pore, which is placed between ordered and disordered obstacles. We probed the static properties of the polymer by calculating the mean square of the radius of gyration and the mean square end-to-end distance of the polymer, and we found that the scaling exponents of both the mean square end-to-end distance h R 2 i and the mean square radius of gyration h R 2 g i as a function of the polymer length N vary with the area fraction of crowding agents, ϕ . The dynamic properties have also been studied by exploring the translocation of the polymer. Our current research shows that the escape time τ increases as ϕ increases. Moreover, in the power-law relation of escape time τ as a function of polymer length N , the scaling exponent ( α ) changes with ϕ . Furthermore, the study has shown that the translocation of the polymer favors the disordered barriers.


Introduction
In polymer physics, the static and dynamic properties of polymers are critical issues. The translocation of biopolymers is primarily focused to cope with the features. Biopolymer translocation is the passage of the polymers through a small aperture. It involves many biological processes, such as the viral injection of linear and cyclic DNA into a host cell, the passage of protein through the channel of the membrane, and the transportation of mRNA through nuclear pores [1]. This implies that, at the scale of molecules, life is controlled by biopolymers: RNA, DNA, and proteins are basic features of biological structure and function [2,3]. Furthermore, the passage of biopolymers through nanopore has a technological applications in rapid sequencing of DNA [4,5], therapy of gene [6], and controlled delivery of drug [7,8].
However, the practical applications of these biopolymers are affected by their physical properties. Mostly, to deal with these properties of polymers, their particular classes are focused. Linear, the polymer architecture to which we are devoted, forms a very interesting class as it has two ends. Hence, it has got our attention.
Even though a few studies have been done with simulation [9][10][11][12][13][14] of this class of biopolymers' translocation through nanopore, the deep level theoretical and simulation investigations of static and dynamic properties in a crowded environment are still underway. For example, the effect of obstacle density on polymer diffusion coefficients has not been studied in previous studies [15,16] nor have the effects of fixed-size obstacles on properties [17]. As it is known, the numerous intracellular and extracellular biological environments are crowded by cells and large molecules densely [18,19]; hence, we believed that it is very impressive to study the properties of polymers in the presence of obstacles (crowding agents).
In the presence of obstacles, Gopinathan and Kim [15] and Cao et al. [16] studied the translocation of polymers without applying external field, which is driven by concentration difference of the obstacles. This means that they studied the translocation of polymers from crowded cis of obstacles density ϕ c > 0 to sparse or uncrowded trans −ϕ t = 0, without applying potential (Δμ = 0). However, the unforced translocation of polymer from sparse or uncrowded cis-(ϕ c = 0) to crowded trans-(ϕ t > 0) is difficult. And, very recently, Gao et al. [20] studied the translocation of linear and ring polymers in a crowded environment of the same density of obstacles at both cis-and trans-sides, without applying electric potential, which is driven by the free energy difference. In our study, we considered a crowded environment of the same concentrations (ϕ c = ϕ t ) and fixed size of obstacles at both sides of the membrane, and we studied the unforced translocation of linear polymers in a crowded environment, which is solely caused by thermal fluctuation. In all simulations, we contemplated excluded volume obstacles and polymers (no monomer-obstacle interaction).
Previously, in order to study polymer properties, numerous simulation methods were used, some of which are often used to simulate polymers, such as the molecular dynamics (MD) [21,22] and Monte Carlo (MC) [23,24] methods. The MC simulation method, which we employed in this study, involves generating and accepting or rejecting possible conformations (states) stochastically [25,26]. To examine the parameters of the static and dynamic properties of the polymer by this method, we employed the bond fluctuation model (BFM).
In this study, we investigate how the density of the crowding agents (ϕ) affects the static and dynamic properties of linear biopolymer translocation from one side of membrane to the other in two dimensions (2D). Moreover, we examined the effects of ϕ on the scaling exponent of mean square end-to-end distance (hR 2 i), mean square radius of gyration (hR g 2 i), and escape time (τ) with the chain length of the polymer N. In addition, we investigated how ϕ influences the probability distribution of escape time and the diffusion of the polymer.

Model and Methods
As we mentioned above, in the simulation method we used, the acceptance or rejection of possible conformations of the polymer is a stochastic process. Therefore, we will require a method for carrying out this task. The bond fluctuation model (BFM) is used in our simulations to complete the objective. For coarse-grained polymer chains, it is an effective lattice Monte Carlo (MC) technique, where each monomer inhabits a certain number of lattice sites exclusively on a square lattice cell [27]. The model's implementation is as follows.
2.1. Bond Fluctuation Model (BFM). BFM, which was proposed by Carmesin and Kremer [27], has been extensively used to investigate the properties of a several variety of polymer systems. The polymer chain is modeled as a string of beads (monomers) placed on a lattice where there is a link (bond) between each monomer. Each of the lattice sites can only accept one monomer to ensure the self-avoiding walk (SAW) restrictions. Each of the beads (monomers) resides at the four vertex points of a square site on the lattice. The monomers were then joined by a predetermined number of bond vectors to their nearest neighbors. The neighbor monomers on a chain must be within a bond distance set, which are in the range of 2 ≤ b ≤ ffiffiffiffiffi 13 p , where b denotes the distance between two consecutive beads. Despite the bond lengths are allowed to fluctuate, they have to be among the set of lengths 2, ffiffi ffi , and ffiffiffiffiffi 13 p . All distances are expressed in lattice size units. The excluded volume effect is guaranteed by the least distance 2, and the maximum bond length, ffiffiffiffiffi 13 p , keeps off the bonds from crossing each other. For a single monomer, there are four possible lattice directions to move, one is chosen randomly with equal probability, the move is accepted if the trial move does not violate the two restrictions (self-avoiding and excluded volume restrictions), and the position is not occupied by an obstacle bead. Hence, throughout the simulations, the local movement of the randomly selected monomer and the attempted displacement in the randomly chosen lattice direction are governed by BFM.
2.2. Simulation Procedure. As we use a lattice for our simulation purposes, the first thing that we do is prepare the two dimensional L × L square lattice. Thus, we prepared a 2D simulation box of 400 × 400 square lattice cells of unit length each with an impenetrable membrane wall in the middle, which was wide enough for the simulation of the polymers considered in our study. And then, we added a pore (width = 6) at the middle of the wall, which is small enough to let only one monomer or two monomers pass through it. The immobile, impenetrable, and fixed size obstacle beads are put on the square lattice with the same area fraction of crowding agents (ϕ) on both the cis-(left) and trans-(right) sides of the membrane. The fraction of surface (volume) occupied by crowding agents is estimated to be 20%-30% (0.2-0.3) of the total surface (volume) or higher [28] .We simulate the polymer, in different values of area fractions of the obstacles, which are ϕ ≤ 0:3 [29] and defined as where N c is number of crowding agents, and A is the total lattice points covered by the obstacles; the factor 4 is due to the fact that each obstacles resides at 4 vertex points of the square lattice. ϕ is a measure of the density of obstacles, in a sense that there is the concentration of the obstacles in the polymer solution. On the trans-side, we placed the disordered obstacles, whereas the ordered obstacles are placed on the cis-side.
To deal with the impact of self-avoidance on the passage of polymers through a pore in the coiled state, simulations of 2 Advances in Polymer Technology dimensions greater than one (d > 1) are necessary. 2D polymers are ideally fitted to this aspiration for two reasons that excluded volume effects which are more ostensive, and shorter computation time is required, comparing with that of three dimensional case. Our simulation is then as follows. To obtain the polymer's initial configuration, we insert the polymer symmetrically at the center of the simulation box, as shown in Figure 1.
The number of monomers (chain length) N is randomly selected and the middle monomer placed in the aperture symmetrically to overcome the entropic barrier. And the excluded volume property is acquainted by considering that each lattice site be allied to one square plainly. The lengths of the links between neighbor monomers are governed by BFM on 2D.
Starting from the initial systematic placement of the polymer, many moves are attained until its equilibration is reached. To attain such an equilibrated configuration, the chain is relaxed by making many local moves and pinching the middle monomer at the center of the pore. The relaxation time or equilibration time is τ equi~N 1+2ν , where ν = 0:75 for 2D. This means that, τ relax = cN 1+2ν , c is a constant given in MC step time unit, as the unit of time in this study is Monte Carlo step (MCS), defining one MC step as N elementary moves. It is counted when the N monomers are made a local move without violating the self-avoiding and excluded volume restrictions. In our simulations, each monomer has relaxed for 10 6 MCS (c = 10 6 MCS). The reason why we used this much MC steps for relaxation was to get the fully relaxed polymer. At t = 0, the middle monomer is permitted to move freely after the equilibration is complete. The relaxed polymer is depicted in Figure 2.
Thus, the end of the simulation is at a time t > 0 when the whole monomers of the polymer are on either side (cis or trans) of the membrane as shown in Figure 3. Throughout this study, this is referred to as the "escape time" (τ). The procedure is repeated a large number of times (5000 runs) for each polymer length N, and the escape time τ that occurred at most is recorded.
In the simulations, the average square end-to-end distance and radius of gyration in square lattice constant unit are calculated as follows. Because the polymer has different conformations at each MC step, we calculated the square end-to-end distance and square radius of gyration at each MCS time interval. Finally, we took the average over the total number of the intervals. End-to-end distance is the distance between the two end monomers, as the name implies. Its square at each MCS interval (R 2 k ) is computed as follows: where k = 1, 2, 3, ..., n (number of MC step intervals), ðx 1 , y 1 Þ, is the position of the first monomer, and ðx N , y N Þ is the position of the last monomer. Then, the average square end-to-end distance (R 2 ave ) is calculated by dividing the total sum of square end-to-end distance of the polymer at each interval for the total number of intervals (n).

Advances in Polymer Technology
We also calculated the average square radius of gyration as follows. First, we calculated the center of mass of the polymer at each MCS time interval ðxcm k , ycm k Þ as x ki , Hence, the radius of gyration (Rg 2 k ) at each MC step interval is Finally, the average square radius of gyration (Rg 2 ave ) is calculated as We repeated this procedure 5000 times (5000 runs) and averaged the results obtained at each run over the   Advances in Polymer Technology number of runs to get hR 2 i and hR 2 g i. The results are discussed in the following section.

Results and Discussion
3.1. Static Properties of Linear Polymer in a Crowded Environment 3.1.1. End-to-End Distance and Radius of Gyration. A polymer's end-to-end distance is a measure of how far it stretches out or how long its chain is. The number of monomers (N) in polymer chains affects this important static property illustrator parameter. And a radius of gyration is the proper quantity to describe all forms of polymers. Like end-to-end distance, it depends on the size of polymers, and ideal polymer models have indicated the power-law scaling relation between the parameters, mean square radius of gyration, and mean square end-to-end distance with chain length N as hR 2 g i~N 1/2 and hR 2 i~N 1/2 , respectively. But in assumption, ideal polymers does not obey excluded volume effect constraint between monomers. However, Flory and Carmesin [30] indicated that the exponent of the scaling relation is ν = 3/ð2 + dÞ, which is dimension-dependent, using clear and efficient approximation on the Rouse model of a SAW polymer chain, where d is dimension. In the plots of the mean square end-to-end distance (hR 2 i) and radius of gyration (hR 2 g i) as a fucntion of chain length N of the polymer (Figures 4 and 5), the interdependence has been shown. In our system, we only consider almost short chain length polymers, due to large amount of computational time needed to simulate very long polymers. The plots in Figure 4 indicate that in the absence of obstacles (free environment), the slope is 1:48 ± 0:01; so, the value of the slope is very close to 2ν (ν = 0:75 in two dimensions) which is in agreement with Flory and Carmesin's [30] scaling exponent of end-to-end distance with degree of polymerization, N. If there are no crowding agents at both sides of the membrane, for unbiased translocation of the considered linear polymer chains, the scaling relation of mean square end-to-end distance as a function of N is hR 2 i~N 2ν . However, as shown in Figure 4, in the presence of obstacles, the scaling exponent changes with the concentration of the obstacles.
And we observed that the end-to-end distance of the polymer increases, as ϕ increases. As the concentration of the obstacles increases, the bead obstacles close to the monomers, this causes the polymers to stretch [31]. Following this, as the polymer stretches, the entropy increases [3]. Therefore, the impact of entropic barrier on the translocation of the polymer increases, as density of the obstacles increases.
And also, the plots of mean square radius of gyration versus N, depicted in Figure 5, manifest the relationship between the two parameters. In the absence of obstacles, the slope = 1:49 ± 0:01 for the polymer is also in agreement with Flory's scaling exponent hR 2 g i~N 2ν . Nonetheless, the nonmonotonic fluctuation in polymer size with obstacle density has been reported in prior studies of polymer characteristics in crowded systems [29,32,33].
In the presence of obstacles, as it can be seen from the figure, for higher ϕ (ϕ > 0), the mean radius of gyration increases, and the universal power law relation of the mean square radius of gyration with chain length, hR 2 g i~N 2ν is violated, due to the conformations of the polymers is twisted by the crowding agents [29]. Therefore, the scaling relation depends on the density of the obstacles.

Advances in Polymer Technology
the probability distribution ðPðτÞÞ of the escape time τ for polymers translocating via the pore in various densities of barriers. Figure 6 depicts the probability distribution of linear polymer of a chain length N = 37 escape time in a crowded environment for ϕ ≤ 0:3.
As we can see from the figure, our simulations reveal that the most probable escape time increases as the environment becomes more crowded, i.e., the peak of the graph shifts to the right, as ϕ increases. And also, as delin-eated in the figure, our simulation results betrayed the histogram of the time is a long tailed probability distribution, as the probability distribution function spoils for the large escape times τ.
In addition, as there is no external force that pushes the polymers to either of the sides, the polymers escape to the side they prefer. From our simulations, as Figure 7 depicts, we also found that the translocation of the polymer prefers the disordered side (trans) line with Karplus and Andrew 0 2e + 07 4e + 07 6e + 07 8e + 07 1e + 08 1.2e + 08 1.4e + 08   Advances in Polymer Technology [20] and in the absence of obstacles, the probability of escaping to either of the sides is equally likely. As we mentioned in the simulation set up, the middle monomer is put inside the pore, and half of the polymer is on the left side, while the other half is on the right side. Unless an external force is applied to push it to the trans-side, the polymer will escape to either side with equal probability. However, in the presence of obstacles, the result shows that as the density of the crowding agent increases, the probability of translocation to trans-increases [20].

Effect of the Obstacles on the Relation of Escape Time and Chain
Length of the Polymer. One of the key concepts in polymer physics is the power law relation between translocation time τ and polymer size N. We obtain and take the most probable values of the escape time from the probability distribution of the escape time of each polymer size and plot them against the sizes N, as delineated in Figure 8. The relationship between the escape time of small linear polymers, like those contemplated in this study, and their chain size N is a power law. The slope of the loglog graph constitutes the scaling exponent α. Therefore, the power law relation is τ~N α . The results of our study manifest, α = 2:49 ± 0:01, for unbiased translocation of the polymers in the absence of obstacles. This result agrees with the Rouse model [34] prediction, which is α = 1 + 2ν, where ν = 0:75 in 2D. Nevertheless, the presence of obstacles shows a difference. In the presence of obstacles, Gopinathan and Kim [15] found that the translocation time decreases, as ϕ c increases. Inconsistent with this, Cao et al. [16] for both unforced (Δμ = 0) and forced (Δμ > 0), they found that the dependence of τ on ϕ c is nonmonotonic. And again, Chen and Lou [17], in the presence of driving force, they found that the escape timeτ decreases, as ϕ increases, and they observed that the scaling relation of τ and chain length N, α, is nonuniversal. Different from these, Karplus and Andrew [20], in their study of the translocation of linear and ring polymers in crowded environments without applying electric potential, they found that the presence of obstacles does not change the translocation time. Nevertheless, our simulation results enabled us to find that the scaling relation varies with the density of the crowding agents, as it can be seen from the figure. And the increment of escape time with ϕ is due to the entropic barrier rising when the concentration of the obstacles increases.

Diffusion of Linear Polymers with Effect of Obstacles.
The dynamic property of the polymer has also been studied by investigating the feature of the diffusion motion of the polymer. The time-dependent mean square displacement (h r 2 ðtÞi) set out the diffusion of the comprehensive system.
It is given by From the Rouse model, the following relations are expected [34]: The diffusion coefficient of polymers (D) is defined as Eq. (9), in the absence of obstacles. This implies that the slope of the graph of the respective mean square In the absence of obstacles, Cao et al. [16] and Ellery and Simpson [35] found that h r 2 ðtÞi~t β , β = 1, means the diffusion is normal diffusion, and in the presence of obstacles, the diffusion changes to subdiffusion (β < 1). In agreement with this, as it can be seen from Figure 9, we observed that the diffusion is Fickian diffusion (normal diffusion, means β = 1) for the absence of obstacles and subdiffusion (β < 1) for the presence of obstacles. The slope of the log-log plot of the mean square displacement vs. time determines the value of β.
We calculated the diffusion coefficient in square lattice constant unit per MCS time using the definition (Eq. (9)) for normal diffusion after computing the mean square displacement as a function of time for each size of polymer. The diffusion coefficient of each chain length N is calculated from the log-log plot of mean square displacement versus time, which is simulated in the absence of obstacles. We displayed the log-log plot of diffusion coefficient D as a function of N, as shown in Figure 10, to investigate the relationship between the diffusion coefficient and the polymer chain length.
In free environment, our simulations result for this relationship that can be put as a power-law scaling expression D~N −0:94±0:01 which is very close to the Rouse model scaling law of the form D~N −1 [34]. As it is depicted in the figure, the diffusion coefficient decreases as the size of the polymer increases.

Conclusion
For both static and dynamical simulations of linear polymers in crowded environments, we used and tested a Monte Carlo technique in our simulation investigation. To do this, we use a bond varying length between nearby monomer molecules in two dimensions (2D) to understand the polymer structures' pragmatic dynamics.
In a crowded environment with the same area fraction of crowding agents, ϕ, we explored polymer diffusion and polymer translocation through a nanoscopic pore. The scaling relations of end-to-end distance and radius of gyration of a linear polymer with chain length N greatly depend on the density of the barriers, ϕ, according to our numerical studies. Furthermore, we looked at how the scaling relations of hRi 2 and hRgi 2 with polymer length N vary with the obstacle density. The escape times τ of the different chain lengths of the polymer for different ϕ are also found from our simulations. And we discovered that τ rises as ϕ rises, resulting in a scaling relationship between τ and polymer length N, and α is not universal. We also found that the polymer prefers the disordered side (trans) to translocate.

Data Availability
The data required to support the findings of this study are included in the manuscript file.