Rational Surface Design for Molecular Dynamics Simulations of

7478

J. Phys. Chem. B 2008, 112, 7478–7488

Rational Surface Design for Molecular Dynamics Simulations of Porous Polymer Adsorbent Media E. Riccardi, J.-C. Wang, and A. I. Liapis* Department of Chemical and Biological Engineering, Missouri UniVersity of Science and Technology, 400 West 11th Street, Rolla, Missouri 65409-1230 ReceiVed: January 4, 2008; ReVised Manuscript ReceiVed: March 7, 2008

The construction and use of nonflat agarose surfaces in a simulation box, together with the employment of criteria for the immobilization of a set of dextran polymer chains on the nonflat agarose surfaces whose mathematical physics is compatible with that of the criteria used for the immobilization of the same set of dextran polymer chains on flat agarose surfaces, are shown to generate, through the use of molecular dynamics simulations whose simulation box has linear dimensions along the lateral directions that are the same when flat and nonflat agarose surfaces are used, dextran porous polymer structures whose pore sizes at the outermost surface and in the vicinity of the outermost surface of the porous medium can be controlled by an indirect manner through the variation of the parameters that characterize the nonflat surface. The use of a nonflat surface for the generation of desired large pores requires only a small or modest increase in the number of solvent molecules in the simulation box, while the use of a flat surface for the construction of the same desired large pores requires significant increases in the size of the linear dimensions of the flat surface. This increases so substantially the number of solvent molecules that the computational loads become intractable. The results in this work show that through the use of nonflat surfaces porous dextran polymer layers having pores of desired sizes can be effectively constructed, and this approach could be used for the design and construction of polymer-based porous adsorbent media that could effectively facilitate the transport and adsorption of an adsorbate biomolecule of interest that must be separated from a mixture of components. A useful definition about the properties that a porous polymer structure must have in order to become, for an adsorbate biomolecule of interest of known molecular size, a useful adsorbent medium, is presented and is used to (1) evaluate the porous polymer structures generated through the employment of different nonflat surface models and (2) determine and select the nonflat surface model from a set of nonflat surface models that is effective in producing promising porous structures. Then a procedure is presented by which a set of porous polymer media is generated through the use of the selected nonflat surface model, and the desired porous structure from this set is determined and could be considered to be used for the transport and immobilization of the selected affinity groups/ligands and the subsequent transport and adsorption of the desired to be separated adsorbate. 1. Introduction In the downstream processing of bioactive macromolecules in the pharmaceutical and biotechnology industries, significant increases in the separation of bioactive macromolecules are realized when ion-exchange chromatography, which utilizes porous adsorbent particles in which the affinity group/ligand is linked to the base matrix of the porous particle via a polymeric extender, is employed.1–14 The mutual interaction of the polymeric extender molecules creates a porous layer on the surface of the pores of the base matrix (stage I11) to provide a large surface area on which the selected affinity groups/ligands are immobilized (stage II11). The adsorption process of the charged biomolecules of interest could be considered to represent stage III11 in the study of the separation of biomolecules by ion-exchange adsorption employing adsorbent particles of the type discussed here. Macroscopic continuum models do not provide appropriate means for studying and understanding the properties of the discrete microscopic porous polymer layers with respect to their * To whom correspondence should be addressed. Telephone: (573) 3414414. Fax: (314) 966-2181. E-mail: [email protected].

three-dimensional (3D) structure characterized by (a) porosity, (b) pore-size distribution, and (c) number of conducting pathways along the direction of net transport (stage I11), as well as with respect to the dynamic effects that items (a)-(c) have on (1) the transport mechanisms and the resulting surface density of the affinity groups/ligands during their attachment onto the activated surface of the pores of the porous polymer layers (stage II11) and (2) the transport and adsorption mechanisms of the desired bioactive macromolecule onto the immobilized affinity groups/ligands (stage III11). On the other hand, molecular dynamics (MD) modeling and simulations provide a discrete microscopic modeling approach11 that could be used to study and understand the dynamic behavior of the porous polymer layers with respect to items (a)-(c) and (1) and (2) above, for stages I, II, and III. Zhang et al.11 modified the M3B model of Molinero and Goddard15 and used the modified M3B model in molecular dynamics studies that constructed porous dextran layers on the flat surface of a base matrix, where the dextran polymeric chains (i.e., polymer extenders) and the flat surface are covered by water. Two different porous polymer layers having 25 and 40 monomers per main polymer chain of dextran, respectively, were constructed,11 and a novel computational

10.1021/jp800078v CCC: $40.75  2008 American Chemical Society Published on Web 06/03/2008

MD Simulations of Porous Polymer Adsorbent Media

Figure 1. (a) Atomistic structure and (b) M3B coarse grain representation of a fragment of the dextran molecule, where beads B1, B4, and B6 of the M3B model correspond to the carbon positions C1, C4, and C6 in the atomistic structure.

method was developed11,16 and used to characterize the 3D porous structures of the porous polymer layers with respect to porosity, pore-size distribution, and number of conducting pathways along the direction of net transport. It is important to mention here that the structural properties exhibited by the porous polymer layers11 having 40 monomers per main polymer chain were for practical reasons more desirable than those of the porous polymer layers having 25 monomers per main polymer chain. The work of Zhang et al.11 demonstrated the efficacy of the MD modeling and simulation approach as a means to construct realistic porous polymer matrixes and study the behavior of such systems. It should be noted here that dextran polymeric extenders have been found in practice3–5,12,13 to provide porous polymer layers whose physicochemical properties have insignificant adverse effects on the bioactivity properties associated with the desired charged adsorbate molecules that are to be separated from a mixture of species by ion-exchange adsorption. Furthermore, the work of Zhang et al.11 on the characterization and analysis of the porous dextran structures has suggested a useful definition for the physical meaning and implications of the pore connectivity of a real porous medium, which is significantly different from the convenient with respect to modeling but artificial in physical meaning pore connectivity parameter employed in pore network models17–23 that have been used to determine the values of the transport parameters employed in macroscopic continuum models. The methodology developed by Zhang et al.11 for the characterization of the 3D structures of real porous media could be used to analyze the experimental data of the structural features of porous media obtained from high-resolution noninvasive24–29 3D methods such as high-resolution optical microscopy.26–29 In the work of Zhang et al.,11 the linear dimensions of the simulation box were 100 Å × 100 Å along the lateral x and y directions that define a flat surface (surface area 10000 Å2) on which 10 dextran chains with 5% side branching and one to two glucose units as side chains were randomly distributed. The dextran chains were completely immersed in water and formed a porous layer whose farthest distance from the base matrix was 100 Å along the z direction normal to the surface. To exhibit proper water density and represent the solution bulk phase,

