Article Cite This: Macromolecules XXXX, XXX, XXX−XXX
pubs.acs.org/Macromolecules
Multiscale Study of the Relationship between Photoisomerization and Mechanical Behavior of Azo-Polymer Based on the CoarseGrained Molecular Dynamics Simulation Junghwan Moon,† Byungjo Kim,† Joonmyung Choi,§ and Maenghyo Cho*,†,‡ Division of Multiscale Mechanical Design, School of Mechanical and Aerospace Engineering, and ‡Institute of Advanced Machines and Design, Seoul National University, Seoul, Republic of Korea § Department of Mechanical Engineering, Hanyang University, Ansan, Gyeonggi-do, Republic of Korea Macromolecules Downloaded from pubs.acs.org by UNIV OF TEXAS AT DALLAS on 03/02/19. For personal use only.
†
S Supporting Information *
ABSTRACT: Cross-linked liquid crystalline polymers (CLCPs) incorporated with photochromic azobenzene moieties can exhibit large and reversible deformations under the exposure of actinic lights. To practically utilize CLCPs in applications such as flexible photoresponsive actuators and soft robots, the user-defined diverse deformation modes have been demonstrated by controlling the liquid crystalline phases and polymer architectures. However, conventional all-atom molecular dynamics (AA MD) simulation is unsuitable for modeling complicated phases and conducting a parametric study of the photomechanical effect in terms of the polymer morphology due to its excessive computational cost. Therefore, a multiscale analysis framework based on coarse-grained molecular dynamics (CG MD) simulations is developed in this study. The mesoscale CLCP network model is constructed to satisfy the structural and thermodynamic properties of the full-atomistic reference. Also, the photoisomerization process is operated by integrating multiscale phenomena as follows. The AA MD reference describes the quantum events by implementing reactive rotational potential for trans-to-cis conversion, and corresponding CG potential sets are derived based on the associated structural changes. Because the mesoscale photoswitching potential is transferable to a wide range of cis-concentrations, the light-induced sequential smectic−nematic−isotropic phase transition is realized. The results show that corruptions of the smectic layer and nematic orientational order reversely influence the overall shape of the network, which causes different photodeformation modes under each transition. Furthermore, the light-induced disordering behaviors are combined with the mechanical deformations of the polymer according to various light and temperature conditions. We expect the presented work contributes to the precise prediction and design of the optomechanical behaviors of photoresponsive devices.
I. INTRODUCTION Photodeformable polymers have promising features compared with other stimuli-responsive materials such as fast responsiveness, remote and wireless control, and utilization of light as an infinite energy source.1,2 In particular, nematic cross-linked liquid crystalline polymers (CLCP) incorporated with photoactive azobenzene derivatives can show reversible large deformations in response to ultraviolet (UV) or visible light irradiation. In a pioneering work, Finkelmann’s group3 introduced remarkable contraction of liquid crystalline elastomer (LCE) film doped with photochromic molecules under UV exposure. The novel optomechanical effect is the result of the cooperative work of the photoresponsive liquid crystals (LCs) and constituent polymer chains.4 Absorption of the photon energy with a certain frequency (365 nm) manipulates the molecular shape of the azobenzene from nematic rod-like (trans-state) to the kinked shape (cis-state). The cis-molecules generated by the photoisomerization process disturb the nematic alignment of the surrounding liquid © XXXX American Chemical Society
crystals (LCs), resulting in light-induced phase transition. The photoisomerization reduces the shape anisotropy of the polymer network along the initial nematic direction, causing the polymer network to exhibit remarkable macroscopic photoshrinkage behavior. Although the photoisomerization of the azo-dyes can be widely utilized in optics, biomedical applications, and photocontrollable micrometer-sized surfaces,5−9 one of the most attractive traits is that the photosensitive CLCPs can convert the light energy into mechanical work. Ikeda et al.10 realized the photobending behavior of the azo-containing polymer, which is enabled by the contractive deformation gradient along the thickness direction of the film. It is due to the rapidly decaying light intensity with respect to the penetration depth. In addition to the bending of the nematic liquid crystalline Received: November 28, 2018 Revised: February 18, 2019
A
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
source. In this study, we expect to discover the mesoscopic mechanism of the macromolecular behavior under the lightdriven phase transition of the smectic LCP. In the research field of the atomistic modeling of the LCs, the hybrid GayBerne/Lennard-Jones interaction potentials have been utilized to characterize the diverse complex structures such as the smectic-A and smectic-B phases.22,29,30 However, the arbitrarily defined potential coefficients cannot explain the subtle structural changes influenced by the types of mesogen, spacer, and backbone chains. Therefore, we derived the CG potential energy sets which can rigorously reproduce the structural and thermodynamic states of the AA MD model. We systematically developed the mesoscale photoreactive potentials for the first time by integrating the physical properties obtained at each different scale. Our previous studies simulated the photochemical process of the fullatomistic azo-LCP by adopting a photoswitchable force field in terms of the appropriate reaction coordinate.19,31 Extending these works, the CG interaction sets are developed, which can capture both light-induced conformational change and the subsequent disordering effect on the surrounding molecules. The gradual phase transition behaviors and mechanical deformations can be examined by applying modified potential to the trans-isomers. While other mesoscale studies have characterized the microstructural difference between two isomers and replaced the photoisomerization effect with the reorientation of the LCs,32,33 these studies are not realistic or not applicable to cases with a wide range of cis-concentrations. Another contribution of this study is to provide the relationship between the LC order parameters and macromolecular shape anisotropy according to internal (initial LC phase and morphology) and external (temperature and photoisomerization ratio) conditions. The mathematical model is defined to couple the evolutions of the translational and orientational orders of the LCs with the overall shape parameter under the smectic−nematic−isotropic transition. Linking the photoinduced effect in mesoscale to the continuum scale stress−strain constitutive relation is an important procedure.18,34 We expect our findings will assist in predicting and designing the complex 3D deformations of the photosensitive actuators and robots with a consideration of the multiscale phenomena.
polymers (LCPs), preprogrammed complex three-dimensional (3D) photodeformations have been achieved.9,11,12 These diverse stimuli-responsive motions can be realized by the control of the local alignments or defects of the LC domains, the use of more complex LC phases (smectic, twisted nematic, etc.), and the architectural design of the polymeric structures.13−15 Combinations of the light-driven deformation modes are promising for the envisaged applications in optomechanical actuators and soft robotic devices.16,17 To accurately predict and design the photomechanical deformations, their multiphysics and multiscale nature should be first understood. To achieve this, diverse scale-physical information needs to be integrated from the quantum event, represented by the photoexcitation, to the macroscopic alternation of the shape. Many researchers have conducted a computational analysis of the light-triggered effects by using quantum mechanics (QM), molecular dynamics (MD), and continuum scale finite element (FE) simulations.18−20 However, there is a lack of a unified solution derived by bridging the physical quantities obtained from the calculations of each scale. Our previous multiscale framework21 showed that the all-atom molecular dynamics (AA MD) simulation can serve as an intermediate scale analysis to provide microscopic molecular shape change in terms of the extent of photochemical reaction. However, the conventional MD simulations with full-atomistic description have obvious limitations in modeling the complex LCP network. First, the AA MD can only exhibit a light-induced nematic-to-isotropic transition and corresponding deformations. While the arrangement of the nematic director requires approximately several nanoseconds, the growth of the smectic phase requires tens or hundreds of nanoseconds.22,23 In terms of the length scale, a few micrometers scale is required to observe convincing twisted nematic structures.22 However, this is beyond the conventional simulation scope of the AA MD due to the excessive computational cost. Another limitation of the conventional MD simulations is the difficulty in considering the complex macromolecular morphology. It has been proven that the optomechanical properties such as the responsiveness, strain, and exerted stress can be tuned by modifying the spacer length and cross-linking density of the azo-LCP network.24,25 Generally, in the case of the side-chain LCP, the spacing between each backbone is about 3−5 nm, and it is comparable to the side length of the unit cell in the AA MD simulation. Accordingly, the full-atomistic computational model is not suitable for exploring the photomechanical effect of the CLCP with diverse LC phases and polymer architectures. In this work, a multiscale optomechanical analysis is conducted by using the coarse-grained molecular dynamics (CG MD) simulations to overcome the limitations of the conventional AA MD. Especially, the simulation scope is extended to the mesoscale phenomena to observe the thermaland light-responsive behaviors of the smectic and the nematic solids. In contrast to the nematic polymers, the smectic LCPs possess both rotational and translational symmetries. Because of the layered structure formed by the self-assembly of the LCs, the highly aligned smectic CLCP is favorable to induce fast photobending deformation and strong mechanical response.26,27 Furthermore, the smectic solid can exhibit novel photomechanical deformations. Yu et al.28 reported sequential switching of the photobending direction: the CLCP first bent in the opposite direction to the light source during the smectic−nematic transition and then bent toward the light
2. METHODS 2.1. Bead Mapping and Formation of LCP Network. The azo-LCP network considered here is composed of the A11AB6 LC monomers and DA11AB LC cross-linkers. Yu and co-workers28 reported that the CLCP film containing these polymer molecules can exhibit a thermotropic smectic A (Sm A −nematic (N)−isotropic (I) transition. Figures 1a,b show the chemical structure and coarse-grained model of the LC monomer and cross-linker. Each coordinate of the CG beads is the group center of mass of the associated segment. First, the trans-azobenzene mesogenic unit is mapped into a three-bead system, which consists of two phenyl (P) beads connected to a centered azo (N) bead. Acrylate polymer chains containing a spacer of undecylene are mapped into four C3 type beads, one CO bead, and one C2 bead. The C3 bead represents the three methylene units. In the case of the bead next to the azobenzene, two methylene units and one oxygen atom are assigned to the bead. One carbon atom and two oxygen atoms in the ester group are denoted together as one CO bead. The C2 bead represents either two methylene units B
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
or a remaining hydrocarbon segment closest to the CO bead. A short alkoxy polymer chain of the LC monomer is mapped into one C3 bead and two C2 beads. Peter et al.33 applied a similar mapping methodology to 8AB8 molecules and proved that it is suitable to decouple intramolecular degrees of freedom. To model the AA MD reference structure, an uncured initial bulk unit cell is constructed with the rigid mesogens oriented to the nematic axis. The molar ratio of the LC monomers and cross-linking agents is 9:1. The un-cross-linked structure is relaxed via sequential steps consisting of conjugate gradient energy minimization, the NVT ensemble for 300 ps at 300 K, and the NPT ensemble for 1.5 ns at 300 K with a time step of 1.0 fs. Then, the molecules are effectively polymerized to form a well-defined network via a two-step cross-linking process, as shown in Figure 2. In the first step, the partially polymerized molecules are generated by primarily linking the diacrylate cross-linkers with surrounding acrylate monomers. This is a crucial procedure to connect polymer backbones and to fix the network structure. Covalent bonding is formed between the reactive atoms as their interatomic distance is the same or less than the designated cutoff. The reaction cutoff distance for the AA MD model is gradually increased from 1.0 to 5.4 Å, with an increment of 0.2 Å. Before shifting the cutoff, an intermediate energetically favorable topology is constructed by applying the energy minimization and NVT runs at 400 K for 50 ps. The first cross-linking simulation proceeds until about 63% of the reactive atoms of the cross-linkers have reacted with each
Figure 1. (a) Chemical structures and (b) CG models of the acrylate LC monomer and LC cross-linker.
Figure 2. AA and CG MD configurations which illustrate how the LC compounds are polymerized into the CLCP structure. The polymeric backbone formed during the cross-linking simulation is highlighted by yellow color in the unit cell. C
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Table 1. Bond Stretching and Angle Bending CG Potential Energy Coefficients (Functional Forms Are Presented in Eqs 3 and 4) bond length
a1
P−N (trans and cis) P−C3 C3−C3 C3−C2 C2−C2 (backbone) C2−C2 (terminal) C3−CO C2−CO bending angle
0.691 0.286 0.114 0.239 0.441 0.209 0.211 0.937
N−P−C3 P−N−P (trans) P−N−P (cis) P−C3−C3 C3−C3−C3 P−C3−C2 C3−C2−C2 C3−C3−CO C3−CO−C2 CO−C2−C2 C2−C2−C2
a1 1.188 1.332 0.066 0.053 0.048 0.048 0.010 0.090 0.032 0.025 0.014
b1 [Å]
l1 [Å]
0.077 0.082 0.117 0.203 0.110 0.323 0.088 0.104 b1 [deg]
3.216 3.797 3.589 2.602 2.609 2.678 3.281 2.416 θ1 [deg]
10.030 9.071 8.724 21.790 6.191 20.010 21.580 11.620 14.990 5.129 15.191
174.9 179.9 95.9 157.0 173.1 144.3 109.3 179.0 121.8 76.27 112.5
a2
b2 [Å]
l2 [Å]
0.095 0.065 0.142 0.334 0.045 0.171
0.164 0.358 0.190 0.125 0.169 0.218
3.698 3.387 2.971 2.379 3.252 3.002
a2
b2 [deg]
θ2 [deg]
0.203
5.722
162.4
0.005 0.006
32.990 3.896
115.7 164.4
0.009 0.032 0.020 0.013 0.105
14.050 35.100 14.010 4.172 19.41
136.1 153.5 146.0 95.7 171.4
a3
b3 [Å]
l3 [Å]
0.092
0.090
3.873
0.027
0.168
2.539
a3
b3 [deg]
θ3 [deg]
0.054
36.800
170.5
0.050
21.130
165.5
0.014
24.850
98.4
are derived to describe the conformational structure of the CLCP network in the atomistic resolution. Meanwhile, the mesoscale intermolecular potential is elaborately optimized to reproduce the structural and thermodynamic properties that are consistent with those of the AA MD reference via the iterative Boltzmann inversion (IBI) method. Because the CG potential is dependent on the state point of the sampled MD model, the choice of reference temperature and the corresponding LC phase is of significant importance.33,35,36 Experimentally, the CLCP presented here exhibits a smecticto-nematic transition at ∼360 K and a nematic-to-isotropic transition at ∼430 K.28 The choice of an all-atom bulk CLCP in the temperature range of the smectic phase (300−360 K) as a target model causes problems in parametrizing the CG potentials because the full-atom MD model cannot rigorously describe the smectic structure. In addition, the radial distribution function (RDF) of the smectic LCP is unsuitable for applying the IBI method because of its translational symmetry. Therefore, we constructed the bonded and nonbonded CG potentials from the nematic ordered CLCP and the isotropic constituent fragments, respectively. The reference temperature for the classical MD structure is 400 K. Although the sampled RDF of the mesogens is obtained from the isotropic liquid, it has been revealed that the short-range correlation between the LCs in the ordered phase does not significantly differ from that of the disordered phase in an appropriate thermodynamic state.33,36 The transferability of the developed mesoscale potential for a wide range of temperatures (300−420 K) and molecular composition is tested by comparing the LC phase behavior of the CG model with that of the experimental result. 2.2.1. Bonded Interactions. The bonded interactions are divided into contributions from relevant conformational quantities: bond length (l) between two adjacent beads and bending angle (θ) between three adjacent beads. The probability distribution functions of each degree of freedom (PCG(l) and PCG(θ)) are obtained by sampling the ensemble trajectories of the sufficiently equilibrated CLCP structure.
other. During the second step, the remaining monomers are polymerized to create long polymer strands with a maximum distance of 7.4 Å. The overall cross-link density of the final model is about 40%. An additional relaxation process is executed by using further minimization, NPT runs for 10 ns at 400 K and 1.0 atm, and NVT runs for 5 ns. Unlike the all-atom model, randomly placing the uncured monomers within the unit cell should be avoided in constructing the initial CG model. It is because excessive polymerization events are required to model a complex network structure, and the polydomain nematic phase can be developed during relaxation. Therefore, a partially cross-linked azo-LCP consisting of 18 acrylate monomers and 2 diacrylate cross-linkers is chosen as the unit molecule. Four smectic layers of the unit molecules (including 1920 mesogenic groups) are placed in an orthorhombic simulation box with the x-axis in the layer normal direction. The CG bulk unit cell consists of 23616 beads, which are equivalent to 165696 atoms. The energetically favorable structure of the Sm A LCP is constructed via stepwise equilibration consisting of geometry optimization and the NPT ensemble for 10 ns at 300 K and 0.1 MPa with a time step of 10.0 fs. The size of the initial unit cell is about 14 nm × 12 nm × 10 nm. During the first crosslinking step, the reactive C2 beads at the terminal of the crosslinkers are connected to both ends of the acrylate backbones. The cutoff distance for a chemical reaction increases from 4.5 to5.5 Å with an increment of 0.5 Å. The energy minimization procedure and NPT runs at 400 K for 3 ns are applied whenever several curing events occur. Next, the backbones of the side-chain LCP are assembled with a maximum cutoff of 9.0 Å. The final cross-link density is the same as that of the AA MD model. Additional relaxation runs are performed with further isobaric NPT runs at 400 K for 10 ns and at 300 K for 120 ns after the completion of the polymerization. 2.2. CG Interaction Potentials. Generally, the total potential energy of the bead system can be separated into bonded and nonbonded energy terms. In this study, the CG potential parameters of the intramolecular degrees of freedom D
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Figure 3. Bond length distributions which are obtained from the isolated 8AB8 molecule and from all-atom and coarse-grained models of the nematic cross-linked liquid crystalline polymer at 400 K.
analytically fitted intramolecular force fields are listed in Table 1. Owing to the difference in the probability of the distance between two C2 beads at the terminal of the LC monomer and acrylate backbone, different sets of structural potential parameters are adopted. Some studies have reported that the CG bonded potential derived by sampling conformation of an isolated LCP molecule cannot appropriately reproduce the ordered structure of the bulk system.35,36 Therefore, we compared the structural distribution of the nematic CLCP network and that of a single 8AB8 in a vacuum. The bead mapping scheme and simulation conditions for sampling the conformation of the isolated 8AB8 molecule can be found in ref 33. Figures 3 and 4 show the bond length and bending angle distributions, respectively. The structural conformation of the CG bead model is in good agreement with that of the parent CLCP with an all-atom description. The deviation of the conformational probability of the CG model from that of the AA MD reference occurs because of the different network formation histories. The CG polymer network is cured in the deep smectic phase. Therefore, the cross-linking process at highly ordered phase forms more aligned polymer network. The deviations do not significantly affect the thermo- and lightinduced phase transition behaviors of the smectic CLCP. It
The detailed simulation conditions for modeling the AA MD reference are given in the Supporting Information. The conformational probabilities are Boltzmann inverted to derive the configurational free energy: UCG(l) = −kBT ln(PCG(l)/l 2)
(1)
UCG(θ ) = −kBT ln(PCG(θ )/sin θ )
(2)
Here, kB is the Boltzmann constant, and T is the reference temperature. The distribution of structural parameters is approximated as a sum of the multicentered Gaussian functions.37 The CG bond and angle potentials in eqs 1 and 2 are then expressed as follows: É ÅÅÄ n 2 |Ñ ÑÑ l i ÅÅ o o y Ñ l − l o o j z ÑÑ o iz UCG(l) = −kBT lnÅÅÅÅ∑ ai expo m −jjj zo }ÑÑ j bi zz o ÅÅ ÑÑ o o o {o (3) n k ~ÑÑÖ ÅÅÇ i = 1
É ÅÄÅ n 2 Ñ ÑÑ l i ÅÅ o o yz | Ñ θ − θ o o j Å o o i zz }ÑÑÑ UCG(θ ) = −kBT lnÅÅÅ∑ ai expo m −jjj z o j b z oÑÑÑ ÅÅ o o ÅÅÇ i = 1 (4) n k i {o ~ÑÑÖ where li and θi are the center positions of the ith Gaussian functions, and ai and bi are constants. The parameters of the
E
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Figure 4. Bending angle distributions which are obtained from the isolated 8AB8 molecule and from all-atom and coarse-grained models of the nematic cross-linked liquid crystalline polymer at 400 K.
are required. It is very difficult to refine such a large number of parameters simultaneously. Therefore, the cross-linked LCP network structure is divided into several small fragments to obtain each potential separately. We constructed three different atomistic fragment systems at 400 K, and the intermolecular interactions of the trans-LCP were then determined based on the RDF and density of each sample. The interaction potentials between the P−P, P−N, and N−N pairs are derived from the trans-azobenzene liquids. The atomistic structural data of the octadecane molecules are used to parametrize the C3−C3 or C2−C2 interactions. Finally, the
should be noted that the structure of the single azo-LCP molecule does not match with the nematic bulk system, especially for the angle between the beads (see Figures 4a,c,d). This is because the molecular alignment between the nematic LCs and the constraint by the cross-links result in a more uncoiled configuration of the polymer chains. 2.2.2. Nonbonded Interactions. The noncovalent CG potentials are optimized to reproduce the target RDF as well as the density of the isotropic molecules. The system considered in this study is composed of five different types of beads, and a total of 15 different nonbonded interaction sets F
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Figure 5. Radial distribution functions which are obtained from the all-atom isotropic fragments of the azo-LCP and from corresponding CG models at 400 K.
(UPC CG(r)), an attractive linear tail function is added to the original potential:38
CO−CO interactions are derived to reproduce the structure of the octyl propanoate liquids. Detailed information about the construction of the fragment unit cell is given in the Supporting Information. The IBI method is an iterative optimization procedure that obtains effective CG nonbonded interactions. In this method, an initial trial potential (UCG,0(r)) can be calculated by the Boltzmann inversion of the target RDF (g(r)) of the reference: UCG,0(r ) = −kBT ln(g (r ))
r yz ji PC UCG (r ) = UCG(r ) − γkBT jjj1 − zzz j rc z{ k
γ can be modified until the pressure of the CG system equals the target pressure (0.1 MPa). We systematically combined the IBI and pressure-correction (PC) steps. Initially, only IBI is applied with β = 1.0 until f target is below 1.0 × 10−4. Then, IBI and PC are implemented at the same time with β = 0.5 and γ = 0.03 until the pressure is about ±0.02 GPa. In the final step, the pressure is elaborately corrected with γ = 10−5−10−3. Figure 5 compares the radial distributions, which are obtained from the atomistic reference and the CG model with the refined potentials. Table 2 lists the deviations of the RDF and density of the CG systems from those of the AA MD unit cells. The precisely refined potential sets derived in this study yield
(5)
where r is the separation distance between the CG units. The CG MD trajectories with the initial guess do not satisfy the target properties. Therefore, the accuracy of the nonbonded energy is gradually improved by using the difference in the potentials of the mean force obtained using the radial distributions generated from the ith iteration (gi(r)) and from the AA MD reference. This step can be expressed as ij g (r ) zy UCG, i + 1(r ) = UCG, i(r ) − βkBT lnjjjj i zzzz k g (r ) {
Table 2. Basic Information on the Unit Cells of the Reference Fragments and the Deviations of the RDF and Density of the CG Models from Those of the Atomistic Simulations
(6)
Here, β is an adjustable constant for stable optimization. We continued the iteration until a target function (f target), which quantifies the deviation of the RDF, reaches below the specific threshold value. It can be calculated using the equation ftarget =
∫0
rc
exp( −r ){gi(r ) − g (r )}2 dr
(8)
isotropic fragments transazobenzene transazobenzene transazobenzene octadecane octadecane octyl propanoate
(7)
where rc is the cutoff distance required to compute the RDF, which is 2.0 nm in this study. Although the structure-based IBI method can find reasonable solutions to reproduce g(r) using several steps, the CG simulations with the refined potential cannot yield the same density as that of the atomistic MD configuration. To compute the pressure-corrected potential G
no. of CG molecules interactions
f target [×10−5]
density (MD) [g/cm3]
density (CG) [g/cm3]
300
P−P
1.119
0.924
0.925
300
P−N
0.569
0.924
0.925
300
N−N
0.615
0.924
0.925
200 200 400
C3−C3 C2−C2 CO−CO
0.047 0.512 10.566
0.662 0.662 0.752
0.664 0.672 0.759
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
can be found in the Supporting Information. In examining the structural properties relevant to the cis-LCP, no substantial changes were detected in the probability distributions, except for the bending angle (θ) between three adjacent P−N−P beads. Figure 6a shows the CG angle potential energy curves of the P−N−P, which are derived from the probability densities of both trans- and cis-states. Before UV excitation, the energy of the trans-azobenzene is minimized at θ = 180°. The internal energy of the rod-like groups jumped to the level of the excited
considerably small relative errors compared with those computed in other publications.35,39 Generally, the CG pair interactions are constructed using two different methods: deriving the numerically tabulated functions or adjusting the parameters of typical analytic forms such as the Lennard-Jones (LJ) or Morse potentials. Although the latter method is simpler and more practical, using the standardized forms causes difficulties in describing the complex anisotropic interactions between the LC molecules. In particular, it has been demonstrated that the common LJ 12−6 potential is too repulsive to represent the CG interactions between the mesogens.33,35 Accordingly, we selected the elaborately constructed numerical forms to describe the pair interaction between the mesogen beads, which is the most dominant in representing the phase transition behaviors. Meanwhile, less important CG interaction sets are simplified by using the LJ 9−6 functional form: 6 l i σ y9 o o iσ y | o o LJ UCG (r ) = εo m2jjj zzz − 3jjj zzz o } o kr{ o kr{ ~ (9) n where ε and σ are the minimum energy and the separation distance between the particles in the most stable state, respectively. That is, we used the IBI- and PC-based tabulated forms to define the noncovalent interactions between the azobenzene chromophores (P−P, P−N, and N−N) and the interactions between the same types of CG beads. The intermolecular potential energy between the different types of CG beads is simplified using the LJ parameters. The Lorentz− Berthelot mixing rules of εij = εiεj and σij = (σi + σj)/2 are used. The interaction parameters between the same types of beads are listed in Table 3. Although exploiting the hybrid CG Table 3. Refined LJ 9−6 Parameters for the Nonbonded Interaction between the Same Types of the CG Beads interaction types
ε [kcal/mol]
σ [Å]
P−P N−N C3−C3 C2−C2 CO−CO
0.476 0.088 0.310 0.275 0.378
6.55 8.35 5.65 5.15 4.45
potential sets generated a density deviation of ∼2.5% compared with the all-atom model of the LCP network, the LC phase behavior is well reproduced, as given in the Results and Discussion section. 2.3. Derivation of Mesoscale Photoswitching Potential. The photoisomerization from the trans- to cis-state causes not only a dramatic change in the molecular shape but also variations of the intermolecular interaction between the mesogens. Accordingly, the mesoscale bonded and nonbonded potentials are switched to reflect the conversion of the conformational geometry and corresponding mesogen−mesogen interactions, respectively. The difference in the geometry between two isomers is reproduced by capturing the intramolecular probability distribution changes in the atomistic model of the azo-LCP. The UV excitation induces a rotation of the dihedral angle of the diazene group (C−NN−C) from 180° to 0°. Therefore, the cis-LCP network in an AA MD reference is modeled by modifying the dihedral angle potential, which was previously proposed to mimic the photon-assisted energy transition pathway.40 The detailed simulation scheme
Figure 6. (a) CG potential energy of two different isomers of the azobenzene derivatives as a function of the angle between the P−N−P beads and variation of the RDFs (b: P−P; c: P−N pairs) with respect to the cis-fractions. H
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Figure 7. CG configurations of the reference unit cells to construct the nonbonded photoswitching potentials using two different ways.
potential sets treating the interactions between the cisC−C molecules (UCG (r)) and between the different isomers UC−C are required. (r)) (UT−C CG CG (r) is developed by applying the IBI and pressure-correction methods to the 100% cisazobenzene unit cell. As shown in Figure 7a, the UT−C CG (r) is derived from the target RDF and density of the mixture with ncis = 50% (e.g., 50% trans-states and 50% cis-states). Even though the proposed strategy explicitly represents the photoactivation of the individual molecule, the transferability over a wide range of conversion ratios should be verified because UT−C CG (r) is derived from the particular state with ncis = 50%. Therefore, the second approach is developed to devise a photoswitching potential that can sensitively reflect changes in the microstructure with respect to ncis. In this method, as shown in Figure 7b, the same interaction sets are applied to both trans- and cis-azobenzenes. Instead, the nonbonded interaction between the azobenzene beads is gradually changed in terms of the cis-population. AA MD unit cells with different molecular compositions (ncis = 0%, 25%, 50%, 75%, and 100%) are constructed, and then the LJ 9−6 potentials (ULJ CG(r,ncis)) are optimized to satisfy the target RDF and density of each system. The LJ parameters are expressed as a function of ncis using interpolation. To sum up, σ and ε are switched whenever the extent of the photochemical reaction varies. The changes in the LJ parameters and a comparison of the light-induced phase behaviors using two different approaches will be discussed in the Results and Discussion section. 2.4. Details of CG MD Simulation Conditions. In this study, the Mesocite module of the Materials Studio package (BIOVIA, Inc.) is used to map atoms to the CG systems, to apply the IBI and PC schemes, and to prepare the coarsegrained LCP network model. During the ensemble simulations, the temperature and pressure are controlled via the Langevin thermostat41 and the Berendsen barostat,42 respectively. The NVT ensemble at 400 K for 8.0 ns is performed to monitor the RDF and the pressure of the fragment-based references at every iteration step to derive the refined CG pair potentials. Here, the time step of the CG MD runs and the cutoff distance
state by altering the CG conformational energy function. The coefficients of the bending angle energy of the cis-molecules are presented in Table 1. The difference in energy between the grounded and photoexcited states is about 76.15 kcal/mol. This is almost the same as the excitation energy induced by UV 365 nm light irradiation (78 kcal/mol). Notably, the switching of the conformational energy associated with the suitable degree of freedom in the mesoscopic simulation clearly describes the photon-assisted hopping between the ground and excited potential energy surfaces. The excited energy is relaxed and minimized at θ = 96° during the CG MD runs to complete the trans-to-cis conversion. The highly bent cis-isomers play a role as an impurity in liquid crystalline ordering. They gradually weaken the energy of the relatively strong interaction between the LCs, and the associated structural variations directly modify the RDF. Figures 6b,c show a comparison of the RDFs between the azobenzene beads, which are extracted from the AA MD trajectories with different fractions of the cis-molecules. As shown in Figure 6b, the first peak of the RDF between the phenyl beads (P−P) decreases by about 6.7% after UV irradiation. The aromatic rings move away to the local minimum which is located between the first and second peaks. By contrast, as seen in Figure 6c, the pair correlation function between the phenyl and azo groups (P−N) changes slightly during the photochemical reaction. No remarkable variation is observed in the RDF of the azo pairs (N−N). Therefore, we can conclude that the microscopic disordering induced by the photoactivated molecules is mainly triggered by the weakened interaction between the aromatic rings. To reproduce this phenomenon, the CG nonbonded potentials are changed in terms of two different variables: (i) a type of azobenzene isomer and (ii) the fraction of the cismolecules (ncis). In the first method, emphasis is placed on distinguishing two isomerizable groups (see Figure 7a). Before a photochemical reaction occurs, the system only needs the nonbonded interaction between the beads of the trans-isomers T−T (r)). As the cis-isomers are generated, additional (UCG I
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules for the nonbonded interaction are set as 10.0 fs and 2.0 nm, respectively. Mesoscale heating and photoisomerization simulations are performed using the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code, which was developed by the Sandia National Laboratory.43 The cross-linked smectic azo-LCP at 300 K is gradually heated to 420 K in increments of 5 K. Isobaric runs for 15 ns are applied to the unit cells at each temperature, and the corresponding heating rate is about 0.3 K/ns. Further long relaxation simulations for 30 ns are conducted at temperatures within the ranges of the experimental clearing points (the smectic-to-nematic transition temperature (TSN = 355−360 K) and the nematic-to-isotropic transition temperature (TNI = 415−420 K)). Three different morphologies of the CLCP models with the same degree of polymerization are utilized to reduce the uncertainty. During the photomechanical simulations, we modulated the photoisomerization ratio to examine sequential light-induced mesophase behavior and deformations. The CG photoswitching potential sets are applied to particular proportions of the photoreactive groups (ncis = 4%, 8%, 16%, 20%, 30%, 50%, and 70%). For each isomerization ratio, we randomly selected four different groups of azobenzenes to convert them into cis-molecules. The structure is then equilibrated by applying the NPT ensemble for 150 ns. The unit cell dimensions are independently varied during the constantpressure simulations to effectively observe the light-activated microscopic strain. The actual photodynamics and reaction pathway of the azobenzene are more complex, and these have been analyzed by using the first-principle calculations. Although the details of the photochemical reaction kinetics are simplified, the mesoscale simulations can treat quasi-static light-induced conformational changes and establish a relationship between the photochemical reaction and the mechanical deformations of macromolecules.
Figure 8. Thermotropic behaviors and corresponding structural variations of the CG CLCP model: changes in (a) the distribution of the center-of-mass positions of the azobenzene molecules along the director (x-direction) and (b) orientational and translational order parameters under temperature elevation.
3. RESULTS AND DISCUSSION 3.1. Thermotropic Behavior. The thermotropic mesophase behavior of the CG CLCP is investigated by characterizing the LC phase using a calculation of specific order parameters. Figure 8a shows the changes in the distribution of the center-of-mass positions of the azobenzene molecules along the smectic layer normal direction during a heating simulation. Four distinct densely packed regions are located periodically along the director. This means that the CG simulation framework can describe the self-assembly of the LC moieties into the smectic layers at T = 300 K. The director of the azobenzene molecules is not tilted from the layer normal direction, which corresponds to the Sm A phase. The smectic layer spacing is about 36.8 Å, which is in good agreement with that of the experimental results (about 37.0 Å) calculated from the X-ray diffraction pattern.28 As shown in Figure 8a, the height of the peaks decreases and the interlayer density gradually increases as the temperature increases. The smectic layers are broken, and the overall LCP network experiences the phase transition into the nematic phase during the heating simulation. The degree of positional ordering of the mesogens can be quantified with a one-dimensional translational order parameter (τ), which is calculated from the following equation: τ=
i r ·n y expjjj2πi k zzz d { k
k
where ⟨ ⟩ denotes the average over all molecules, rk is the center-of-mass position vector of the kth azobenzene molecule, d is the periodic layer spacing, and n is the director of the liquid crystals, which coincides with a unit vector perpendicular to the smectic plane. Meanwhile, the alignment of the mesogens along the nematic direction is evaluated by calculating the scalar orientational order parameter (s), which can be derived as follows: s=
3 cos2 θ − 1 2
(11)
where θ is the angle between the molecular axis of the moieties and the nematic axis. The molecular axis of the azobenzenes in the CG MD system is defined as a line segment connecting two phenyl (P) beads. Here, we distinguish the three liquid crystalline phases of the azo-containing polymer using these two variables: (i) smectic A (τ > 0, s > 0), (ii) nematic (τ = 0, s > 0), and (iii) isotropic (τ = s = 0). Figure 8b shows the evolutions of the order parameters during temperature changes. At T = 300 K, both the positional and orientational order parameters are above zero (τ = 0.51 and s = 0.67), implying that the LC molecules are arranged in
(10) J
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules the smectic phase. As the temperature increases to 360 K, the LC phase scarcely changes. However, when the system is heated above 360 K, τ primarily decreases, and s decreases slightly to 0.55−0.61. Although self-organized layers gradually vanish because of kinetic thermal energy, the rigid mesogens are still aligned to the long-range direction, and the overall system is transformed into the nematic structure. Finally, the CG LCP structure becomes the isotropic phase without translational and rotational symmetries at a temperature above 410 K. As seen in Figure 8b, the order parameters vary smoothly without discontinuity under phase transitions. This weak first-order transition is consistent with the experimental results of the liquid crystalline elastomers44 and arises from the quenched disorder because of the adjoining chains. The crossover from the smectic to the nematic phase evolves more smoothly than that of the N−I transition. This is attributed to the relatively greater energy required to disrupt the selfassembled layer. The characteristic phase transition temperatures observed in this study (TSN = 360−370 K and TNI = 410−420 K) are comparable to the real clearing points of the CLCP (TSN = 360 K and TNI = 430 K).28 Accordingly, it is concluded that the mesoscale potential sets derived in this study yield not only the structural and thermodynamic properties at the given reference state but also the thermotropic mesophase behavior satisfactorily. Figure 9 clearly shows the thermotropic structural changes in the CG model. The beads representing the trans-azobenzene groups are indicated by a light blue color. 3.2. Photomechanical Behavior. 3.2.1. Light-Induced Structural Changes. Figure 10 shows the variations of characteristic LJ 9−6 nonbonded interaction parameters as the photoisomerization ratio gradually increases. The red dotted curve is plotted by fitting the data points as a nonlinear regression function. For the interaction between the aromatic hydrocarbons (P−P), both σ and ε decline as the photochemical reaction proceeds (Figures 10a,b). Meanwhile, the energy well depth of the CG interaction between the phenyl and azo beads (P−N) increases (Figure 10d). A noteworthy finding is that the energy depth of the interparticle interaction between the P beads decreases to be comparable to that of the interaction between the P and N beads. The trans-molecules exhibit a stronger noncovalent interaction between the unsaturated organic groups (ε = 0.476 kcal/mol) than that between the phenyl rings and the azo group (ε = 0.263 kcal/ mol). However, the photon-assisted generation of the kinked molecules significantly weakens the interaction energy of P−P, and it becomes indistinguishable from the attraction energy of P−N. The light-induced structural changes simulated by using two different CG nonbonded interaction energy sets, which were discussed earlier, are compared. Figure 11 shows the phase transition behavior of the CG azo-containing polymer triggered by gradual photoisomerization. Similar to the thermotropic behavior, the molecular shape change of the photosensitive groups induces a sequential Sm A−N−I phase transition. At ncis = 20%, τ decreases to 0.28−0.31 and s decreases to 0.42−0.49, which implies the Sm A−N transition. By increasing ncis up to 50%, the entire LCP network becomes isotropic without any symmetric ordering. Figure 12 shows the configurations of the azo-LCP network with different cisfractions. Only 4% of the bent-shaped azobenzenes do not severely influence the smectic structures. However, the CG unit cell with ncis = 20% shows that the smectic layers are
Figure 9. Coarse-grained configurations of the azo-containing LCP network models at (a) T = 300 K (Sm A), (b) T = 370 K (N), and (c) T = 420 K (I).
primarily disassembled with the retained alignment of the LCs. As seen in Figure 12c, the ordering of the mesogens is clearly broken when the concentration of the cis-molecules increases up to 50%. Another remarkable feature of the results is the similar trends in the variations of both the translational and orientational orders obtained by implementing the two different CG potential sets (Figure 11). Accordingly, we can claim that UT−C CG (r) calculated by targeting the RDF of the unit cell with ncis = 50% can describe both the local corruption of the mesogenic order and the phase behaviors with a wide range K
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
Figure 10. CG LJ 9−6 pair interaction parameters (σ and ε) for the beads associated with the azobenzene molecules as a function of the photoisomerization ratio: (a, b) P−P interaction; (c, d) P−N interaction.
direction as the photoisomerization ratio gradually increases. At low concentrations (ncis = 4−16%), the unit cell length along the initial alignment direction increases, and the nominal tensile strains are about 1.8−6.8%. On the other hand, the azopolymer exhibits a monotonically increasing compressive strain from 4.6% to 20.4% under higher cis-population conditions (ncis = 30−70%). As shown in Figure 11, the light-induced stretching at low ncis dominantly originates from the Sm A−N transition. The uniaxial shrinkage is induced by the decrease in the orientational order during the N−I transition. The distinguished shape change of the smectic solid is consistent with that of the experimental findings.28,45 As shown in Table 4, the compressive photostrain of the all-atom CLCP at ncis = 20% is larger than that of the CG model because the AA MD model only exhibits the contraction behavior under the destabilization of the nematic phase. However, under higher ncis (50−70%), the strain of the mesoscale model has a good agreement with that of the atomistic reference. Contrary to the switching of the deformation modes of the Sm A polymer, the nematic CLCP with s = 0.55 and τ = 0.09 at T = 390 K exhibits only compressive behavior regardless of the fraction of the cis-state (see Table 4). In particular, the nematic polymer network at higher temperature can show extremely large shrinkage up to 18.8% with the introduction of a small amount of photoactivated molecules. This is in good agreement with observations that the light-induced N−I transition of the monodomain nematic polymer doped with only 9−18 mol % of the photosensitive compounds results in
Figure 11. Changes in the translational and orientational order parameters of the azo-containing LCP as a function of the photoisomerization ratio (ncis). The data points are obtained by implementing two different CG photoswitching potential sets at T = 300 K.
of ncis. All of the following results are generated by implementing these mesoscale photoswitching potential sets. 3.2.2. LC Phase-Dependent Photomechanical Deformations. The liquid crystalline phase transition and structural changes upon UV irradiation are accompanied by the shape change of the cross-linked polymers. Table 4 shows the photostrains of the CG and AA MD models along the LC director axis with different proportions of the photoactivated molecules. The sequential photomechanical responses of the Sm A CLCP at T = 300 K show a change in the deformation L
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
elements and thermal kinetic energy provide a sufficient internal mobility for the photoactivated molecules to disorder the surrounding LCs. The compressive strain dramatically increases up to 37.6% with 70% of the cis-azobenzenes, and it is about twice that of the Sm A solids at T = 300 K. The attractive feature demonstrated in this study is that the photomechanical deformation modes can be manipulated by adjusting the operation temperature and corresponding LC structures. The uniaxial expansion during the Sm A−N transition is attributed to the macromolecular shape change of the crosslinked polymer network. Actually, since the side-chain mesogens are assembled to form the layered structure, the backbone polymer chains are compactly packed between the layers. The growth of the smectic structure orients the main chains perpendicular to the director axis, and the entire network becomes more oblate. Hence, the photoinduced breakdown of the smectic layers rather increases the shape anisotropy of the backbones. To quantify the conformational change during the photoisomerization, a shape parameter (r) is calculated, which effectively captures the degree of the shape anisotropy of the macromolecules. First, the radius of the gyration tensor components in terms of the Cartesian coordinates (α and β) are measured using the following equation: N
Gαβ =
1 ∑ (ri ,α − R α)(ri ,β − Rβ) N i=1
(12)
Here, ri,α represents the position coordinate of bead i with respect to the α-axis, Rα is the position of the backbone center of mass, and N is the number of beads. Then, r is defined as the square of the averaged inertia radii ratio, which can be expressed as follows: r=
RP R⊥
2
=
2
Gxx ( Gyy +
Gzz )/2
(13)
where RΡ and R⊥ are the radii parallel and perpendicular to the nematic direction, respectively. Figure 13 shows the change of r under the photomechanical responses. In the initial Sm A state, the shape parameter of the backbone is less than unity, which indicates an oblate polymeric conformation. Under the Sm A−N transition (ncis = 4−16%), the shape anisotropy along the initially oriented direction increases by up to 23.6%.
Figure 12. Sequential light-induced structural changes of the LCP network model at T = 300 K with ncis = (a) 4%, (b) 20%, and (c) 50%.
Table 4. Comparison of the Photostrains of the CG and AA MD Models photostrain (%) ncis (%) 4 8 16 20 30 50 70
T = 300 K (CG) 6.84 4.51 1.79 −0.48 −4.56 −12.11 −20.38
± ± ± ± ± ± ±
1.21 1.22 1.19 1.55 1.40 1.50 1.62
T = 390 K (CG) −18.05 −19.07 −24.09 −24.41 −31.27 −35.42 −37.63
± ± ± ± ± ± ±
3.08 3.04 3.83 3.48 3.84 4.89 5.13
T = 300 K (AA)
−4.21 ± 1.98 −11.65 ± 2.43 −17.10 ± 3.52
considerable deformations.4 Despite the low concentrations of the cis-state, the kinked molecules slowly but steadily collapse the molecular ordering and induce a corresponding shape change. Another reason for the amplified mechanical actions at higher temperature is because the increases of free volume
Figure 13. Change in the backbone shape parameters of the CG CLCP model at T = 300 K. M
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules However, after the first transition is completed, the photonassisted structural change does not generate any elongation of the main chains. Thus, we can conclude that the unique photomechanical deformation of the Sm A CLCP is independent of the decrease in the nematic orientational order but is mainly affected by the removal of the positional symmetry of the LCs. 3.2.3. Relation between Macromolecular Shape and LC Order Parameters. To predict the large-scale behavior of the photoresponsive polymer, it is essential to reveal the coupling mechanism between the evolution of the microscopic LC structure in response to the photon energy and corresponding macromolecular shape change. In particular, the shape parameter successfully reflects the effect of photoisomerization on the mechanical response of the cross-linked polymer network. The optomechanical constitutive equation of the finite element analysis has previously been described by the neoclassical elastic free energy, which is expressed as a function of the shape parameter of the polymer.18,34 In this study, the mesoscale simulation framework is utilized to link the LC symmetry and shape parameters in terms of various external conditions (temperature and photoisomerization ratio). Notably, it is useful to parametrize the photomechanical response as a function of ncis because the continuum scale simulations require a deformation profile in terms of the effective light intensity and corresponding concentration of the cis-molecules at each material point. Figure 14 shows the light-induced disordering behavior of the LCP network with different external temperature conditions. The figure shows that similar to the response at T = 300 K, the positional symmetry primarily decreases, and then the rotational symmetry disappears. However, the operating temperature clearly affects the photoinduced phase and microstructural evolutions. As seen in Figure 14a, while the translational order parameter is almost the same as the initial LC state of T = 300 K, it diminishes even more rapidly at T = 330 K. While τ decreases to 0.2 with ncis = 20−30% at room temperature, only about 10% of the cis-molecules are required to reach the same value at higher temperature. This tendency becomes more pronounced for both translational and orientational order parameters with increasing temperature. This indicates that thermally induced variations in the macromolecular morphology render the overall system more sensitive to the light illumination. The proportional disruption of the LC symmetry with an increase in ncis and its temperature dependency are analogous to the classical dilution model,4,46 which assumes that the characteristic phase transition temperatures (TSN and TNI) decrease linearly with an increase in the number of cis-molecules. In the dilution model, every increment of the photoisomerization ratio is equivalent to an increment of the temperature because both phenomena result in the decline of the order parameters. The dilution model utilizes only one curve of the evolution of the order parameters (Figure 8b) by simply mapping the effect of the light illumination to that of the temperature increase. However, Figure 14 shows that the multiple τ values form a pair with the same s by changing the operating temperature. This indicates that the dependence of the photoinduced phase behavior on the external condition cannot be simply described by only one temperature-order parameter curve. Indeed, the presented simulation is capable of reflecting the complex photoresponsive microstructures that cannot be captured by the classical model.
Figure 14. Changes in the translational and orientational order parameters with respect to the variation of the photoisomerization ratio at (a) T = 330 K, (b) T = 360 K, and (c) T = 390 K.
Generally, the shape parameter of the nematic polymer has a linear relationship with the uniaxial orientational order parameter, which can be expressed as r = 1 + αs. The coefficient α is equal to 3 for freely jointed polymers.46 In the case of the smectic CLCP, an additional term that represents the effect of τ should be considered. As shown in Table 4, while the corruption of the nematic ordering results in the contraction of the polymer conformation along the director, the disassembly of the smectic layers causes the elongation. In addition, the major cause of the expansion during the Sm−N transition is the increase in the shape anisotropy of the backbone chains. As seen in Figure 13, in the initial state, the increase of the shape parameter is very dramatic, with only a small amount of decrease in τ. On the other hand, as τ decreases, the variation becomes progressively smaller. Accordingly, we define the increment of r induced by the reduction of the positional symmetry as a function of N
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
simulation results. α, a linear coefficient for the orientational order, ranges from 2.98 to 3.42, which is very close to that of the simple nematic LCP (α = 3). Both α and η become larger with an increase in temperature, which implies that more photomechanical deformations can be generated by the same extent of variations in the LC structures at higher temperature. Hence, the calculated coefficients in this study can successfully describe the modification of the polymeric morphology by increasing the external temperature and the corresponding responsiveness of the photomechanical behavior. Figure 15 shows the light-triggered evolutions of r at each temperature. The data points (black squares) are obtained by averaging the shape anisotropy of the overall CG CLCP molecules using eq 13. They show good agreement with the fitted curve (red dashed lines) predicted using eq 14. At the Sm A phase (Figures 15a,b), r increases at low ncis, and then the expansion mode switches to contractive deformation as the photochemical reaction proceeds. On the other hand, Figures 15c,d show that photoinduced elongation behavior cannot be observed at a temperature equal to or greater than the smectic−nematic transition temperature. Yu et al.28 reported that when the temperature was lower than 90 °C, the smectic polymer film bent away from the light source because of the dominant expansion at the surface, and the bending direction then switched. As the temperature reached above 90 °C, the polymer film only bent toward the light source because of the N−I transition-induced shrinkage. Our multiscale- and multi-
exponential terms. The final expression of r(ncis,T) in terms of the photoisomerization ratio and the temperature is r(ncis , T ) = 1 + α(T )s(ncis , T ) − [exp{η(T )τ(ncis , T )} − 1] (14)
At τ = 0, the effect of the smectic layer disappears, and the equation becomes r = 1 + αs, which is the formula for the nematic polymer. As τ gradually increases, the main chain aligns perpendicular to the oriented direction, which reduces r. The shape parameter changes exponentially in response to the variation of τ. Here, α and η are the coefficients that are related to the material properties and are expressed as a function of T to consider the mechanical response of the macromolecular network at each operation temperature. Table 5 lists the coefficient values obtained by fitting r in terms of the order parameters computed from the CG MD Table 5. Fitted Coefficient Values of Eq 14 Which Define a Coupled Relation between Photoinduced Phase Transition and Mechanical Deformation Behaviors temperature [K]
α(T)
η(T)
300 330 360 390
2.977 2.969 3.302 3.415
0.909 1.190 1.391
Figure 15. Evolutions of the shape parameter of the CG CLCP molecules with varying the extent of the photoisomerization at (a) T = 300 K, (b) T = 330 K, (c) T = 360 K, and (d) T = 390 K. O
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules
cis-population, which is proportional to local light intensity. The shape parameter can be utilized to perform the continuum scale optomechanical constitutive equation to accurately predict the large-scale deformations. Thus, we plan to extend the multiscale framework to a continuum scale analysis by considering the polymeric morphology parameters such as the cross-link density and spacer length between the backbone and mesogens. We therefore expect that the outcome presented here can be used to identify the microstructure-dependent photomechanical effect and precisely design the macroscale deformations of the photoresponsive smart materials.
physics-integrated simulation results are expected to be utilized to analyze this type of complex 3D photomechanical deformation with diverse LC phases.
4. CONCLUSIONS In this work, a comprehensive multiscale simulation was conducted to examine the thermo- and photoresponsive behaviors of the azo-containing polymer with diverse LC phases. In particular, we describe the complex structure of the LC polymer network and the global motions triggered by photochemical reaction by using the CG simulation techniques. The number of degrees of freedom is efficiently reduced by mapping the repeating units into the CG beads to enlarge the time and length scales of the conventional MD simulation. The elaborately developed intramolecular and intermolecular potentials of the CG system are derived by primarily targeting the polymeric conformation and the radial distribution functions of the AA MD reference model, respectively. In addition, the CG photoswitching potential sets which can exhibit mesoscopic polymeric structural changes driven by the quantum event are constructed for the first time. The light-induced manipulations of the CLCP are classified according to the individual nature of the photoisomerization: (i) the molecular shape change during the trans-to-cis conversion and (ii) the role of the activated molecules as an impurity to break the LC symmetry. Especially, the azobenzene unit is represented by three connective beads to enable molecular bending under illumination, which cannot be realized by using anisotropic Gay-Berne particles. The disordering effect of the kinked molecules on the surrounding mesogens is reflected by incorporating the nonbonded interaction between the trans- and cis-states. The pair interaction is calibrated by applying the iterative Boltzmann inversion and pressure correction methods to the atomistic unit cell in which two different microstates coexist. The potential sets bridging the quantum-mechanical event and the intermolecular interaction are not only capable of capturing the local distortion of the ordering generated by individual cisazobenzenes but also transferrable to various photoisomerization ratios. The mesoscale simulations can successfully predict thermotropic and photoinduced transitions between the three specific LC phases (Sm A, N, and I). The photomechanical deformations induced by the phase behaviors are also investigated in detail. In particular, two different deformation modes appear: the expansion during the Sm A−N transition and the contraction during the N−I transition. The unusual elongation behavior is attributed to the transformation of the main chain into a more prolate conformation as the smectic layers are broken. The coupled relationship between the detailed photoinduced changes of the LC state and the mechanical deformations is identified by using the analytic equation, which defines the macromolecular shape parameter as a function of the LC order parameters. The developed model satisfactorily reproduces the effects of the extent of the isomerization and external temperature on the photomechanical response of the CLCP network. The presented mesoscale description of the photoisomerization can be extensively applied to explore not only the photomechanics of the azo-polymer but also the mechanism of the molecular switch that is usually implemented into biomacromolecules to modify their functionality. In addition, the molecular shape anisotropy is examined as a function of the
■
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.8b02535. Construction of all-atom molecular dynamics reference model and atomistic description of the photoisomerization (PDF) CG model with 4% of the cis-molecules (AVI) CG model with 20% of the cis-molecules (AVI) CG model with 50% of the cis-molecules (AVI)
■
AUTHOR INFORMATION
Corresponding Author
*(M.C.) Tel +82-02-880-1693, Fax +82-02-886-1693, e-mail
[email protected]. ORCID
Maenghyo Cho: 0000-0003-3942-9261 Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS This work was supported by a grant from the National Research Foundation of Korea (NRF) funded by the Korea government (MSIP) (Grant 2012R1A3A2048841).
■
REFERENCES
(1) Jiang, H.; Kelch, S.; Lendlein, A. Polymers move in response to light. Adv. Mater. 2006, 18 (11), 1471−1475. (2) Chatani, S.; Kloxin, C. J.; Bowman, C. N. The power of light in polymer science: photochemical processes to manipulate polymer formation, structure, and properties. Polym. Chem. 2014, 5 (7), 2187− 2201. (3) Finkelmann, H.; Nishikawa, E.; Pereira, G.; Warner, M. A new opto-mechanical effect in solids. Phys. Rev. Lett. 2001, 87 (1), 015501. (4) Hogan, P.; Tajbakhsh, A.; Terentjev, E. UV manipulation of order and macroscopic shape in nematic elastomers. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 2002, 65 (4), 041720. (5) Volgraf, M.; Gorostiza, P.; Numano, R.; Kramer, R. H.; Isacoff, E. Y.; Trauner, D. Allosteric control of an ionotropic glutamate receptor with an optical switch. Nat. Chem. Biol. 2006, 2 (1), 47. (6) Ikeda, T.; Tsutsumi, O. Optical switching and image storage by means of azobenzene liquid-crystal films. Science 1995, 268 (5219), 1873−1875. (7) Fukuda, T.; Matsuda, H.; Shiraga, T.; Kimura, T.; Kato, M.; Viswanathan, N. K.; Kumar, J.; Tripathy, S. K. Photofabrication of surface relief grating on films of azobenzene polymer with different dye functionalization. Macromolecules 2000, 33 (11), 4220−4225. (8) Zhang, F.; Zarrine-Afsar, A.; Al-Abdul-Wahid, M. S.; Prosser, R. S.; Davidson, A. R.; Woolley, G. A. Structure-based approach to the photocontrol of protein folding. J. Am. Chem. Soc. 2009, 131 (6), 2283−2289.
P
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX
Article
Macromolecules (9) McConney, M. E.; Martinez, A.; Tondiglia, V. P.; Lee, K. M.; Langley, D.; Smalyukh, I. I.; White, T. J. Topography from topology: photoinduced surface features generated in liquid crystal polymer networks. Adv. Mater. 2013, 25 (41), 5880−5885. (10) Yu, Y.; Nakano, M.; Ikeda, T. Photomechanics: directed bending of a polymer film by light. Nature 2003, 425 (6954), 145. (11) Serak, S.; Tabiryan, N.; Vergara, R.; White, T. J.; Vaia, R. A.; Bunning, T. J. Liquid crystalline polymer cantilever oscillators fueled by light. Soft Matter 2010, 6 (4), 779−783. (12) Iamsaard, S.; Aßhoff, S. J.; Matt, B.; Kudernac, T.; Cornelissen, J. J.; Fletcher, S. P.; Katsonis, N. Conversion of light into macroscopic helical motion. Nat. Chem. 2014, 6 (3), 229. (13) Kowalski, B. A.; Guin, T. C.; Auguste, A. D.; Godman, N. P.; White, T. J. Pixelated polymers: directed self assembly of liquid crystalline polymer networks. ACS Macro Lett. 2017, 6, 436. (14) White, T. J.; Broer, D. J. Programmable and adaptive mechanics with liquid crystal polymer networks and elastomers. Nat. Mater. 2015, 14 (11), 1087. (15) Natansohn, A.; Rochon, P. Photoinduced motions in azocontaining polymers. Chem. Rev. 2002, 102 (11), 4139−4176. (16) Yamada, M.; Kondo, M.; Mamiya, J. i.; Yu, Y.; Kinoshita, M.; Barrett, C. J.; Ikeda, T. Photomobile polymer materials: towards lightdriven plastic motors. Angew. Chem. 2008, 120 (27), 5064−5066. (17) Palagi, S.; Mark, A. G.; Reigh, S. Y.; Melde, K.; Qiu, T.; Zeng, H.; Parmeggiani, C.; Martella, D.; Sanchez-Castillo, A.; Kapernaum, N.; et al. Structured light enables biomimetic swimming and versatile locomotion of photoresponsive soft microrobots. Nat. Mater. 2016, 15 (6), 647. (18) Lin, Y.; Jin, L.; Huo, Y. Quasi-soft opto-mechanical behavior of photochromic liquid crystal elastomer: Linearized stress−strain relations and finite element simulations. Int. J. Solids Struct. 2012, 49 (18), 2668−2680. (19) Choi, J.; Chung, H.; Yun, J.-H.; Cho, M. Photo-isomerization effect of the azobenzene chain on the opto-mechanical behavior of nematic polymer: A molecular dynamics study. Appl. Phys. Lett. 2014, 105 (22), 221906. (20) Böckmann, M.; Doltsinis, N. L.; Marx, D. Nonadiabatic hybrid quantum and molecular mechanic simulations of azobenzene photoswitching in bulk liquid environment. J. Phys. Chem. A 2010, 114 (2), 745−754. (21) Chung, H.; Choi, J.; Yun, J.-H.; Cho, M. Nonlinear photomechanics of nematic networks: upscaling microscopic behaviour to macroscopic deformation. Sci. Rep. 2016, 6, 20026. (22) Wilson, M. R. Progress in computer simulations of liquid crystals. Int. Rev. Phys. Chem. 2005, 24 (3−4), 421−455. (23) Berardi, R.; Muccioli, L.; Zannoni, C. Can nematic transitions be predicted by atomistic simulations? A computational study of the odd−even effect. ChemPhysChem 2004, 5 (1), 104−111. (24) Sánchez-Ferrer, A.; Finkelmann, H. Opto-mechanical effect in photoactive nematic main-chain liquid-crystalline elastomers. Soft Matter 2013, 9 (18), 4621−4627. (25) Garcia-Amorós, J.; Finkelmann, H.; Velasco, D. Influence of the photo-active azo cross-linker spacer on the opto-mechanics of polysiloxane elastomer actuators. J. Mater. Chem. 2011, 21 (4), 1094−1101. (26) Yu, Y.; Maeda, T.; Mamiya, J. i.; Ikeda, T. Photomechanical effects of ferroelectric liquid-crystalline elastomers containing azobenzene chromophores. Angew. Chem. 2007, 119 (6), 899−901. (27) Lansac, Y.; Glaser, M. A.; Clark, N. A.; Lavrentovich, O. D. Photocontrolled nanophase segregation in a liquid-crystal solvent. Nature 1999, 398 (6722), 54. (28) Zhang, Y.; Xu, J.; Cheng, F.; Yin, R.; Yen, C.-C.; Yu, Y. Photoinduced bending behavior of crosslinked liquid-crystalline polymer films with a long spacer. J. Mater. Chem. 2010, 20 (34), 7123−7130. (29) Zannoni, C. Molecular design and computer simulations of novel mesophases. J. Mater. Chem. 2001, 11 (11), 2637−2646.
(30) Wilson, M. R. Molecular dynamics simulations of flexible liquid crystal molecules using a Gay-Berne/Lennard-Jones model. J. Chem. Phys. 1997, 107 (20), 8654−8663. (31) Moon, J.; Choi, J.; Cho, M. Opto-mechanical behavior and interfacial characteristics of crosslinked liquid crystalline polymer composites with carbon nanotube fillers. Carbon 2017, 121, 181−192. (32) Ilnytskyi, J. M.; Neher, D.; Saphiannikova, M. Opposite photoinduced deformations in azobenzene-containing polymers with different molecular architecture: molecular dynamics study. J. Chem. Phys. 2011, 135 (4), 044901. (33) Peter, C.; Delle Site, L.; Kremer, K. Classical simulations from the atomistic to the mesoscale and back: coarse graining an azobenzene liquid crystal. Soft Matter 2008, 4 (4), 859−869. (34) Chung, H.; Park, J.; Cho, M. Numerical study of light-induced phase behavior of smectic solids. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 2016, 94 (4), 042707. (35) Zhang, J.; Su, J.; Ma, Y.; Guo, H. Coarse-grained molecular dynamics simulations of the phase behavior of the 4-cyano-4′pentylbiphenyl liquid crystal system. J. Phys. Chem. B 2012, 116 (7), 2075−2089. (36) Mukherjee, B.; Delle Site, L.; Kremer, K.; Peter, C. Derivation of coarse grained models for multiscale simulation of liquid crystalline phase transitions. J. Phys. Chem. B 2012, 116 (29), 8474−8484. (37) Milano, G.; Goudeau, S.; Müller-Plathe, F. Multicentered Gaussian-based potentials for coarse-grained polymer simulations: Linking atomistic and mesoscopic scales. J. Polym. Sci., Part B: Polym. Phys. 2005, 43 (8), 871−885. (38) Reith, D.; Pütz, M.; Müller-Plathe, F. Deriving effective mesoscale potentials from atomistic simulations. J. Comput. Chem. 2003, 24 (13), 1624−1636. (39) Ju, S.-P.; Huang, S.-C.; Lin, K.-H.; Chen, H.-Y.; Shen, T.-K. Prediction of Optical and Dielectric Properties of 4-Cyano-4pentylbiphenyl Liquid Crystals by Molecular Dynamics Simulation, Coarse-Grained Dynamics Simulation, and Density Functional Theory Calculation. J. Phys. Chem. C 2016, 120 (26), 14277−14288. (40) Heinz, H.; Vaia, R.; Koerner, H.; Farmer, B. Photoisomerization of azobenzene grafted to layered silicates: simulation and experimental challenges. Chem. Mater. 2008, 20 (20), 6444− 6456. (41) Schneider, T.; Stoll, E. Molecular-dynamics study of a threedimensional one-component model for distortive phase transitions. Phys. Rev. B: Condens. Matter Mater. Phys. 1978, 17 (3), 1302. (42) Berendsen, H. J.; Postma, J. v.; van Gunsteren, W. F.; DiNola, A.; Haak, J. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81 (8), 3684−3690. (43) Plimpton, S.; Crozier, P.; Thompson, A. LAMMPS-large-scale atomic/molecular massively parallel simulator. Sandia National Laboratories 2007, 18, 43. (44) Selinger, J. V.; Jeon, H. G.; Ratna, B. Isotropic-nematic transition in liquid-crystalline elastomers. Phys. Rev. Lett. 2002, 89 (22), 225701. (45) Aßfalg, N.; Finkelmann, H. A smectic A liquid single crystal elastomer (LSCE): Phase behavior and mechanical anisotropy. Macromol. Chem. Phys. 2001, 202 (6), 794−800. (46) Warner, M.; Terentjev, E. M. Liquid Crystal Elastomers; OUP Oxford: 2003; Vol. 120.
Q
DOI: 10.1021/acs.macromol.8b02535 Macromolecules XXXX, XXX, XXX−XXX