1. Introduction
Aromatic heterocyclic polyimides and their composites attract a great attention during the last years as perspective membrane materials [
1,
2,
3,
4,
5,
6]. The presence of the polar groups and water permeability make these polymers and composites on their base the promising candidates for preparation of membranes for dehydration of different types of organic solvents [
2,
3,
7,
8,
9]. Matrimid
® is one of the commercially available amorphous polyimide that demonstrates excellent thermal stability, good mechanical properties, high glass transition temperature and superior chemical stability. These properties make Matrimid
® a promising polymeric material for separation industry.
Permeability of penetrant during pervaporation depends on two processes: sorption and diffusion, both of which are affected by the polymer-penetrant interactions. Diffusion of water in the polymers depends on temperature, packing of polymer chains and polymer–water interaction. It is known that water molecules may present in polymer in different states. It was shown by differential scanning calorimetry investigations for poly (vinyl alcohol) [
10] and sodium poly (styrenesulfonate) [
11] that three states of water exist: (i) nonfreezing water (bound water), (ii) freezable bound water and (iii) free water.
Additionally, the polymer can swells under influence of penetrant, which leads to the decreasing of density of polymer chains and increasing of water mobility in membrane. Computer simulation of polymers and their permeability for penetrants with molecular dynamics (MD) method is an effective tool for the deeper understanding of the mechanism of molecular mobility in membranes and allows proposing adequate models for further investigations [
12,
13,
14].
As it was pointed by Nguyen et al. [
15], the clustering of penetrant molecules occurs if they exhibit stronger affinity to each other than to the polymer and mobility of polymer chains provokes local volume expansion. The clusterization process influences the sorption, diffusion and permeation properties of membranes. Clustering leads to association of molecules into groups consisting of two or three species, which diffuse as the single unit. Usually, it leads to the decreasing of mobility of penetrant in membrane [
16]. It was shown by Sun et al. [
17] that the diffusion of organic vapours in the poly (dimethyl siloxane) membrane is complicated by the Langmuir-type penetrant immobilization and formation of clusters. The authors pointed out that these two processes (immobilisation on sorption centres and cluster formation) can be more important than additional free-volume, which appears due to swelling of polymer under influence of penetrant. Zimm-Lundberg model is often applied for determination of activity of sorbate at which clusterization occurs [
18]. Another model, which allows simultaneous calculation of amount of sorbate in the cluster and bonded to the sorption centres, is the model of Laatikainen-Lindstrom [
19]. Nevertheless, both these methods are indirect tools for the finding conditions of clusterization on the base of sorption data; both methods have limitations, which arise from initial postulates used for the deriving of model equations.
Summarizing literature data, it can be derived that the sorption on specific sorption centres and formation of clusters reduce the mobility of penetrant inside the polymer, while swelling of macromolecules often leads to the increasing of the diffusion rate. As a result, theoretical analysis of the behaviour of penetrant in polyheteroarylene membrane will be useful for understanding its separation performance at different composition of feed mixture.
The simulation of penetrant behaviour inside the polyheteroarylenes can give evidences for the proposed mechanisms of slowing of the penetrant diffusion in polymeric membranes. The aim of our work was the analysis the concentration dependence of mobility of water molecules in the polyheteroarylene Matrimid
®5218 using the molecular dynamics simulation. The choosing of Matrimid
®5218 as the object is primarily connected with the growing demand for membrane technology. As mentioned above, membranes on the base of Matrimid
® have stable technological characteristics that are typical for the membranes on the basis of polyheteroarylenes. In spite Matrimid
® earlier was applied mainly for gas separation [
20,
21] many recent works are connected with its application for separation of liquid mixtures, in pervaporation or nanofiltration processes [
22,
23]. Nevertheless, the consideration of peculiarities of liquids permeation is mostly descriptive, that is, quite understandable due to practical directions of these investigations. For example, hydrophilicity of Matrimid is explained by the preferential interaction between the water molecules and imide groups in polymer (results from contact angle determination [
22]). The analysis of molecular origin of the different penetrant interactions with Matrimid is given in [
23]. In the interesting work of Stanford et al. [
24] the sorption and diffusion of penetrants in Matrimid is considered. In particular, the diffusional behavior of water is in contrast to the alcohols: the diffusion coefficient slightly decreases as a function of activity. Stanford et al. [
24] explained it by the clustering behavior of water at higher activities. Our study is an attempt to obtain some additional information on the mobility of penetrants in Matrimide, as membrane material, especially for water.
2. Materials and Methods
The material selected for this study is poly (5 amino-1-(4′aminophenyl)-1,3-trimethyl indane), trade name Matrimid
®5218. Chemical structure of Matrimid
®5218 repeating unit is presented in
Figure 1. On the first stage of computer simulation we calculated the electronic structure of dimer of Matrimid
®. We used density functional theory (DFT) method realized in DMol
3 [
25,
26] module of Materials Studio package (version 7.0) with Perdew-Burke-Ernzerhof (PBE) functional [
27]. The charges on atoms were estimated with Mulliken scheme. The dimer was considered instead of monomer in order to find the charges of atoms near the bond between monomers.
The dimers were combined into octamers of Matrimid
® to simulate the structure of polymer. The cube cell with periodic boundary condition was fulfilled with 8 octamers using amorphous cell module. The initial density was set 0.1 g·cm
−3. Forcite module from materials studio was applied for classical MD simulation. The choice to use a short-range force field was made in favour of Universal Force Field (UFF) [
28] with calculated charges on atoms. The UFF gives the best results for polymer density in comparison with other force fields (COMPASS II, DREIDING). The calculated ensemble was isothermal-isobaric ensemble (NPT) with Nose thermostat and Berendsen barostat. To create the model of polymer, the cell was prepared at
T = 1000 K and equilibrated during 1 ns using MD simulations in canonical ensemble (NVT). The autocorrelation function of distances between the ends of octamers was used for analysing of quality of equilibration. Next, the temperature was decreased stepwise with 100 K per step up to 300 K. At each step the 1 ns of MD simulation was performed for equilibration at NPT regime (
p = 1 atm). As a result, the relaxed structure of polymer was obtained with density 1.175 g·cm
−3, which is close to the experimental value 1.2 g·cm
−3. The resulted cell size was 37 Å × 37 Å × 37 Å. The cell was filled with water molecules at
T = 300 K using sorption module of materials studio and fixed loading regime. The numbers of water molecules per cell were 10, 15, 20, 30, 40 and 60, which correspond to water concentrations of 0.5 wt.%, 0.75 wt.%, 1 wt.%, 1.5 wt.%, 2 wt.% and 3 wt.%. These concentrations were chosen according to the experimental data on sorption capacity of Matrimid
® respectively to water vapour [
16]. The larger cell with size 55 Å × 55 Å × 55 Å containing 27 octamers and 135 water molecules (2 wt.%) was also constructed with the same methodology and used for MD simulations in order to evaluate the influence of cell size onto the obtained results. The mean square displacement (MSD) was calculated using standard technique in Forcite module from materials studio.
3. Results
Mean square displacement curves for water molecules at different concentrations in polymer are given in
Figure 2a. The small amount of water does not provoke the conjugation of water molecules in clusters. The disposition of MSD curves in
Figure 2a does not follow the concentrations due to interaction not only water molecules between each other, but also with water clusters which are formed at higher water content. In order to verify the time region, at which the analysis of MSD curves will be performed the curves had been also plotted in a log-log scale (
Figure 2c). For all experiments the slope of curves was lower than 1 for the times up to 2.5 ns, which is the limit of performed simulation. Thus, the subdiffusion regime of water movement was observed in all cases and the analysis of results is presented in terms of the following equation:
where
t is the time, α is an exponent and
D* is the coefficient, which will be considered as the value proportional to the water mobility. In order to evaluate the region of stability of α, the derivative of log-log plot of MSD curve was plotted for MD in the cell containing 27 polymer chains (see
Figure 2b). It is seen that the approximately constant value of α is observed up to
t = 600 ps, and gradual increase is observed at higher
t. The subdiffusion nature of water mobility can be caused by bonding of water with the polar groups of the polymer such as C=O groups and by the hampering of water penetration through the channels between polymer chains. The value of
t = 100 ps was chosen as the starting time for the analysis of water mobility on the base of the observation that before this time some of MSD curves in log-log scale demonstrate inflection point. This is more clearly seen in
Figure 2c for the MSD curves in log-log scale.
Based on the discussion above, the MSD curves in log-log scale in the range of 100 to 600 ps were fitted with linear equation (red lines in
Figure 2d) and the parameters of subdiffusion were estimated using the logarithm of Equation (1): log
10(MSD) = log
10(
D*) + αlog
10(
t). The exponent value α, shown as blue curve in
Figure 2d, demonstrates the increasing in the concentration range of 0.5 wt.% to 1.5 wt.% and reaching of the approximately constant value 0.43 at higher concentrations. This can be attributed to the sorption of water molecules onto the active sorption centres of polymer (C=O groups). Thus, 0.5% can be considered as strongly bonded water. Limited amount of accessible sorption centres determines the maximal amount of strongly bonded water and increasing of water content leads to the fast increasing of α.
As it is seen in
Figure 2d, the values of α for concentrations >1.5 wt.% are equal to 0.43, which demonstrates the hampering of water diffusion. Additionally, the value of parameter
D* decreases with increasing of water content. In addition to the data for 8 polymer chains, the single result for 27 chains is also presented in
Figure 2d. The good agreement between parameters obtained for cells with eight and 27 polymer chains confirms the correct choice of optimal number (eight) of simulated chains. From
Figure 2d it can be also proposed that the mobility of water in the polymer is restricted by the additional factor, which is directly proportional to the concentration. The explanation for this phenomenon was suggested earlier on the base of experimental data [
24]: the nontypical behaviour of water diffusion in Matrimid
® can be connected with clusterization of water inside the volume of polymer. Clusterization leads to the increasing of size of penetrant units, which results in decreasing of probability of their penetration between macromolecular chains. This suggestion makes it interesting to analyse the radial distribution function (g(
r)) of distance between oxygen atoms of water molecules at different concentrations.
In
Figure 3, the values of g(
r) are presented for simulations with 0.5 wt.%, 1.0 wt.%, 1.5 wt.% and 3.0 wt.% of water, respectively. The characteristic of the distance between oxygen atoms for water molecules, connected with formation of hydrogen bonds, is approximately 3 Å. Thus, the peak with the maximum of approximately 3 Å in the g(
r) function may be considered as a measure of the probability of the formation of hydrogen bonds between water molecules inside the polymer. It can be seen that at the lowest water concentration (
Figure 3a) this peak is very small, which reflects high distance and weak interaction between water molecules. This may be connected with strong bonding of water molecules onto the polar groups of the polymer. To check this assumption, the radial distribution function (RDF) of water around oxygen atoms in carbonyl groups of polyheteroarylene were plotted in
Figure 4a,b for the water content in polymer 0.5 wt.% and 1.5 wt.%, respectively. The difference between 1.0 wt.% and 1.5 wt.% is concerned with the beginning of creation of water clusters. The statistical values of g(
r) reflect the behaviour of atoms with additional types of bonding. The lower concentration does not provide the stability of water clusters in time. It can be seen, that these RDF are very similar to each other, which proves sorption of water on C=O groups of polymer. It can be also mentioned that adhesion to sorption centres gives significant impact to the hampering of diffusion only at very low concentration of water in polymer. Although the amount of sorption centres in the cell used for this study was 256 (4 C=O groups per monomer multiplied with eight octamers in cell), they can bond only ~10 water molecules. It allows proposing that the most part of active sorption centres is not accessible for the interaction with water because of dense chain packing of polyheteroarylene.
From
Figure 3a–c it is seen that the intensity of peak near 3 Å increases when water content growth from 0.5 till 1.5. This result demonstrates the increasing of probability of water–water interaction and the formation of clusters. These clusters exist up to maximal investigated content of water, which can be derived from significant similarity of RDF presented in
Figure 3c,d.
The mean cluster size (MCS) of sorbate inside the polymer is one of the significant parameter for understanding of penetrant state in the membrane. The MCS can be calculated using the model equations of Zimm-Lundberg or Laatikainen-Lindstrom on the base of sorption data. Simulations allow calculating MCS directly from the distribution of water inside the polymer. For this purpose, the positions of oxygen atoms of water in the moment of time
t = 500 ps were extracted from the trajectories and were analyzed. The exact water molecule was considered as belonging to exact cluster if the distance between water and at least one water molecule in the cluster was lower than
r = 5 Å. This value of r was chosen on the base of water–water RDF functions presented in
Figure 3. The value 5 Å is the point where peak of g(
r) is ending. MCS was calculated as
where
Ni is the number of clusters containing
i water molecules. The results are presented in
Table 1.
As it presented in
Table 1, the appearance of clusters at water content >0.75% is observed. The increasing of water concentration leads to the increasing of number sizes of clusters. This is expressed in the increasing of MCS from 1, for water contents 0.75 wt.% and 0.5 wt.%, to 1.4, for the simulation with 3 wt.% of water. The mean cluster sizes calculated for cells with eight and 27 polymer chains at 2 wt.% of water were 1.33 and 1.31 respectively, which demonstrates the independency of obtained MCS on the cell size. These results give evidence for the suggestion that decreasing of mobility of water at concentrations >1 wt.% is connected with adhesion of water molecules to each other and formation of clusters. This conclusion on the role of cluster behaviour is in a good agreement with experimental results [
24] for Matrimid and similar conclusions on water diffusion in other polymers [
29].