J. Phys. Chem. B, Vol. 112, No. 25, 2008 7479 30000 water molecules were employed in the simulation system to obtain a water phase that not only covers the dextran layer but also extends 30 Å beyond the dextran layer. The MD simulations performed by Zhang et al.11 produced porous polymer structures whose effective pore sizes in the region located farthest from the surface of the base matrix would play a significantly more important role than the pore sizes in the region located next to the surface. In can be further reasoned that only when the effective sizes of the pores of the porous dextran region located farthest from the surface of the base matrix are substantially larger than the effective radii of the affinity groups/ligands and adsorbate molecules of interest could the dextran layer be considered as a potential candidate for immobilizing the adsorption sites (stage II) and adsorbing the charged adsorbate (stage III). A sensible approach to determining such effective molecular radii is through the experimentally measured hydrodynamic radius. Since in real practice the molecular size of the desired to be separated adsorbate, as well as the molecular size of the selected affinity group/ligand that would effectively adsorb and separate the desired adsorbate from a mixture of species, is known, the approximate effective size of the pores of the dextran layer, especially those in the region of the dextran layer located furthest from the surface of the base matrix and which could provide effective transport and adsorption for the adsorbate molecules, could be estimated in an a priori manner.22,30–34 However, in order to construct by MD modeling and simulations porous dextran layers with desired pore sizes for the effective separation of large adsorbate biomolecules by ion-exchange adsorption, the linear dimensions of the simulation box would have to increase so significantly that tremendous amounts of computational time would be required to carry out the MD simulations and render intractable the construction by MD simulations of porous dextran layers that have structural features similar to those used for effective ion-exchange adsorption in practice. In this work, a novel approach is presented with respect to the construction and treatment of the nonflat geometry of the surface of the base matrix on which the dextran polymer chains are attached. The variation of the controlled degree of curvature of the nonflat surfaces introduced in the present work is used as a means to control, in an indirect manner, the size of the pores of the porous dextran layers that are constructed by MD simulations and resulted from the interactions of the dextran chains randomly distributed on the nonflat surface. By varying the degree of curvature of the nonflat surface, different structures of porous dextran layers can be obtained and the characteristics of their porous structures are studied and presented in this work; these results represent data for stage I. In addition, the dependencies of the desired pore sizes on the type of the nonflat surface and the degree of nonlinearity of the nonflat surface are also analyzed in order to achieve useful guidance for future work that would study stages II and III. 2. System Formulation, Nonflat Surfaces, and Simulation Models and Methods In this work, 10 dextran polymer chains with 40 monomers per main polymer chain and 5% side branching, and having one end attached to and immobilized at randomly selected locations on the surface of the base matrix (i.e., agarose), are employed in the construction of the porous dextran polymer layer by MD simulation. The complexity of the system studied in this work is such that fully atomistic models would render the simulations intractable because of the computational times that would be required. Therefore, in order to make computa-

7480 J. Phys. Chem. B, Vol. 112, No. 25, 2008

Riccardi et al.

TABLE 1: Bending and Torsion Potential Parameters (Equations 1 and 2) angle type

Kθ (kcal/mol/rad2)

θ0 (deg)

161′ 616′ 461′ 416′

75 150 100 110

140.3 127.8 149.9 122.4

torsion type

B1 (kcal/mol)

φ01 (deg)

616′1′ 416′1′ 461′4′ 461′6′ 6416′ 4616′ 1461′ 4161′ 61′6′1′′

1.0 1.0 3.6 4.3 15.7 30.0 30.0 42.0 15.0

-85 -50 152 -145 90 -80 -55 45 140

For computational efficiency and consistency, it is reasonable to coarse grain one water molecule into one particle with proper size and interaction characteristics.11,15 The nonbonded like and unlike interactions between the coarse-grain beads of dextran and the water molecules are considered to be described by the Morse potential11,15

B2 φ02 B3 φ03 (kcal/mol) (deg) (kcal/mol) (deg) 1.8 2.1 2.5 2.8 4.3

-40 80 -85 110 -72

4.0

50

tionally tractable the simulations for our system and still retain physically relevant and meaningful information regarding the structure of the porous dextran polymer layers, a modified11 coarse-grain model based on the M3B model of Molinero and Goddard15 as shown in Figure 1 is employed for the representation of the dextran polymer chains. In general, the interaction between the coarse-grain beads of the dextran polymer chains can be expressed as a sum of valence and nonbonded terms,11,15 Ut ) Uv + Unb, where Ut represents the total interaction energy, Uv denotes the valence potential energy for bonds within the same dextran polymer chain, and Unb represents the nonbonded potential energy between different dextran polymer chains. For the bonds between two consecutive coarse-grain beads in the dextran, the SHAKE algorithm35,36 is used to constrain their bond lengths to equilibrium values adopted from previous studies.11,15 The valence potential energy, Uv, is therefore a sum of three-body interactions associated with the bending angle between three connected beads, Ub, and four-body interactions associated with the torsion angle between four connected beads, Utor. The coarse-grain angle bending motion is considered to be harmonic,11,15 and described by the following potential:

1 Ub(θ) ) Kθ(θi - θ0)2 2

(1)

where θi and θ0 denote instantaneous and equilibrium values of a bending angle, respectively, and Kθ is a coarse-grain bending force constant. For the four-body valence interactions, each coarse-grain torsion is considered11,15 to be represented by a sum of shifted dihedral functions: 3

Utor(φ) )

∑ 21 Bj[1 + cos(jφ - φj0)]

(2)

j)1

where the coarse-grain torsional angle φ ) φinkl denotes the angle between the planes formed by ink and nkl (i, n, k, and l are B1, B4, and B6) beads, and each term has a different integral periodicity j, barrier Bj, and phase φ0j . The values of Kθ, B1, B2, and B3 have been reported in previous studies.11,15 While the parameter values for θ0, φ01, φ02, and φ03 in our previous work11 provide one feasible representation, they are further optimized by using an approach37 based on the MM3 force field and a 10-monomer dextran chain. The new values used in this work are presented in Table 1. The dextran polymer chains in our simulation system are completely immersed in water as in the actual separation system.

{[ (

(

Unb(rij) ) D0,ij exp -0.5Rij

)) ] }

2 rij -1 -1 -1 R0,ij

(3)

where rij is the distance between a coarse-grain dextran bead and/or a water molecule i and j, R0,ij is the distance of minimum energy D0,ij, and Rij is the steepness parameter of the interaction potential. The values of the parameters D0,ij, R0,ij, and Rij have been fitted,11 and they are different for dextran beads (B1, B4, B6) in terminal or internal monomers. For unlike pair interactions, the mixing rule used is D0,ij ) (D0,iD0,j)1/2, R0,ij ) (R0,iR0,j)1/ 2, and Rij ) (1/2)(Ri + Rj).11,15 With such potentials, Zhang et al.11 have shown that the spatial distribution of the beads in the polymer chains composed of R-D(1f4)-linked or R -D(1f6)linked glucans is similar. The rational design of the structure of a nonflat surface where the degree of the surface nonlinearity could be appropriately varied produces the following three benefits: (i) It permits an increase in the effective surface area without increasing (within certain limits) the linear dimensions along the lateral x and y directions of the simulation box. Furthermore, by generating only a minor to a small increase in the height of the simulation box along the normal z direction, the introduction of a nonflat surface allows the simulation of a larger system to become computationally tractable because the increase in the number of water (solvent) molecules required for the implementation of a large system that employs a nonflat surface is significantly much smaller than the increase in the number of water molecules that is required for the same large system when a flat surface is employed, since in the latter case the linear dimensions of the simulation box along the lateral x and y directions would have to be increased significantly in order to accommodate the size of the system of interest to be simulated. (ii) It permits indirect control of the size of the pores of the porous dextran polymer layer formed by the interaction of the immobilized dextran polymer chains on the nonflat surface of agarose, so the degree of the surface nonlinearity influences in a kind of a priori manner the construction by MD simulations of pores whose sizes could accommodate the transport and adsorption of the desired adsorbate molecules in stage III. (iii) It provides the possibility to an investigator of simulating the behavior of a desired surface concentration of immobilized dextran polymer chains without having to modify the linear dimensions of the simulation box. It is important to mention here that Busini et al.38 have performed MD simulations employing a fully atomistic model for the agarose, but they have used very few agarose chains in their simulations and, thus, their data cannot provide a representative surface description for the agarose surface of our system on which the dextran polymer chains are immobilized. Since agarose is composed of similar monomers as dextran, the modified coarse-grain M3B11,15 model introduced above is used to develop agarose surface models for our simulations. For the agarose-dextran system considered in this work, the nonflat surface ought to be described by a continuous potential function to prevent water molecules from leaking out of the system and to provide the presence of sufficient space to fully accommodate the dextran polymer chains in their equilibrium structure on the surface and without interfering with (1) the criterion of randomness with regard to the selection (through


J. Phys. Chem. B, Vol. 112, No. 25, 2008 7481

Figure 2. Depictions of (a) nonflat surface model A, (b) nonflat surface model B, and (c) the immobilized ligand and adsorbate molecule at the pore surface located close to the outer surface of the porous dextran layer.

the use of a random number generator) of points on the agarose surface where one end of the polymer chain is immobilized, considering that the size of the surface area and the nonlinearity of the surface allow the selected set of immobilized polymer chains to form a porous polymer layer over the surface, and (2) the criterion for the direction of positioning the dextran polymer chain on the surface grid, which requires that the polymer chain is attached along the direction of the normal vector to the plane that is tangent to the nonflat surface at the randomly selected point on the surface where the polymer chain is fixed. This latter criterion for the direction of positioning a polymer chain at a randomly selected point on the nonflat surface is employed due to the fact that, when the surface of the simulation box is flat, the immobilized polymer chain’s direction is along the z axis of the simulation box that is normal to the flat surface defined by the lateral directions x and y of the simulation box.11 Among several possible nonflat surface structures that can satisfy the two above-mentioned criteria, two different surface models have been selected and employed in this work to represent the agarose surface as shown in Figure 2a. In surface model A, the nonflat surface is considered to have sinusoidal features. To realize this model with computational efficiency, we modify a surface model of Steele39 as follows to treat the interactions between the trigonometric surface and the particles that represent water molecules and dextran M3B beads:

[

( ) () ()

UisA(x, y, z) ) 2πεis

σis R0s

2

2 σis 5 zs

10

-

σis 4 zs

( )(

3

[

√2

R0s zs 0.61 R0s + σis σis √2 σis

( 2πb x) + cos( 2πb y)]

zs ) z - a 2 + cos

)] 3

of the dextran porous polymer layer, the structural details and the values of R0s and ε0s of the agarose substrate are not of major concern in this work. The parameters a and b determine the amplitude and wavelength of the surface trigonometric waves, respectively. In this work, we consider trigonometric waves in square patterns in a simulation box whose x and y linear dimensions are 100 Å, respectively. As a result, the nonflat surface model in eq 4a not only meets the criteria set forth above but also permits the extent of the surface nonlinearity to be varied by changing the amplitude and/or number of trigonometric waves. However, since those first dextran monomers to be immobilized have a fixed finite size and cannot overlap with the nonflat surface, certain restrictions should be followed concerning the values of a and b in order to not lose the advantages of the nonflat surface model A. For this purpose, we consider the projected size of one dextran monomer on the xy plane, depicted as δ in Figure 2a, to be smaller than the distance between the two inflection points around the valley of the trigonometric wave, which is equal to half the wavelength. As a result, the restriction followed in this work is b > 2δ. Another possible approach to creating nonflat surfaces is to embed relatively large spheres on a flat surface in place of trigonometric waves. Shown as surface model B in Figure 2b, this nonflat surface conceptually is a simpler one than surface model A and its interaction with various particles is considered to be a sum of two contributions:

UisB(x, y, z) ) Uis(z) + (4a)

(4b)

Here, (x,y,z) represents the coordinate of a particle and R0s is the nearest distance between surface atoms and is taken to be the same as the diameter of a dextran monomer based on the M3B model. The cross-interaction parameters σis and εis are evaluated as σis ) (1/2)(σi + R0s) and εi0 ) (εiε0s)1/2, where σi and εi are the diameter and interaction strength of particle i, respectively, and ε0s ) max(εB1,εB4,εB6) is the interaction strength between agarose monomers. It is worth pointing out that the lowest points on this surface correspond to zs ) z and coincide with a flat surface at z ) 0. Since no sufficient information is available to guide the construction of a coarsegrain surface model and our main interest is in the upper part

[

∑ Uik(x, y, z)

(5a)

k

( ) () ()

Uis(z) ) 2πεis

σis R0s

2

2 σis 5 z

10

-

σis 4 z

√2 R0s z 0.61 R0s 3 + σis σis √2 σis

( )( ) [( ) ( ) ]

Uik(x, y, z) ) 4εik

σik rik

12

-

σik rik

6

3

]

(5b)

(5c)

where Uis denotes the interaction with a flat surface, which is the same as that employed in the work of Zhang et al.11 and is dependent only on z, the distance from the flat surface. Uik denotes the interaction with the surface-embedded spheres indexed by k, and rik is the distance between a particle and the center of a sphere. The interaction model used here is the common Lennard-Jones 12-6 potential. For comparison with

7482 J. Phys. Chem. B, Vol. 112, No. 25, 2008 surface model A, we consider one case with 3 × 3 spheres uniformly embedded in the flat surface and in contact with each other. In this case, all the spheres have the same diameter σsp ) 33.333 and, thus, the same σik ) (σi + σsp)/2; each sphere is apparently equivalent to a few glucose monomers. To estimate εik, the interaction of a particle with M3B beads packed into a sphere is examined and coarse grained to obtain εik ≈ 5(εiε0s)1/2. Between the two nonflat surfaces, surface model B actually requires a longer time to implement than surface model A in our simulations over the same number of time steps. It is worth noting that, with 10 dextran chains, a minimum of nine concave surface segments (“hills” as per Figure 2) have been used in the simulation box for both surface models in order to minimize any artificial ordering effects due to the periodic boundary conditions. It is also worth emphasizing here that surface models A and B do not necessarily represent an accurate description of the geometry of the agarose surface, which is actually unknown. In fact, the accurate knowledge of the geometry of the agarose surface is not a strict requisite for the construction of the porous dextran polymer layers and their characterization, especially in the regions of the porous polymer structure located furthest from the surface which represent the most important porous regions with respect to transport and adsorption of the desired adsorbate biomolecules (stage III), as shown in Figure 2c. Because the dextran chains are attached perpendicularly to the local surface, they tend to collapse into the convex regions (regions I), leaving the concave regions (regions II) sparsely populated to form larger pores. The two nonflat surfaces thus allow the use of wavelength and amplitude in surface model A and sphere diameter in surface model B to control the formation of large pores that are most critical to biomolecular adsorption and transport. It is important to mention here that while the concept of the nonflat agarose surface is used as a theoretical tool with which a methodology is developed and used to produce the three important benefits discussed in items (i)-(iii) above, it is conceivable that practitioners3–5,12,13 who are engaged in the construction of agarose media for use as base matrices on which dextran polymer chains are immobilized could be interested in constructing in practice nonflat agarose surfaces with nonlinearities similar to those employed in the MD simulations and have resulted, as will be shown in this work, in producing dextran porous polymer layers that have desirable pore structures for use in the development of effective adsorbent media. The simulations in this work are carried out using the leapfrog integration algorithm together with a Gaussian thermostat method40 for temperature control. Each dextran chain has two side branches, one with single monomer and one with two monomers, that are randomly bonded to the main chain through B4 beads between the third and the 37th monomers. The dextran chains are first equilibrated without water, and they exhibit a one-dimensional helical structure. At the start of every simulation, these dextran chains are attached to the flat or nonflat surfaces at randomly selected points with random rotation around their long molecular axes and the B1-B6 bonds of their first monomers along the direction normal to the local surface. In case a dextran chain overlaps with the surface or existing dextran chains, it is reattached to a different random point. The water phase contains 32000 water molecules and is equilibrated separately. Afterwards, it is combined with the dextran-agarose surface system and the water molecules overlapping with dextran and surface are relocated to the top of the water phase with a downward velocity. The complete system has periodic boundary conditions applied in the x and y directions and is reequilibrated

Riccardi et al.

Figure 3. (a) Averaged number of pore openings and (b) averaged pore volume as a function of distance z from the surface for the porous dextran layer on the flat surface.

under 25 °C with dextran chains frozen during the first stage to allow a large integration time step of 12 fs to expedite equilibration. After the first stage, every particle is allowed to move except the first two monomers of the dextran chains in order to maintain their perpendicular configurations. The time step is switched back to 0.01 fs and gradually increases to 8 fs after having completed 10 ps of simulation time, and then the time step remains at 8 fs until the end of the simulation runs that span from 600 ps to 10 ns. Furthermore, it is important to report here that all the systems studied in this work and presented in Figures 3–12 attained equilibration (dynamic equilibration) in simulation times that are less than 600 ps. The results presented in Figures 3–9, 11, and 12 describe the structural properties of the dextran porous polymer layers obtained at 600 ps. It should be noted that employing before the start of the MD simulation of the complete system the simulation strategy with respect to the water molecules and the dextran polymer chains, and which was discussed above, has led to shorter simulation times for equilibration. Also, it is important to note here that the dextran polymer chains that interact with each other to form the porous polymer layer are not totally stationary or static at their equilibrated state because the dextran polymer chains are flexible and find themselves in a balanced form of movement that deals with bodies that remain at apparent rest (dynamic equilibrium). To facilitate the removal of unwanted influence from the initial condition, an annealing procedure is used where the system temperature is raised to 80 °C and remains at the temperature over an extended period of time before cooling to 25 °C again to start a production run. 3. Results and Discussion As demonstrated in our previous study,11 the analyses of the number of pore openings and pore volume together provide adequate characterization and evaluation of the porous medium constructed by discrete particles. For the studies undertaken in this work, a number of molecular dynamics simulations have


Figure 4. (a) Averaged number of pore openings and (b) averaged pore volume as a function of distance from z ) 0 (agarose surface) to the outermost surface of the porous dextran layer on the nonflat surface model A with a ) 6.249 Å and b ) 33.333 Å.

Figure 5. (a) Averaged number of pore openings and (b) averaged pore volume as a function of distance from z ) 0 (agarose surface) to the outermost surface of the porous dextran layer on the nonflat surface model B with nine spheres and dextran chains immobilized only on the surface of the nine spheres.

been performed that consider the flat surface and the nonflat surfaces with different wavelength and amplitude for surface model A and one sphere diameter for surface model B. It should be noted here that the major differences between the work of Zhang et al.11 and the work reported here are as follows: this works introduces and uses (a) both flat and nonflat surfaces in


Figure 6. (a) Averaged number of pore openings and (b) averaged pore volume as a function of distance from z ) 0 (agarose surface) to the outermost surface of the porous dextran layer on the nonflat surface model B with nine spheres and dextran chains immobilized on both the sphere surface and the flat parts between spheres.

Figure 7. (a) Number of pore openings and (b) pore volume as a function of distance from z ) 0 (agarose surface) for the outermost surface of the porous dextran layer of the most representative realization of the system employing the flat surface.

the simulation box while Zhang et al.11 used only flat surfaces in the simulation box, (b) appropriate potential functions for the nonflat surfaces employed in the simulations, and (c) new improved values for the parameters shown in Table 1 of this work and which were determined by using an approach37 based on the MM3 force field. The results in Figures 3–6 present the

7484 J. Phys. Chem. B, Vol. 112, No. 25, 2008


Figure 9. (a) Inside side view of the porous structure of the dextran layer of the system shown in Figure 4 and (b) inside side view of the porous structure of the dextran layer of the system shown in Figure 8.

distributions of the number of pore openings and pore volume averaged from multiple realizations and categorized by different ranges of pore radii along the distance z from the surface of the

Riccardi et al.

Figure 10. (a) Averaged number of pore openings and (b) averaged pore volume as a function of distance from z ) 0 (agarose surface) to the outermost surface of the porous dextran layer on the nonflat surface model A with a ) 4.688 Å and b ) 25.000 Å. The simulation time is 2 ns.

base matrix, respectively, for (i) the flat surface, (ii) the nonflat surface model A with a ) 6.249 Å and b ) 33.333 Å that imply nine surface trigonometric waves and the ratio 4a/b being equal to 0.75, (iii) the nonflat surface model B with nine spheres of diameter 33.333 Å and the 10 dextran polymer chains immobilized only on the surface of the spheres, and (iv) the nonflat surface model B with nine spheres of diameter 33.333 Å similar to case (iii) but with the dextran polymer chains allowed to be immobilized not only on the surface of the spheres but also on the flat parts existing between the spheres. In all cases, the linear dimensions along the lateral x and y directions are equal to 100 Å. It is worth mentioning here that the increment of 5 Å that defines the largest difference in the magnitudes of the pore radii considered in each size range indicated in Figures 3–6 (as well as in the figures discussed later on) was selected on the basis of the conventional consideration where, in histograms representing pore volumes for different pore sizes, the increment employed between pore sizes is often taken to be equal to 5 Å so that the information on the histogram is not cluttered; of course, a different magnitude for this increment could have been selected. By comparing the results in Figures 3–6, it is observed that the nonflat surface models A and B provide a larger number of large pores in the range 20-25 Å and a larger pore volume in this range of pore sizes than those obtained when the flat surface is used. Furthermore, the pore volume available at the region located close to the upper surface of the porous polymer layer comprised from the pores of size 20-25 Å is significantly larger in the systems employing the nonflat surfaces than that obtained from using the flat surface. A better insight can be ascertained from Figure 7, where the results of the most representative realization for the case of the flat surface are presented. It can be observed in Figure 7 that the pores in the size range 20-25 Å do not start from the outermost surface of the porous polymer layer but rather from z ≈ 66.5 Å; thus there

MD Simulations of Porous Polymer Adsorbent Media are no pore openings and pore volume available from these large pores for the entrance and effective transport and adsorption of biomolecules at the outermost surface and for a distance of about 25 Å from the outermost surface into the porous polymer structure. In general, for effective biomolecular transport and adsorption, it would be required to have pores whose diameters are at least 2-3 times larger than the hydrodynamic diameter of the adsorbate molecules and which are available at the outermost surface of the porous polymer layer. While many proteins have hydrodynamic diameters that are as large in magnitude as the diameters of the large pores encountered in the dextran porous polymer layers of the systems in Figures 3–6, and, therefore, these porous media would not be appropriate for facilitating the adsorption and transport of such proteins, it is apparent that the pore sizes in the systems of Figures 3–6 could accommodate biomolecules such as peptides (e.g., desmopressin) whose hydrodynamic diameters are about 1/3 or less than the diameter of the largest pores of the systems in Figures 3–6. Furthermore, the purpose of the current work is to show the proof of concept which indicates that the employment of nonflat surfaces in the simulation box can provide an indirect method for controlling the size and pore-size distribution of the pores in the dextran porous polymer layers, and thus, if the hydrodynamic diameter of the adsorbate biomolecule of interest (peptide or protein) is known, then one could always use the methodology presented in this work with appropriate modifications in the length of the linear dimensions of the simulation box and in the number of immobilized on the agarose surface dextran polymer chains, to construct dextran porous polymer structures whose large pores located at the outermost surface of the porous layer are at least 2-3 times larger than the hydrodynamic diameter of the adsorbate biomolecule. It is important to mention here that the most representative realizations for surface model A having the values of a and b reported in Figure 4 and surface model B for the case described in item (iv) above and whose averaged results are presented in Figure 6 have the large pores of size 20-25 Å located at the outermost surface of the porous polymer layer, as is desired to be the case. Furthermore, another important criterion in assessing the structure of a porous medium is the pore connectivity11 that represents a measure of the number of conducting pathways that facilitate the transport along the direction of net transport (the z direction in our systems) in the porous medium; the magnitude of this number of conducting pathways is related to the average total number of pore openings of different pore radii across lateral planes perpendicular to the direction of net transport, and to the average length along the direction of net transport of the pores originating from the outermost surface of the porous medium and ending at the point where a discontinuity (a discontinuity denotes that the number of pore openings is equal to zero for a given range of pore radii) in a significant number of pore openings of different pore radii occurs. In general, the more the distribution curves in Figures 3a, 4a, 5a, 6a, and 7a intersect each other and are continuous along z, the higher pore connectivity the porous polymer medium has. Therefore, Figures 3a, 4a, 5a, 6a, and 7a clearly indicate that the pore connectivity of the larger pores represented by the size ranges of 15-20 and 20-25 Å is significantly higher in the porous structures formed through the employment of the nonflat surface models A and B than that obtained when the flat surface is used. In summary, the results in Figures 3–7 and the above discussion have clearly shown that, with the same linear dimensions along the x and y directions, the nonflat surface models A and B produce (1) a larger number of pore openings



of large size, (2) a larger pore volume associated with the pores of large size at the outermost surface and in the region of the porous polymer layer located close to the outermost surface so that large in size adsorbate molecules could enter and be transported and adsorbed in the porous polymer layer, and (3) a higher pore connectivity for the larger in size pores in the porous medium than the flat surface for the same number of immobilized dextran polymer chains. A comparison of the results in Figures 4–6 indicates that the pore connectivity of the larger pores represented by the size ranges of 15-20 and 20-25 Å obtained from the nonflat surface model A is significantly higher than those obtained from the nonflat surface model B and, furthermore, the pore volume associated with the largest pore size range (20-25 Å) in the system employing the nonflat surface A (as per Figure 4) is continuous along z for the whole length of the porous polymer layer and is larger in magnitude than that obtained from the system in Figure 5 that employs the nonflat surface B, while the pore volume associated with the same pore size range (20-25 Å) in the system shown in Figure 6, where the nonflat surface model B is employed, is discontinuous along z. A discontinuity11 in a pore volume distribution curve shown in Figures 3b, 4b, 5b, 6b, and 7b is defined as a position in the z direction at which the value of the pore volume becomes equal to 0. Since biomolecules to be transported and adsorbed in a porous polymer layer enter from the outermost surface, it would be most desirable for any given curve in Figures 4–7 to begin from the largest value of z located at the outermost surface of the porous polymer layer and end at the smallest possible value of z. From this point of view, a discontinuity in a curve in Figures 4–7 is indicative of the relevant pores to be dead-ended and the adsorbate biomolecules likely to be trapped there. As can be understood and shown in Figures 6 and 7, the undesirable discontinuities in the number of pore openings and pore volume coincide with each other. In Figures 3b, 4b, 5b, 6b, and 7b, the pore volume distributions that are continuous have shoulders

7486 J. Phys. Chem. B, Vol. 112, No. 25, 2008 on the corresponding curves which reflect the fluctuations along the z direction of the number of pore openings in the relevant size ranges. By comparing the two porous structures obtained when using the nonflat surface model B and presented in Figures 5 and 6, the curve of pore volume having the largest pore size range of 20-25 Å in Figure 5b is continuous and has a number of shoulders while the curve in Figure 6b has no shoulders and is undesirably discontinuous, although the total number of pore openings in the size range of 20-25 Å in Figure 6a is greater than that in Figure 5a. Therefore, it could be stated that the larger the number of shoulders on a given curve, the more probable that the pore volume distribution is continuous. In Figure 6b, the fact that the pore volume curve having a pore size range of 20-25 Å has two smooth parts without shoulders and separated by a discontinuity indicates an undesirable case where only two regions in the porous layer have the large pores but they are disconnected. Two additional important variables of interest in the characterization of porous media are those associated with the pore surface area and its distribution along z, the direction of net transport. The results in Figures 4b, 5b, and 6b can be considered to provide a relative measure with respect to the pore surface areas and their corresponding distributions along the direction of net transport for the porous structures obtained by using the nonflat surface models A and B. It is indicated by the results in these figures that, for the pores of size range 20-25 Å, the porous structure obtained from the employment of the nonflat surface model A has a relatively larger and better distributed pore surface area along the z direction than those obtained from the use of the nonflat surface model B. The analysis and discussion presented above regarding the results in Figures 4–6 have thus indicated that the use of the nonflat surface model A produces porous polymer structures having significantly more desirable features than those obtained from the employment of the nonflat surface model B. Therefore, we judiciously focus on using the nonflat surface model A to continue this work. In Figure 8 the average number of pore openings and pore volume as a function of z determined from multiple realizations are presented for a system that employs the nonflat surface model A with a ) 4.688 Å and b ) 25 Å (nine waves and 4a/b ) 0.75); the systems of Figures 4 and 8 have the same surface area of 15278.61 Å2 and also the same linear dimensions along the lateral x and y directions that are equal to 100 Å, respectively. By comparing the results in Figures 4 and 8, it can be observed that while the larger pores having sizes in the ranges 15-20 and 20-25 Å are characterized by continuity along z from the outermost surface of the porous polymer layer all the way to a region close to the agarose surface, the pore volume associated with the largest in size pores (20-25 Å) of the system in Figure 8b is larger over the total length along z than that of the system in Figure 4b. Also, the results in Figure 8a indicate that the pore connectivity of the pores corresponding to the large pore size ranges of 15-20 and 20-25 Å is comparable to the pore connectivity of the system in Figure 4a along the z direction, but the pore connectivity of the system in Figure 8a is found to be higher in the upper region of the porous polymer layer, which represents the operational region of interest for transport and adsorption of adsorbate molecules. Furthermore, the results in Figures 4b and 8b discussed above indicate that the porous polymer layer of the system in Figure 8 has a relatively larger along the z direction pore surface area than the system in Figure 4, and, thus, the use of the nonflat surface model A with a ) 4.688 Å and b ) 25.000 Å provides a porous dextran polymer layer that has significantly better and more

Riccardi et al.


desirable features for (a) the transport and immobilization of affinity groups/ligands in stage II, and (b) the transport and adsorption of desired adsorbate biomolecules in stage III, than the porous dextran polymer layer obtained from using the nonflat surface model A with a ) 6.249 Å and b ) 33.333 Å, even though the two systems have the same surface area. This difference in porous structures occurs because the two sets of a and b values have different effects on the distribution and the relative deviation of the direction of the immobilized dextran polymer chains on the agarose surface with respect to the z direction through the criteria, discussed in section 2 of this work, that are required to be followed for determining the positioning of the 10 dextran polymer chains on the agarose surface. As implied by eq 4b, the positioning direction along which the dextran polymer chain is immobilized on the surface is a stronger function of the wavelength parameter b than the amplitude parameter a. Another illustration that can indicate visual differences between the porous structures of the systems in Figures 4 and 8 could be obtained by looking at the inside side-view visualizations shown in Figures 9a and 9b of the porous structures of the systems presented in Figures 4 and 8, respectively. These visualizations are snapshots taken at 600 ps from the most representative realization of the systems in Figures 4 and 8. It should be noted here that the visualizations presented in Figure 9 are showing the inside side view of the dextran porous structures at a single position of depth and could be different at other depth positions. Furthermore, due to the periodic boundary conditions employed in the MD simulations, it is worth mentioning that the visualizations in Figure 9 are repetitive in the porous structures for the same position of depth. It is important to mention here that the side-view visualizations can only provide a relative visual qualitative structural difference between the porous media of the systems in Figures 4 and 8, while it was shown above that only the approach represented by the data in Figures 4a, 4b, 8a, and 8b can provide useful quantitative information for evaluating and comparing the porous

MD Simulations of Porous Polymer Adsorbent Media structures formed by the interacting dextran polymer chains. In Figure 10, the averaged from multiple realizations structural properties of the system considered in Figure 8 are presented when the simulation time is equal to 2 ns, while the data in Figure 8 represent results obtained from a simulation time of 600 ps. It can be observed that the structural characteristics of the systems in Figures 8 and 10 are equivalent and indicate that equilibration (dynamic equilibrium) had been attained at 600 ps. The analysis and discussion of the results presented in the above paragraphs has led to the selection of the functional form of the nonflat surface that provides pore structures having desirable properties (i.e., surface model A provides desirable pore structures). After the functional form of the nonflat surface has been appropriately determined, then for a selected (i) set of values for the linear dimensions along the lateral x and y directions, (ii) number of polymer chains, and (iii) set of values for the parameters a and b that satisfies b > 2δ in order to prevent overlapping of the immobilized polymer chains with the surface, as discussed in section 2, and the criterion of random distribution of the selected number of polymer chains, different porous polymer structures could be constructed by MD simulations through the variation of the wavelength b while keeping the surface area the same as that obtained from the initially selected values of a and b. This computational procedure should be repeated for additional values of the amplitude a to explore other achievable porous structures under different surface areas. It should be mentioned here that the selection of the values in items (i)-(iii) above should involve judicious choices whose goal would be to enable an investigator to construct porous polymer structures where the largest pores would have diameters that are at least 2-3 times larger than the hydrodynamic diameter of the adsorbate molecule of interest. The different porous polymer structures that would be obtained from the MD simulation procedure described above could be studied and compared by examining their porous structures through the presentation and analyses of their data, as discussed in this work (cf. Figures 4 and 8). From this set of generated porous structures, the porous polymer medium that exhibits high pore connectivity and continuous pore volume distribution for the large diameter pores that start from the outermost surface of the porous polymer structure should then be selected as a desirable porous polymer structure to be employed in the transport and immobilization of the affinity groups/ligands (stage II) and in the subsequent transport and adsorption of the desired adsorbate molecule (stage III). In Figures 11 and 12, the results determined from multiple realizations for porous polymer structures that use the nonflat surface model A with a ) 5.625 Å, b ) 25.000 Å, and a ) 6.563 Å, b ) 25.000 Å, are represented, respectively; the surface areas for the systems in Figures 11 and 12 are 17053.35 Å2 and 18927.37 Å2, respectively, while the linear dimensions along the lateral x and y directions are still 100 Å, respectively. The number of dextran polymer chains immobilized on the surface of these two systems is equal to 10, as was the case for the systems in Figures 3–10. It can be observed from the results in Figures 11 and 12 that the increased surface areas, when compared with the surface areas of the systems in Figures 4 and 8, lead to the formation of pores whose diameters are in the range of 25-30 Å, which were not obtained in the systems in Figures 4 and 8. Furthermore, it is important to note that the system in Figure 12 has also pores with diameters in the range 20-25 Å while the system in Figure 11 has no such pores. Also, the results in Figures 11a and 12a suggest that the system

J. Phys. Chem. B, Vol. 112, No. 25, 2008 7487 in Figure 12a has a higher pore connectivity for the large in size pores than the system in Figure 11a, and also the data in Figures 11b and 12b indicate that the surface area of the pore structure in the region close to the outermost surface of the porous polymer layer generated by the large pores having diameters in the size ranges of 20-25 and 25-30 Å is larger for the system in Figure 12b than that in Figure 11b. 4. Conclusions and Remarks It was shown that, by employing nonflat surfaces in the simulation box and using a set of criteria for the immobilization of the dextran polymer chains on the nonflat surfaces that are compatible in mathematical physics with those employed for the immobilization of the same polymer chains on a flat surface having the same linear lateral dimensions as the nonflat surface, porous dextran polymer structures could be constructed on model agarose surfaces that have larger pores at the outermost surface of the porous polymer structure than the pore sizes obtained when the flat surface is used. This interesting and useful result obtained from the employment of the nonflat surfaces is of practical importance because most adsorbate biomolecules of interest to be separated by adsorption have rather large hydrodynamic diameters and, thus, require for their transport and adsorption large pores starting at the outermost surface of the porous polymer structure and extending along the direction z that represents the direction of net conduction for the mass flux of adsorbate molecules. Compared to the simulation that uses a flat surface, the use of the nonflat surfaces requires only a small increase in the number of solvent molecules and, therefore, the desired result of constructing porous media whose pore sizes can be controlled indirectly by the variation of the parameters that characterize the nonflat surface can be accomplished without substantial increase in computational loads. It is worth mentioning again here that while the concept of using a nonflat agarose surface in the simulation box has been shown to provide a powerful theoretical tool that can be used to control in an indirect way the size of the pores and the pore-size distribution of the dextran porous polymer structure immobilized on the agarose surface so that practically desirable porous structures could be constructed for a given adsorbate of interest, it is conceivable that practitioners3–5,12,13 who are involved in the construction of agarose media for use as base matrixes on which dextran polymer chains are immobilized could be interested in constructing different nonflat agarose surfaces which could have the potential to generate dextran porous polymer structures that have desirable pore sizes, pore-size distributions, and pore connectivities for a given adsorbate species of interest. By defining that a porous polymer structure is desirable if it has (a) suitable in size large pores that start at the outermost surface of the porous medium and can accommodate effectively the transport and adsorption of the adsorbate molecules, (b) high pore connectivity for the desired in size large pores along the z direction that represents the direction of net conduction for the mass flux of the adsorbate molecules, and (c) continuous pore volume distributions for the desired large pores along the z direction so that there is a relatively large and properly distributed surface area in the porous polymer structure associated with the desired large pores, it was shown that the use of the nonflat trigonometric surface model A generated substantially more desirable porous dextran polymer structures than those obtained when the nonflat surface model B with surfaceembedded spheres was used, and, thus, the nonflat surface model A was considered to be more suitable for further use in stages

7488 J. Phys. Chem. B, Vol. 112, No. 25, 2008 II and III. Furthermore, after the comparison (based on items (a)-(c) above) of different dextran porous polymer structures generated from a set of different nonflat surface models leads to the selection of the proper nonflat surface model, it is important to construct by the selected nonflat surface model the porous structure that would best facilitate the transport and adsorption of the adsorbate biomolecule of interest whose hydrodynamic diameter is known. This latter requirement could be realized through the procedure presented in this work and which generates a set of porous polymer structures using different values for the parameters that characterize the selected nonflat surface model, and then by employing the criteria (items (a)-(c) above) that establish what constitutes, for a given adsorbate molecule of interest, a desirable porous polymer medium, the porous polymer structure that could best facilitate the transport of the desired adsorbate molecule is identified and selected to be employed in the two subsequent stages involving (i) the transport and immobilization of the affinity groups/ligands (stage II) and (ii) the transport and adsorption of the adsorbate species of interest (stage III). References and Notes (1) Liapis, A. I.; Grimes, B. A.; Lacki, K.; Neretnierks, I. J. Chromatogr., A 2001, 921, 135. (2) Grimes, B. A.; Liapis, A. I. J. Colloid Interface Sci. 2002, 248, 504. (3) Berg, H.; Hansson, H.; Kagedal, L. Adsorption/Separation Method and a Medium for Adsorption/Separation. U.S. Patent 6,428,707 B1, Aug 6, 2002. (4) Johansson, B.-L.; Belew, M.; Eriksson, S.; Glad, G.; Lind, O.; Maloisel, J.-L.; Norman, N. J. Chromatogr., A 2003, 1016, 21. (5) Johansson, B.-L.; Belew, M.; Eriksson, S.; Glad, G.; Lind, O.; Maloisel, J-L.; Norman, N. J. Chromatogr., A 2003, 1016, 35. (6) Zhang, X.; Grimes, B. A.; Wang, J.-C.; Lacki, K. M.; Liapis, A. I. J. Colloid Interface Sci. 2004, 273, 22. (7) Zhang, X.; Wang, J.-C.; Lacki, K. M.; Liapis, A. I. J. Colloid Interface Sci. 2004, 277, 483. (8) Zhang, X.; Wang, J.-C.; Lacki, K. M.; Liapis, A. I. J. Colloid Interface Sci. 2005, 290, 373. (9) Liapis, A. I.; Grimes, B. A. J. Sep. Sci. 2005, 28, 1909. (10) Liapis, A. I. Ind. Eng. Chem. Res. 2005, 44, 5380.

Riccardi et al. (11) Zhang, X.; Wang, J.-C.; Lacki, K. M.; Liapis, A. I. J. Phys. Chem. B 2005, 109, 21028. (12) Harinarayan, C.; Mueller, J.; Ljunglo¨f, A.; Fahrner, R.; Van Alstine, J. M.; van Reis, R. Biotechnol. Bioeng. 2006, 95, 775. (13) Ljunglo¨f, A.; Lacki, K. M.; Mueller, J.; Harinarayan, C.; van Reis, R.; Fahrner, R.; Van Alstine, J. M. Biotechnol. Bioeng. 2007, 96, 515. (14) Liapis, A. I.; Grimes, B. A. J. Sep. Sci. 2007, 30, 648. (15) Molinero, V.; Goddard, W. A. J. Phys. Chem. B 2004, 108, 1414. (16) Zhang, X. Ph.D. Dissertation, Department of Chemical and Biological Engineering, University of MissourisRolla, Rolla, MO, 2005. (17) Rege, S. D.; Fogler, H. S. AIChE J. 1988, 34, 1761. (18) Sahimi, M.; Jue, V. Phys. ReV. Lett. 1989, 62, 629. (19) Petropoulos, J. H.; Petrou, J. K.; Liapis, A. I. Ind. Eng. Chem. Res. 1991, 30, 1281. (20) Hollewand, M. P.; Gladden, L. F. Chem. Eng. Sci. 1992, 47, 1761. (21) Stauffer, D.; Aharony, A. Introduction to Percolation Theory, 2nd ed.; Taylor and Francis: London, 1992. (22) Meyers, J. J.; Liapis, A. I. J. Chromatogr., A 1998, 827, 197. (23) Bryntesson, L. M. J. Chromatogr., A 2002, 945, 103. (24) Latour, L. L.; Kleinberg, R. L.; Mitra, P. P.; Sotak, C. H. J. Magn. Reson. 1995, 112, 83. (25) Lori, N. F.; Conturo, T. E.; Bihan, D. L. J. Magn. Reson. 2003, 165, 185. (26) Egner, A.; Jordas, S.; Hell, S. W. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 3370. (27) Levi, V.; Ruan, Q.; Gratton, E. Biophys. J. 2005, 88, 2919. (28) Levi, V.; Ruan, Q.; Plutz, M.; Belmont, A. S.; Gratton, E. Biophys. J. 2005, 89, 4275. (29) Ragan, T.; Huang, H.; So, P.; Gratton, E. J. Fluoresc. 2006, 16, 325. (30) Renkin, E. M. J. Gen. Physiol. 1954, 38, 225. (31) Brenner, H.; Gaydos, L. J. J. Colloid Interface Sci. 1977, 58, 312. (32) Mason, E. A.; Wendt, R. P.; Bresler, E. H. J. Membr. Sci. 1980, 6, 283. (33) Petropoulos, J. H.; Liapis, A. I.; Kolliopoulos, N. P.; Petrou, J. K.; Kanellopoulos, N. K. Bioseparation 1990, 1, 69. (34) Dullien, F. A. L. Porous Media: Fluid Transport and Pore Structure, 2nd ed.; Academic Press: New York, 1992. (35) Rapaport, D. C. The Art of Molecular Dynamic Simulation; Cambridge University Press: New York, 1995. (36) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, U.K., 1987. (37) Bohne, A.; Lang, E.; von der Lieth, C. W. J. Mol. Model. 1998, 4, 33. (38) Busini, V.; Moiani, D.; Moscatelli, D.; Zamolo, L.; Cavallotti, C. J. Phys. Chem. B 2006, 110, 23564. (39) Steele, W. A. Surf. Sci. 1973, 36, 317. (40) Brown, D.; Clarke, J. H. R. Mol. Phys. 1984, 51, 1243.

JP800078V

Rational Surface Design for Molecular Dynamics Simulations of

Recommend Documents