Simulation of Stimuli-responsive Polymer Networks

The structure and material properties of polymer networks can depend sensitively on changes in the environment. There is a great deal of progress in the development of stimuli-responsive hydrogels for applications like sensors, self-repairing materials or actuators. Biocompatible, smart hydrogels can be used for applications, such as controlled drug delivery and release, or for artificial muscles. Numerical studies have been performed on different length scales and levels of details. Macroscopic theories that describe the network systems with the help of continuous fields are suited to study effects like the stimuli-induced deformation of hydrogels on large scales. In this article, we discuss various macroscopic approaches and describe, in more detail, our phase field model, which allows the calculation of the hydrogel dynamics with the help of a free energy that considers physical and chemical impacts. On a mesoscopic level, polymer systems can be modeled with the help of the self-consistent field theory, which includes the interactions, connectivity, and the entropy of the polymer chains, and does not depend on constitutive equations. We present our recent extension of the method that allows the study of the formation of nano domains in reversibly crosslinked block copolymer networks. Molecular simulations of polymer networks allow the investigation of the behavior of specific systems on a microscopic scale. As an example for microscopic modeling of stimuli sensitive polymer networks, we present our Monte Carlo simulations of a filament network system with crosslinkers.


Introduction
If you put a gummy bear into a glass of water, the next day you will find the same bear grown in volume by more than a factor of five.(The effect depends on the choice of Gummy bears.)This is a very simple example of a polymer network that is strongly responsive to the stimuli of the environment (in this case humidity).Different kinds of smart polymer networks are sensitive to different stimuli, such as temperature changes [1,2], illumination [3,4], or properties of the solvent, such as the pH value [1,[5][6][7], the ion concentration [7,8], or the chemical composition [5,[8][9][10].Some materials are responsive to more than one kind of stimuli [11,12].Reversible networks can be switched back and forth between different states by turning the stimulus on and off.Smart polymer networks can be used for various applications, such as chemical sensors [8,13], biosensors [13,14], actuators [7,15], artificial muscles [5], or for drug transport and release [1,6,16].For the development and optimization of these devices, it is crucial to investigate the network behavior, both experimentally and theoretically.In this article, we focus on numerical methods for studying intelligent polymer networks and describe modeling techniques on different scales, ranging from macroscopic methods to molecular modeling.
There are two important types of smart polymer networks: (a) Covalently crosslinked networks that swell or shrink if exposed to a stimulus, and (b) networks with a variable topology, in which crosslinks can be broken or created by external stimuli.
Both types of smart polymer networks can change their density, their mechanical properties, and their permeability in response to changes of the environment.The majority of applications use chemically crosslinked polymer networks (see the application references in the first paragraph), however, networks with stimuli-sensitive crosslinks are also used for the development of sensors [9] and for controlled drug delivery and release [5,17,18].One important aspect of reversibly crosslinked polymer networks is their ability to rebuild broken bonds and to rearrange in a completely different topology and shape.This way, the polymer network can heal if it has been plastically deformed or broken [19].On the other hand, a fixed network topology has its own advantages: Even strongly swollen, the gummy bear still keeps its shape.
This article is dedicated to numerical studies of smart polymer networks that are sensitive to external stimuli.Modeling methods are discussed on different levels of details, from macroscopic theories to molecular modeling.The introduction to the different techniques is combined with more detailed descriptions of three numerical methods that have been developed by the authors.We think that they are particularly interesting examples of numerical studies of intelligent polymer networks that show novel routes for polymer network modeling [20][21][22][23][24][25].We use these examples to describe important aspects of modeling methods on different levels of detail.For a comprehensive description of our numerical studies, we refer to the original articles [20][21][22][23][24][25].
We start with macroscopic models, in which all material properties are described by continuous fields.Here, we discuss various modeling principles, before we explain, in more detail, a phase field model for chemo-sensitive hydrogels that was recently presented by one of the authors [20].Afterwards, we give a short introduction to self-consistent field theory.The starting point of the self-consistent field theory is the exact statistically physical description of a polymer model system so that constitutive equations are not required.It has been used, for many decades, for the studying of polymer melts.We elucidate an extension of the method that allows modeling reversibly crosslinked polymer networks with switchable microphases [ and present, as an example, a Monte Carlo study of a network formed by stiff polymers and physical crosslinkers that, depending on external paramet properties [24,25].

Macroscopic Models, Finite Element
Stimulated by external conditions solvent, the volume of hydrogels can change strongly.For the investigation of hydrogel swelling on a macroscopic scale, various macroscopic models have been developed.Overviews on macroscopic hydrogel models have been presented in The first model, developed by Flory polymer network.There are various extensions of this statistical approach [ state of a hydrogel in a given environment can be calculated by minimizing the free energy system.Relevant contributions to the free energy are the interaction and the entropy of polymer chains, the entropy of ions in the hydrogel and hydrogel models have been used, for exam The equilibrium ion distribution in the hydrogel and in the bath can be obtained from the Donnan equilibrium and the postulation of electric neutrality, as discussed below.
An important model for hydrogels is the multiphase mixture model, in which the system consists of two or more phases.In the biphasic theory, the system is divided into an incompressible polymer phase and an incompressible fluid phase that is in contact with a fluid reservoir [ extended to three-phase (or multiphase) mixture models that include a polymer phase, a solvent phase, and one or several ionic phases, formed by the mobile ions [ the theory of porous media, has been developed by Freiboth Macroscopic hydrogel theories are frequently based on conservation laws and constitutive equations [39,40,41,43].If the ion, solvent, and polymer phases are incompressible, mass conservation leads to: where and are the volume fractions of the polymer network and the solvent phase, while are the volume fractions of the corresponding concentrations .In the absence of chemical reactions, the conce conservation law: where is the µth component of the local velocity of phase convention for Greek indices that appear twice.
Frequently, electric neutrality is postulated for the hydrogel charged ions are attached to the polymer network with a concentration networks with switchable microphases [21][22][23].Finally, we discuss molecular modeling techniques and present, as an example, a Monte Carlo study of a network formed by stiff polymers and physical crosslinkers that, depending on external parameters, can change its structural and topological lement and Phase Field Methods Stimulated by external conditions, such as a change of the pH or the chemical composition of the hydrogels can change strongly.For the investigation of hydrogel swelling on a macroscopic scale, various macroscopic models have been developed.Overviews on macroscopic hydrogel models have been presented in [26][27][28].
The first model, developed by Flory and Rehner [29,30], is based on the entropy changes in a swelling polymer network.There are various extensions of this statistical approach [31 state of a hydrogel in a given environment can be calculated by minimizing the free energy system.Relevant contributions to the free energy are the interaction and the entropy of polymer chains, the entropy of ions in the hydrogel and in the bath, and the Coulomb interactions.Thermodynamic hydrogel models have been used, for example, by Caykara et al. [35], and Kramarenko The equilibrium ion distribution in the hydrogel and in the bath can be obtained from the Donnan equilibrium and the postulation of electric neutrality, as discussed below.
drogels is the multiphase mixture model, in which the system consists of two or more phases.In the biphasic theory, the system is divided into an incompressible polymer phase and an incompressible fluid phase that is in contact with a fluid reservoir [37] phase (or multiphase) mixture models that include a polymer phase, a solvent phase, and one or several ionic phases, formed by the mobile ions [38][39][40][41].A multiphase approach, based on been developed by Freiboth et al. [42].

Macroscopic hydrogel theories are frequently based on conservation laws and constitutive
].If the ion, solvent, and polymer phases are incompressible, mass conservation are the volume fractions of the polymer network and the solvent phase, while are the volume fractions of the ion phases.The volume fractions .In the absence of chemical reactions, the conce component of the local velocity of phase .In this article, we use the summation convention for Greek indices that appear twice.
Frequently, electric neutrality is postulated for the hydrogel system [41,43 charged ions are attached to the polymer network with a concentration , electric neutrality requires ].Finally, we discuss molecular modeling techniques and present, as an example, a Monte Carlo study of a network formed by stiff polymers and physical ers, can change its structural and topological a change of the pH or the chemical composition of the hydrogels can change strongly.For the investigation of hydrogel swelling on a macroscopic scale, various macroscopic models have been developed.Overviews on macroscopic and Rehner [29,30], is based on the entropy changes in a swelling [31][32][33][34].The equilibrium state of a hydrogel in a given environment can be calculated by minimizing the free energy of the total system.Relevant contributions to the free energy are the interaction and the entropy of polymer chains, and the Coulomb interactions.Thermodynamic Kramarenko et al. [36].The equilibrium ion distribution in the hydrogel and in the bath can be obtained from the Donnan drogels is the multiphase mixture model, in which the system consists of two or more phases.In the biphasic theory, the system is divided into an incompressible polymer phase ].The method has been phase (or multiphase) mixture models that include a polymer phase, a solvent Here, and are the valence number and the concentration of fixed ions, attached to the polymer network and is the valence number of the mobile ion species the system is studied with the help of momentum equations [ phases, the forces induced by the chemical potential of solvent potentials of the ions, are balanced by the friction forces between the different phases that flow with different velocities.The equations are of the form with chemical/electrochemical potentials is the mass of a solvent molecule ( coefficients between phases and where are the components of the stress tensor.Furthermore, several co used: Usually, the chemical potential is expressed by a linear combination of the ion concentrations and the hydrostatic pressure, while the electrochemical potentials of the ions can be expressed as a function of the chemical activity term, which is proportional to and a term that is proportional to the electric potential [ between the stress tensor of the polymer network and its elastic strain tensor If the ion concentrations in the solvent are small, so that volume fractions of the mobile ions .Then, the volume fraction of water molecules is deter fraction and does not have to be considered explicitly.
The flow of ions can be treated more generally by using the Poisson equation instead of the electric neutrality condition.Models that include the Poisson equation have be Wallmersperger et al. [47], Li et al hydrogels induced by chemical and electric stimuli.
Following the method used in [ with the Poisson equation.The time dependence of the concentration be described by the Nernst-Planck equation: with the Faraday constant temperature .Ions of species have diffusion constants to the continuity equation Equation concentration gradient and the electric potential non-ideal mixing behavior of the solvent, extended versions of Eq include a term for the chemical activity or other source terms for the ions [ are the valence number and the concentration of fixed ions, attached to the polymer s the valence number of the mobile ion species .In several models, the dynamics of the system is studied with the help of momentum equations [37][38][39][40][41]43,46].For the solvent and the ion phases, the forces induced by the chemical potential of solvent molecules, and the electrochemical are balanced by the friction forces between the different phases that flow with different velocities.The equations are of the form: with chemical/electrochemical potentials .The mass densities are given by is the mass of a solvent molecule ( ) or an ion ( ).The constants and .For the polymer phase, one has: are the components of the stress tensor.Furthermore, several constitutive equations are used: Usually, the chemical potential is expressed by a linear combination of the ion concentrations and the hydrostatic pressure, while the electrochemical potentials of the ions can be expressed as a vity term, which is proportional to with an activity coefficient and a term that is proportional to the electric potential [39].Furthermore, a linear or nonlinear relation between the stress tensor of the polymer network and its elastic strain tensor is considered [ If the ion concentrations in the solvent are small, so that and are distinctly larger than the volume fractions of the mobile ions ( ), Equation (1) can be approximated by .Then, the volume fraction of water molecules is determined by the polymer volume fraction and does not have to be considered explicitly.
The flow of ions can be treated more generally by using the Poisson equation instead of the electric neutrality condition.Models that include the Poisson equation have been used, for example, by et al. [48], and Chen et al. [39], to study the swelling and shrinking of hydrogels induced by chemical and electric stimuli.
Following the method used in [48], here we discuss a typical set of equations used in combination with the Poisson equation.The time dependence of the concentration of mobile ions of type Planck equation: , the gas constant have diffusion constants .The Nernst-Planck equation corresponds uation (2) under the assumption that the ion flow is determined by the concentration gradient and the electric potential .For systems that show chemical reactions or a ideal mixing behavior of the solvent, extended versions of Equation ( 5) include a term for the chemical activity or other source terms for the ions [49 (3) are the valence number and the concentration of fixed ions, attached to the polymer .In several models, the dynamics of 6].For the solvent and the ion and the electrochemical are balanced by the friction forces between the different phases that flow with (4) are given by , where ).The constants are friction nstitutive equations are used: Usually, the chemical potential is expressed by a linear combination of the ion concentrations and the hydrostatic pressure, while the electrochemical potentials of the ions can be expressed as a with an activity coefficient , ].Furthermore, a linear or nonlinear relation is considered [28].are distinctly larger than the can be approximated by mined by the polymer volume The flow of ions can be treated more generally by using the Poisson equation instead of the electric en used, for example, by to study the swelling and shrinking of ions used in combination of mobile ions of type can (5) , and the Planck equation corresponds under the assumption that the ion flow is determined by the show chemical reactions or a have been used, which 49].In [39], a convection term , with a flow velocity potential obeys the Poisson equation with the vacuum dielectric constant depends on the differences respective concentrations at the external boundary of the hydrogel The osmotic pressure causes a deformation non-linear deformations can be considered by finite strain theory [ small, they fulfill the linear equation where and are components of the full pressure tensor in the hydrogel region and the unit tensor, respectively.The coefficients hydrogel.The Green-Lagrangian strain tensor For small deformations, the third term in the brackets can be neglected.Solving Eq numerically, requires appropriate bo the corresponding surface region is fixed.For the remaining surface region with a surface normal vector pointing outward and a surface traction vector Various types of stimuli-induced swelling and deformation of hydrogels have been investigated with transport theories.A numerical study of temperature For hydrogel spheres in water, the radial polymer density dist of swelling time [38].Cheng et al.

non-uniform fixed charge density [
A system of particular interest, which is frequently investigated experiment hydrogel strip in water between two electrodes.The electric field from the electrodes is perpendicular to the strip area.The strip is fixated at one end electric field, leads to a bending of the non which can be used for sensors and actuators in biological and other aqueous environments reproduced and analyzed in detail in various numerical studies [ Equations ( 6)-( 9) are typically solved with finite element methods.Alternatively, meshless methods have been used [39,40,50,[53][54][55][56], which are conceptually more complex but avoid the definition and adaptation of a mesh grid during the simulation.In finite element and meshless methods, boundary conditions have to be defined on the hydrogel surface with a flow velocity , is added to the right side of Equation obeys the Poisson equation: with the vacuum dielectric constant and the relative dielectric constant between the ion concentrations in the hydrogel and the at the external boundary of the hydrogel: The osmotic pressure causes a deformation of the initial geometry of the hydrogel.considered by finite strain theory [45,46,48].If the deformations are small, they fulfill the linear equation: are components of the full pressure tensor in the hydrogel region and the unit tensor, respectively.The coefficients and are the first and second Lamé coefficients of the Lagrangian strain tensor is given by: For small deformations, the third term in the brackets can be neglected.Solving Eq numerically, requires appropriate boundary conditions.If the hydrogel is firmly attached to a substrate, the corresponding surface region is fixed.For the remaining surface region, pointing outward and a surface traction vector induced swelling and deformation of hydrogels have been investigated with transport theories.A numerical study of temperature-dependent swelling is presented in [ For hydrogel spheres in water, the radial polymer density distribution has been calculated as a function et al. have studied the curling effect of articular cartilage with uniform fixed charge density [46].A system of particular interest, which is frequently investigated experiment hydrogel strip in water between two electrodes.The electric field from the electrodes is perpendicular to the strip area.The strip is fixated at one end, or in the middle.The diffusion of ions ds to a bending of the non-uniformly swelling strip.The deformation of the strip can be used for sensors and actuators in biological and other aqueous environments reproduced and analyzed in detail in various numerical studies [39,43,52].
are typically solved with finite element methods.Alternatively, meshless methods ], which are conceptually more complex but avoid the definition and adaptation of a mesh grid during the simulation.In finite element and meshless methods, boundary conditions have to be defined on the hydrogel surface, which typically moves at each time Equation (5).The electric (6) .The osmotic pressure in the hydrogel and the For small deformations, the third term in the brackets can be neglected.Solving Equations ( 6)-( 9), undary conditions.If the hydrogel is firmly attached to a substrate, , one has: (10) .induced swelling and deformation of hydrogels have been investigated dependent swelling is presented in [50].ribution has been calculated as a function have studied the curling effect of articular cartilage with A system of particular interest, which is frequently investigated experimentally [51], consists of a hydrogel strip in water between two electrodes.The electric field from the electrodes is perpendicular or in the middle.The diffusion of ions, induced by the uniformly swelling strip.The deformation of the strip, can be used for sensors and actuators in biological and other aqueous environments, could be are typically solved with finite element methods.Alternatively, meshless methods ], which are conceptually more complex but avoid the definition and adaptation of a mesh grid during the simulation.In finite element and meshless methods, boundary which typically moves at each time step so that the hydrogel surface has to be traced during the swelling process.Recently, we have developed a phase field model, which avoids the precise localization of the phase boundaries and can be easily coupled with other phase fields that may descri

Phase Field Theory
The phase field method (PF method in the following) is generally used to study multiphase systems with interfacial dynamics.It allows following the dynamics of the system interfaces inherently co to all relevant physical and chemical driving fields.In recent years, there has been a growing number of articles that apply the PF method on polymer systems, polymer melt or a polymer solution (see, fo that has been developed in our group [20] the concept of the PF method and describe how it has been applied on polymer systems.Our motivation for the model's development was to introduce a new modeling concept that can be easily adapted to different network properties and external conditions.
The PF method, which is used to describe the non based on a free energy functional of the form: where is the free energy density of the bulk considers spatial inhomogeneities standard Ginzburg-Landau free energy with coefficients and .An overview over phase field methods and in [62].Usually, one has an extended region in spa another region in which the phase is absent and region of finite width, in which the phase field decreases smoothly.If the phase field is described by a Cahn Here, the flux density is given by with a variational derivative polymer blend, is often taken as the concentration of one component and the bulk free energy density is the Margules activity from the Flory noise term is added to the Cahn the hydrogel surface has to be traced during the swelling process.Recently, we have developed a phase field model, which avoids the precise localization of the phase boundaries and can be easily coupled with other phase fields that may describe additional phases in the system.
The phase field method (PF method in the following) is generally used to study multiphase systems with interfacial dynamics.It allows following the dynamics of the system interfaces inherently co to all relevant physical and chemical driving fields.In recent years, there has been a growing number of articles that apply the PF method on polymer systems, i.e., systems in which at least one phase is a polymer melt or a polymer solution (see, for example, [57][58][59][60][61]).A phase field model for a hydrogel that has been developed in our group [20] is elucidated in detail below.To that end, we first introduce the concept of the PF method and describe how it has been applied on polymer systems.Our ation for the model's development was to introduce a new modeling concept that can be easily adapted to different network properties and external conditions.
The PF method, which is used to describe the non-equilibrium dynamics of phase boundaries, is .If there is only one phase field variable is the free energy density of the bulk and usually a double well potential.The second part considers spatial inhomogeneities.A simple example of the free energy in Landau free energy: . An overview over phase field methods and their ].Usually, one has an extended region in space in which the phase field variable is another region in which the phase is absent and .Both regions are separated by an interface region of finite width, in which the phase field decreases smoothly.If is conserved, the dynamics of eld is described by a Cahn-Hilliard equation [63]: Here, the flux density is given by: , and a mobility , which, in general, may depend on is often taken as the concentration of one component and the bulk free energy is the Margules activity from the Flory-Huggins theory [64,65].In many cases, a Gaussian is added to the Cahn-Hilliard equation so that the equation has the hydrogel surface has to be traced during the swelling process.Recently, we have developed a phase field model, which avoids the precise localization of the phase boundaries and can be easily coupled The phase field method (PF method in the following) is generally used to study multiphase systems with interfacial dynamics.It allows following the dynamics of the system interfaces inherently coupled to all relevant physical and chemical driving fields.In recent years, there has been a growing number systems in which at least one phase is a ]).A phase field model for a hydrogel is elucidated in detail below.To that end, we first introduce the concept of the PF method and describe how it has been applied on polymer systems.Our ation for the model's development was to introduce a new modeling concept that can be easily equilibrium dynamics of phase boundaries, is , the free energy is (11) usually a double well potential.The second part .A simple example of the free energy in Equation ( 11) is the (12) their applications are given ce in which the phase field variable is and .Both regions are separated by an interface is conserved, the dynamics of , which, in general, may depend on .In a fluid is often taken as the concentration of one component and the bulk free energy ].In many cases, a Gaussian Hilliard equation so that the equation has a Langevin form: (15) Considering a flow velocity equation with a noise term [66]: The fluid dynamics can be described by a continuity equation and a (advanced) Navier equation, which is given by: where is the mass density, is the hydrostatic pressure, considers the surface tension from the interfaces [ gravity, while leading to viscoelastic behavior [ model [68], which has been used, viscoelastic droplets [69] and the coalescence of polymer drops and interfaces [ The connection of the Cahn-Hilliard and the Navier nomenclature of Hohenberg and Halperin [ polymer phase is described as a Newtonian fluid.In principle, all terms on the right hand side (beside ) can be coupled to the phase field variable the composition.The model H has been used to study spinodal decomposition with and without external velocity fields and Rayleigh slightly different Cahn-Hilliard free energy has been used to study the shape of fluid films on a dewetting substrate [73].
The dynamics of non-conserved phase field variables equation [74]: where is the mobility coefficient and term.In various studies, the Allen polymer crystals [60,75,76].In the following example, it is used to characterize the extension of the hydrogel [20].

Phase Field Model of a Hydrogel
Recently, we have developed a PF model for studying the swelling behavior of a hydrogel as the ion concentration in the surrounding is changed described to demonstrate the concept the hydrogel is characterized by a phase field variable Equation ( 18), with a free energy [20] , Equation ( 15) can be extended to the convective Cahn The fluid dynamics can be described by a continuity equation and a (advanced) Navier is the hydrostatic pressure, is the viscous stress tensor considers the surface tension from the interfaces [59].The term represents external forces denotes the elastic response of polymers to a flow field, leading to viscoelastic behavior [67].An alternative description of viscoelasticity is the Oldroy , together with the PF method, to study the viscoelastic droplets [69] and the coalescence of polymer drops and interfaces [ Hilliard and the Navier-Stokes equation is called "Model H nomenclature of Hohenberg and Halperin [71].In many cases, the part polymer phase is described as a Newtonian fluid.In principle, all terms on the right hand side (beside ) can be coupled to the phase field variable , as the involved material parameters may depend on the composition.The model H has been used to study spinodal decomposition with and without external velocity fields and Rayleigh-Taylor instabilities [57,59,72].The described method with Hilliard free energy has been used to study the shape of fluid films on a conserved phase field variables is described by the Allen is the mobility coefficient and consists again of a bulk free energy density and a gradient term.In various studies, the Allen-Cahn equation has been used to investigate the growth of 76].In the following example, it is used to characterize the extension of the

Hydrogel
Recently, we have developed a PF model for studying the swelling behavior of a hydrogel as the ion concentration in the surrounding is changed [20].In the following, some details are described to demonstrate the concept.The system is sketched in Figure the hydrogel is characterized by a phase field variable .The dynamics of energy [20] 49 can be extended to the convective Cahn-Hilliard (16) The fluid dynamics can be described by a continuity equation and a (advanced) Navier-Stokes (17) is the viscous stress tensor, and represents external forces, such as denotes the elastic response of polymers to a flow field, of viscoelasticity is the Oldroyd-B to study the elongation and burst of viscoelastic droplets [69] and the coalescence of polymer drops and interfaces [61,70].
Stokes equation is called "Model H" in the is neglected so that the polymer phase is described as a Newtonian fluid.In principle, all terms on the right hand side (beside the involved material parameters may depend on the composition.The model H has been used to study spinodal decomposition with and without ].The described method with Hilliard free energy has been used to study the shape of fluid films on a is described by the Allen-Cahn (18) nsists again of a bulk free energy density and a gradient Cahn equation has been used to investigate the growth of 76].In the following example, it is used to characterize the extension of the Recently, we have developed a PF model for studying the swelling behavior of a hydrogel [20].In the following, some details are ure 1.The extension of .The dynamics of is controlled by (19) The second term is a double well potential as used in Landau theory.The second and third term provide an interface region of finite width. of the mixture, respectively.In alloy systems, the coefficients the surface energy and the interface width.For the function , and , so that the derivative of the region.The term is the osmotic pressure at the initial equilibrium configuration from which the displacement measured.The osmotic pressure The elastic energy density: includes the difference of the first Lamé coefficients in the hydrogel region bath region .
Figure 1.Sketch of a hydrogel attached to the network, and mobile cations and anions.In the phase field method, all system components are represented by continuous field variables.The extension of the hydrogel is characterized by a phase field variable region (light green area).At the interface to the bath, the phase field decreases smoothly to [20].
The bath solution is assumed to be a shear mechanical equilibrium for the strain tensor With defined by: The second term is a double well potential as used in Landau theory.The second and third term provide an interface region of finite width.and are the potential height and the molar volume of the mixture, respectively.In alloy systems, the coefficients and the surface energy and the interface width.For the function , so that the derivative of the first part of is only finite in the interface is the difference between the osmotic pressure and is the osmotic pressure at the initial equilibrium configuration from which the displacement is defined as [20]: includes the difference of the first Lamé coefficients in the hydrogel region Sketch of a hydrogel system with a polymer network, a surrounding bath, ions and mobile cations and anions.In the phase field method, all system components are represented by continuous field variables.The extension of the by a phase field variable, which is inside the hydrogel region (light green area).At the interface to the bath, the phase field decreases smoothly to The bath solution is assumed to be a shear-free liquid.Analogously to Equation mechanical equilibrium for the strain tensor by an equation: and .The concentration of fixed charges is 50 The second term is a double well potential as used in Landau theory.The second and third term are the potential height and the molar volume can be extracted from one has , is only finite in the interface is the difference between the osmotic pressure and , where is the osmotic pressure at the initial equilibrium configuration from which the displacement is (20) (21) and in the surrounding system with a polymer network, a surrounding bath, ions and mobile cations and anions.In the phase field method, all system components are represented by continuous field variables.The extension of the inside the hydrogel region (light green area).At the interface to the bath, the phase field decreases smoothly to Equation (8), we ensure (22) .The concentration of fixed charges is where and are the concentration of fixed ions and the phase field variable at the initial equilibrium state.Note that of fixed charges is zero inside the bath region while it is a finite constant inside the hyd The fields , , , , and are coupled via In the simulation, a hydrogel with negative fixed ions ( Na + ( ) and Cl − ( hydrogel is .In the bath, the Na .The initial values of Donnan equilibrium and electric neutrality conditions At the boundary of the system, the displacement and the electric potential is kept the concentration of mobile ions is continuously changed with time.The time step is chosen to be .In Figure 2, profiles of the phase field (Figure 2a) the gel is strongly swollen and shows a rather quadratic shape, due to the quadratic boundary of the simulation cell.At large ion concentration region becomes more circular in order to minimize the interface energy.are the concentration of fixed ions and the phase field variable at the initial at all times.In the bath region, the concentration of fixed charges is zero inside the bath region while it is a finite constant inside the hyd are coupled via Equations ( 19)-( 23  Our new phase field model is an alternative approach for studying the dynamics of hydrogel swelling induced by the ion concentration of the bath.With the chosen decrease or increase rate for the boundary ion concentration, the hydrogel swells and shrinks without significant hysteresis.The method goes without explicit boundary conditions at the hydrogel and can be easily adapted to other hydrogel systems.For example, one can study a hydrogel shell that includes an active agent.Then, a concentration field of drug molecules can be added to the given model in order to study the drug release kineti

Self-Consistent Field Theory
The macroscopic theories of polymer networks depend on phenomenological, constitutive equations.Microscopic aspects of the polymer system can be studied with molecular simulations but they are restricted to small time and length scales.An alternative method for polymer modeling is the self-consistent field theory.Here, polymers are represented by continuous lines.They interact by steric and attractive interactions, where the latter is weighted by the Flory parameter.account the connectivity of polymer chains, the fact that two spatially separated monomers are much more correlated if they belong to the same polymer.For the given model system, the theory starts from  [20].In general, the equilibrium gel volume is monotonously shrinking with increasing ion concentration in the bath.The result is in qualitative agreement with experiments, in which the pH is kept constant and the ion concentration is varied [ Our new phase field model is an alternative approach for studying the dynamics of hydrogel duced by the ion concentration of the bath.With the chosen decrease or increase rate for the boundary ion concentration, the hydrogel swells and shrinks without significant hysteresis.The method goes without explicit boundary conditions at the hydrogel and can be easily adapted to other hydrogel systems.For example, one can study a hydrogel shell that includes an active agent.Then, a concentration field of drug molecules can be added to the given model in order to study the drug release kinetics.
The macroscopic theories of polymer networks depend on phenomenological, constitutive equations.Microscopic aspects of the polymer system can be studied with molecular simulations but and length scales.An alternative method for polymer modeling is the consistent field theory.Here, polymers are represented by continuous lines.They interact by steric and attractive interactions, where the latter is weighted by the Flory parameter.account the connectivity of polymer chains, the fact that two spatially separated monomers are much more correlated if they belong to the same polymer.For the given model system, the theory starts from ].In general, the equilibrium gel volume is monotonously ation in the bath.The result is in qualitative agreement with experiments, in which the pH is kept constant and the ion concentration is varied [35,77].
Our new phase field model is an alternative approach for studying the dynamics of hydrogel duced by the ion concentration of the bath.With the chosen decrease or increase rate for the boundary ion concentration, the hydrogel swells and shrinks without significant hysteresis.The method goes without explicit boundary conditions at the hydrogel surface and can be easily adapted to other hydrogel systems.For example, one can study a hydrogel shell that includes an active agent.Then, a concentration field of drug molecules can be added to the given The macroscopic theories of polymer networks depend on phenomenological, constitutive equations.Microscopic aspects of the polymer system can be studied with molecular simulations but and length scales.An alternative method for polymer modeling is the consistent field theory.Here, polymers are represented by continuous lines.They interact by steric and attractive interactions, where the latter is weighted by the Flory parameter.The method takes into account the connectivity of polymer chains, the fact that two spatially separated monomers are much more correlated if they belong to the same polymer.For the given model system, the theory starts from an exact partition function.Then, the equilibrium properties are typically determined with the help of a saddle-point approximation.The self and Helfand [78,79], is a well-established method for polymer chains without Polymer networks with permanent random crosslinks theory [85][86][87], but the application of the method in numerical studies is difficult, since the fields are defined in the replica space.Recently, which allows studying smart polymer networks with reversible crosslinks [ application of the method for a blend of reversibly crosslinked and B homopolymers, is described.The monomers that the correspnding homopolymers would demix macroscopically in the absence of the copolymers.In the AB + A + crosslinks, the system has been studied in detail theoretically [ Many aspects of the model can be transferred to other systems with amphiphilic dissolved in a mixture of incompatib non-covalent crosslinks, as in [92].
The model system resembles an incompressible polymer blend with a constant overall density and three types of polymers monomers and a length , (ii) polymers of type , which are symmetric diblock copolymers system stoichiometry is determined by th and the relative A fraction of homopolymers We use paths with All chains have a Gaussian shape distribution that can be divided into the density from the A part of the copolymers.Analogously, one has B monomers.
Without crosslinks, the partition function of the with path integrals over all configurations of the useful mathematical tools.They can be considered as the integral that considers all configurations L/N.In practice, most path integrals cannot be solved numerically.One exception is the path integral of a single Gaussian chain.The basic idea of the self such that only path integrals of single Gaussian chains have to be solved.In interactions between polymers are considered by hen, the equilibrium properties are typically determined with the help of a point approximation.The self-consistent field theory, which has been introduced by Edwards established method for polymer chains without Polymer networks with permanent random crosslinks have been studied with a replica field ], but the application of the method in numerical studies is difficult, since the fields are defined in the replica space.Recently, we have developed an extended selfwhich allows studying smart polymer networks with reversible crosslinks [21application of the method for a blend of reversibly crosslinked ABdiblock copolymers is described.The monomers A and B are assumed to have a low compatibility so that the correspnding homopolymers would demix macroscopically in the absence of the + B mixture, various nanostructures may form.Without reversible the system has been studied in detail theoretically [81,88,89] and experimentally [ Many aspects of the model can be transferred to other systems with amphiphilic dissolved in a mixture of incompatible solvents.We assume that the AB copolymers can bind with ].The model system resembles an incompressible polymer blend with a constant overall density .The system consists of (i) homopolymers of type homopolymers of type b with B monomers and a length which are symmetric diblock copolymers AB that have a length of system stoichiometry is determined by the volume fraction of copolymers: fraction of homopolymers: to parametrize the configuration of the All chains have a Gaussian shape distribution .A and B monomers have a local density that can be divided into the density from type homopolymers and the density part of the copolymers.Analogously, one has Without crosslinks, the partition function of the system is given by: over all configurations of the th polymer of type useful mathematical tools.They can be considered as the limit for integral that considers all configurations of a polymer with fixed length L and .In practice, most path integrals cannot be solved numerically.One exception is the path integral of a single Gaussian chain.The basic idea of the self-consistent field theory is to transform such that only path integrals of single Gaussian chains have to be solved.In Eq interactions between polymers are considered by: 53 hen, the equilibrium properties are typically determined with the help of a consistent field theory, which has been introduced by Edwards established method for polymer chains without crosslinks [80][81][82][83][84].
been studied with a replica field ], but the application of the method in numerical studies is difficult, since the fields are -consistent field theory, -23].In Section 2.1, an diblock copolymers, mixed with A are assumed to have a low compatibility so that the correspnding homopolymers would demix macroscopically in the absence of the AB orm.Without reversible ] and experimentally [90,91].Many aspects of the model can be transferred to other systems with amphiphilic AB copolymers copolymers can bind with The model system resembles an incompressible polymer blend with a constant overall density homopolymers of type with A monomers and a length , (iii) that have a length of .The (25) (26) to parametrize the configuration of the th polymer of type .monomers have a local density homopolymers and the density for the density of (27) th polymer of type .Path integrals are of a 3N-dimensional and N atoms of diameter .In practice, most path integrals cannot be solved numerically.One exception is the path integral of a consistent field theory is to transform Equation ( 27) Equation ( 27), non-bonded The Flory-Huggins interaction in the exponent includes a Flory parameter considers the incompressibility of the polymer blend.So far, the ansatz corresponds to that of Edwards and Helfand [78,79].We have developed an extended model, which considers the formation of reversible crosslinks.The result is the first self The method is described in detail in [21].In the following, we explain the basic concept.

Self-Consistent Field Theory for Reversibly Crosslinked Polymer Networks
We apply our extended self-consistent field theory on an reversible crosslinks between the copolymers.We assume that crosslinks between adjacent monomer pairs A + A, A + B, or B + B form with a probability proportional to respectively.Here, is the crosslink strength and types of the crosslinked monomers.

Setting
, we use the crosslink strength to specify the crosslink properties.For convenience, the crosslink strength is represented by the parameter .With the help of the Hubbard-Stratonovich (HS) transformation a free energy can be derived, which can be solved with a saddle-point approximation.We present the resulting free energy expression with dimensionless quantities, indicated by a tilde.Dimensionless lengths are defined as gyration radius of a free Gaussian polymer of length , while we use dimensionless free energies of the form system with dimensions.We use field free energy can be divided into thr with a Flory-Huggins part: a crosslink contribution: and a free energy term: which considers the entropy of the incompressible polymer blend for given monomer densities.The auxiliary fields and transformation.The quantity fluctuating in the auxiliary fields (see [ Huggins interaction in the exponent includes a Flory parameter considers the incompressibility of the polymer blend.So far, the ansatz corresponds to that of Edwards 79].We have developed an extended model, which considers the formation of reversible first self-consistent field theory for reversibly crosslinked polymer networks.The method is described in detail in [21].In the following, we explain the basic concept.

Consistent Field Theory for Reversibly Crosslinked Polymer Networks
consistent field theory on an AB + A + B polymer blend that may form reversible crosslinks between the copolymers.We assume that crosslinks between adjacent monomer form with a probability proportional to is the crosslink strength and and are weighting factors, considering the of the crosslinked monomers.The partition function can then be written as [ , we use the crosslink strength and the crosslink asymmetry to specify the crosslink properties.For convenience, the crosslink strength is represented by the Stratonovich (HS) transformation a free energy can be derived, which point approximation.We present the resulting free energy expression with dimensionless quantities, indicated by a tilde.Dimensionless lengths are defined as of a free Gaussian polymer of length .Dimensionless d , while we use dimensionless free energies of the form dimensions.We use for the crosslink strength.The dimensionless mean field free energy can be divided into three terms [21]: which considers the entropy of the incompressible polymer blend for given average monomer densities.
with and are introduced in the HS is the partition function of a single polymer chain of type fluctuating in the auxiliary fields (see [20]).It can be calculated numerically by solving modified 54 (28) Huggins interaction in the exponent includes a Flory parameter .The delta function considers the incompressibility of the polymer blend.So far, the ansatz corresponds to that of Edwards 79].We have developed an extended model, which considers the formation of reversible crosslinked polymer networks.The method is described in detail in [21].In the following, we explain the basic concept.polymer blend that may form reversible crosslinks between the copolymers.We assume that crosslinks between adjacent monomer , , or , are weighting factors, considering the partition function can then be written as [21][22][23]: (29) and the crosslink asymmetry to specify the crosslink properties.For convenience, the crosslink strength is represented by the Stratonovich (HS) transformation a free energy can be derived, which point approximation.We present the resulting free energy expression with dimensionless quantities, indicated by a tilde.Dimensionless lengths are defined as with the .Dimensionless densities are given by for a for the crosslink strength.The dimensionless mean (30) (31) (32) (33) which considers the entropy of the incompressible polymer blend for given average monomer densities.are introduced in the HS is the partition function of a single polymer chain of type ]).It can be calculated numerically by solving modified diffusion equations.In the saddle p the desired monomer densities energy with respect to , The method has been used to calculate a phase diagram as a function of the relative adhesion strength , the crosslink asymmetry investigated system has homopolymers of lengths and a volume fraction can form a lamellar phase or a hexagonal phase monomer densitities in Figure 4.The phase diagram has been determined by study box length ratio of .The geometry supports the formation of a hexagonal as well as a lamellar phase.We have varied lowest free energy density.We find that for stable.Otherwise, the hexagonal phase is stable if phase boundaries between the hexagonal and the lamellar phase are shown for various the homopolymers in the range of The self-consistent field theory study shows that the reversibly crosslinked polymer network in an blend switches from the hexagonal to the lamellar phase and back if the stoichiometry of the system or the crosslink properties are changed.The latter switching method can be achieved with crosslinks that depend on external parameters illumination [3].Note that in the hexagonal structure, the minority pha domains, while in the lamellar phase the stripes are extended on a large length scale.This may be used for drug release applications.diffusion equations.In the saddle point approximation, the auxiliary fields can be determined by considering that the variations of the free , and must vanish.An algorithm is described in [ s been used to calculate a phase diagram as a function of the relative adhesion , the crosslink asymmetry , and the relative fraction of homopolymers The investigated system has homopolymers of lengths with and a volume fraction of copolymers.In the studied parameter range, the system can form a lamellar phase or a hexagonal phase, as shown in the sketches and the calculated 4. The phase diagram has been determined by study .The geometry supports the formation of a hexagonal as well as a in the range of and determined the phase with the We find that for and , only lam stable.Otherwise, the hexagonal phase is stable if and are suitably small.In Fig phase boundaries between the hexagonal and the lamellar phase are shown for various .consistent field theory study shows that the reversibly crosslinked polymer network in an blend switches from the hexagonal phase to the lamellar phase and back if of the system or the crosslink properties are changed.The latter switching method can achieved with crosslinks that depend on external parameters, such as Note that in the hexagonal structure, the minority phase is separated into small domains, while in the lamellar phase the stripes are extended on a large length scale.This may be used and, finally, can be determined by considering that the variations of the free must vanish.An algorithm is described in [21].s been used to calculate a phase diagram as a function of the relative adhesion fraction of homopolymers The , a Flory parameter of of copolymers.In the studied parameter range, the system as shown in the sketches and the calculated 4. The phase diagram has been determined by studying unit cells with a .The geometry supports the formation of a hexagonal as well as a and determined the phase with the , only lamellar structures are are suitably small.In Figure 5, the phase boundaries between the hexagonal and the lamellar phase are shown for various fractions of consistent field theory study shows that the reversibly crosslinked polymer network in an blend switches from the hexagonal phase to the lamellar phase and back if of the system or the crosslink properties are changed.The latter switching method can such as the temperature [2] or se is separated into small domains, while in the lamellar phase the stripes are extended on a large length scale.This may be used

Monte Carlo and Molecular D
Properties of polymer networks can be studied on a microscopic level with the help of simulations and molecular dynamic simulations, in which molecular details of the system can be investigated explicitly.Detailed descriptions of molecular simulation techniques can be found in [ and [94].Special simulation methods for polym Carlo (MC) simulations, small changes of the molecule configurations are performed randomly at each step and the new configuration is accepted with a probability depending on the energy change.The acceptance criterion produces a sequence of configurations with a Boltzmann distribution.Other ensembles, such as a grand-canonical ensemble or a constant allowing changes of system parameters like the number of molecules or molecular dynamics (MD) simulations, the positions and velocities of the atoms are obtained from equations of motion that depend on external and interatomic forces.The method allows direct investigation of the dynamical prope dynamic simulations that include a random noise term and diffusive particle dynamic (DPD)

Dynamic Simulations of Polymer Networks
Properties of polymer networks can be studied on a microscopic level with the help of simulations and molecular dynamic simulations, in which molecular details of the system can be investigated explicitly.Detailed descriptions of molecular simulation techniques can be found in [ ]. Special simulation methods for polymer systems are described in [95 Carlo (MC) simulations, small changes of the molecule configurations are performed randomly at each step and the new configuration is accepted with a probability depending on the energy change.The tance criterion produces a sequence of configurations with a Boltzmann distribution.Other canonical ensemble or a constant-pressure ensemble can be sampled by changes of system parameters like the number of molecules or the volume of the system.In molecular dynamics (MD) simulations, the positions and velocities of the atoms are obtained from equations of motion that depend on external and interatomic forces.The method allows direct investigation of the dynamical properties of the system.Similar simulation techniques are Brownian dynamic simulations that include a random noise term and diffusive particle dynamic (DPD) 56 polymer blend: The colored curves denote the ) and the lamellar (dark blue), (cyan) [22].Phase behavior of consistent field theory study, Journal of 10.1557/jmr.2013.315.
Properties of polymer networks can be studied on a microscopic level with the help of Monte Carlo simulations and molecular dynamic simulations, in which molecular details of the system can be investigated explicitly.Detailed descriptions of molecular simulation techniques can be found in [93] in [95].In canonical Monte Carlo (MC) simulations, small changes of the molecule configurations are performed randomly at each step and the new configuration is accepted with a probability depending on the energy change.The tance criterion produces a sequence of configurations with a Boltzmann distribution.Other pressure ensemble can be sampled by the volume of the system.In molecular dynamics (MD) simulations, the positions and velocities of the atoms are obtained from equations of motion that depend on external and interatomic forces.The method allows direct rties of the system.Similar simulation techniques are Brownian dynamic simulations that include a random noise term and diffusive particle dynamic (DPD) simulations, in which the dynamics of interpenetrating spheres resemble the flow behavior of different species in the system (see [93,94 colloids from stimuli-responsive microgel capsules [ studied with Brownian dynamics [97].Compared to field details but are restricted to much smaller time and length scales.In multiscale simulations, material parameters obtained in molecular simulations are used in field on larger scales.
Atomistic simulations of polymers are very time consuming.Consequently, only a few atomistic simulations of polymer networks have been performed [98 reactions in epoxy networks, atomistic simulations simulation effort can be reduced by using a coarse In this case, a polymer is approximated by a chain of beads, where each bead represents a monomer or a sequence of neighboring monomers.A model potential for a polymer chain is often composed of a bond length potential , which depends on the distance of neighboring beads, a bending potential , which considers the angle between neighboring bond vectors torsion of the polymer chain.Then, the potential energy of a linear polymer with where is the distance vector between the centers of mass of neighbor atoms angle between and , while .This type of potentials has been used in Gibbs ensemble MD simulations o polymer networks [104].In a more coarse that is distinctly longer than the persistence length of the polymer.In this case, the angular dependent parts of the polymer potential can be neglected.The remaining bond length potential is typically modeled by a harmonic spring or a The resulting chains are rather flexible and can be used to represent suitably long polymers on the respective length scale.If the persistence length is comparable with the polymer length, which is often the case for short biopolymers, a wormlike chain bond-length and the bending potential but neglects the torsion [ polymers, one can also neglect the bending and bond polymer by a rigid rod.While this approach neglects many internal degrees of freedom, it enables very efficient sampling of the configuration space.simulations, using the example of our filament network simulations

Monte Carlo Simulation of Filament
We have used this method by studying a system of rigid filaments that can form filament networks with the help of reversibly binding crosslinkers [ short actin biopolymers and myosin crosslinkers.The system can form various structures dissolved actin filaments, separated bundles of actin filaments bundle structures.In living cells, the filament network is a dynamic system.Using ATP, the actin simulations, in which the dynamics of interpenetrating spheres resemble the flow behavior of different 94]).DPD simulations have been used to investigate the release of responsive microgel capsules [96].Polyelectrolytes in electric fields have been dynamics [97].Compared to field theories, molecular simulations provide more to much smaller time and length scales.In multiscale simulations, material parameters obtained in molecular simulations are used in field-based or other coarse Atomistic simulations of polymers are very time consuming.Consequently, only a few atomistic simulations of polymer networks have been performed [98][99][100][101].In a multiscale simulation of curing reactions in epoxy networks, atomistic simulations have been combined with DPD simulations [102].The simulation effort can be reduced by using a coarse-grained model, such as the united atom model [103].In this case, a polymer is approximated by a chain of beads, where each bead represents a monomer or a sequence of neighboring monomers.A model potential for a polymer chain is often composed of a , which depends on the distance of neighboring beads, a bending potential , which considers the angle between neighboring bond vectors, and a potential torsion of the polymer chain.Then, the potential energy of a linear polymer with is the distance vector between the centers of mass of neighbor atoms , while is the angle between the vectors .This type of potentials has been used in Gibbs ensemble MD simulations o ].In a more coarse-grained approach, each bead can represent a chain sequence ly longer than the persistence length of the polymer.In this case, the angular dependent parts of the polymer potential can be neglected.The remaining bond length potential is typically modeled by a harmonic spring or a finitely extensible nonlinear elastic (FENE) The resulting chains are rather flexible and can be used to represent suitably long polymers on the respective length scale.If the persistence length is comparable with the polymer length, which is often the case for short biopolymers, a wormlike chain model is frequently used, which includes the length and the bending potential but neglects the torsion [108,109].In the case of very stiff polymers, one can also neglect the bending and bond-length fluctuations and represent the whole rigid rod.While this approach neglects many internal degrees of freedom, it enables very efficient sampling of the configuration space.We elucidate the method and facilities of rigid rod our filament network simulations [24,25].

ilament Networks with Reversibly Binding Crosslinkers
We have used this method by studying a system of rigid filaments that can form filament networks with the help of reversibly binding crosslinkers [24].One example for such systems are mixtures of short actin biopolymers and myosin crosslinkers.The system can form various structures dissolved actin filaments, separated bundles of actin filaments, and filament networks with and without ving cells, the filament network is a dynamic system.Using ATP, the actin .In a multiscale simulation of curing have been combined with DPD simulations [102].The the united atom model [103].In this case, a polymer is approximated by a chain of beads, where each bead represents a monomer or a sequence of neighboring monomers.A model potential for a polymer chain is often composed of a , which depends on the distance of neighboring beads, a bending potential , and a potential for the local torsion of the polymer chain.Then, the potential energy of a linear polymer with beads is given by: (34) is the distance vector between the centers of mass of neighbor atoms and , is the is the angle between the vectors and .This type of potentials has been used in Gibbs ensemble MD simulations of swelling grained approach, each bead can represent a chain sequence ly longer than the persistence length of the polymer.In this case, the angular dependent parts of the polymer potential can be neglected.The remaining bond length potential is typically (FENE) potential [105][106][107].The resulting chains are rather flexible and can be used to represent suitably long polymers on the respective length scale.If the persistence length is comparable with the polymer length, which is often model is frequently used, which includes the ].In the case of very stiff length fluctuations and represent the whole rigid rod.While this approach neglects many internal degrees of freedom, it enables very the method and facilities of rigid rod rosslinkers We have used this method by studying a system of rigid filaments that can form filament networks r such systems are mixtures of short actin biopolymers and myosin crosslinkers.The system can form various structures, such as and filament networks with and without ving cells, the filament network is a dynamic system.Using ATP, the actin polymers can change their length and the myosin crosslinkers can perform directed walks along the actin filaments.In vitro, in the absence of ATP, the system is passivated and the fixed length.Filaments with lengths in the range of the persistence length can be created with the help of suitable proteins.The resulting biomimetic hydrogels are of great interest the structure of the network depends sensitively on the filament density, the crosslinker and the temperature, and can be altered reversibly, by mechanical, chemical We have studied the system with the help of Monte Carlo simulations.The filaments are represented by long hard spherocylinders, while the crosslinkers are modeled by short spherocylinders with specifically binding ends.A snapsho have a diameter D and the cylindrical part has a length of crosslinkers, respectively.The configuration of a a unit vector that is parallel to the cylinder axis.Steric interactions between spherocylinders are considered by a potential Each crosslinker, , has adhesive sites at the both ends, with which it can physically bind to a neighboring filament, .The definition of the adhesion potential is based on the shortest distances between the axis of filament, , and the adhesive sites of the crosslinker with .The centers of the adhesive sites are localized on the away from the crosslinker's center of mass.We use sites lie inside the hemispheres of the crosslinkers.
The adhesion is mimicked by a square where is the potential depth.The condition filaments (which have a monodisperse length distribution).With a potential range crosslinker can only bind to the filament if the surface of the spherocylinders ar and if the crosslinker axis is almost perpendicular to the filament surface.The total energy is given by The filament hydrogel is modeled with the help of a Metropolis algorithm.At each step, position and/or orientation of a filament or a crosslinker is changed randomly by a small amount configuration is accepted with a probability change of the configurational energy.
In the simulations, we have studied the occurrence of filament bundles and the pe of the system.A cluster is defined as a set of filaments that are interconnected by crosslinker bonds.A filament bundle is a cluster in which the filaments are suitably well aligned.We have analyzed the average structure of the filament hydrogel in phase diagrams as a function of the filament volume fraction , the crosslinker-filament ratio which more than of the filaments are part of a bundle, is called the polymers can change their length and the myosin crosslinkers can perform directed walks along the , in the absence of ATP, the system is passivated and the fixed length.Filaments with lengths in the range of the persistence length can be created with the help of suitable proteins.The resulting biomimetic hydrogels are of great interest as they are biocompatible and network depends sensitively on the filament density, the crosslinker and the temperature, and can be altered reversibly, by mechanical, chemical, or temperature stimuli.
We have studied the system with the help of Monte Carlo simulations.The filaments are represented by long hard spherocylinders, while the crosslinkers are modeled by short spherocylinders with specifically binding ends.A snapshot of the system is given in Figure 6a.The spherocylinders and the cylindrical part has a length of and crosslinkers, respectively.The configuration of a filament is determined by the center of mass o the cylinder axis.Steric interactions between spherocylinders that gets infinite if the two spherocylinders overlap.has adhesive sites at the both ends, with which it can physically bind to a .The definition of the adhesion potential is based on the shortest distances and the adhesive sites of the crosslinker, : . The centers of the adhesive sites are localized on the crosslinker axis a distance away from the crosslinker's center of mass.We use and sites lie inside the hemispheres of the crosslinkers.
The adhesion is mimicked by a square-well potential: is the potential depth.The condition denotes that crosslinkers can only adhere to filaments (which have a monodisperse length distribution).With a potential range crosslinker can only bind to the filament if the surface of the spherocylinders ar and if the crosslinker axis is almost perpendicular to the filament surface.The total energy is given by The filament hydrogel is modeled with the help of a Metropolis algorithm.At each step, position ment or a crosslinker is changed randomly by a small amount configuration is accepted with a probability change of the configurational energy.
In the simulations, we have studied the occurrence of filament bundles and the pe of the system.A cluster is defined as a set of filaments that are interconnected by crosslinker bonds.A filament bundle is a cluster in which the filaments are suitably well aligned.We have analyzed the ent hydrogel in phase diagrams as a function of the filament volume filament ratio , and the adhesion strength .The phase diagram region, in of the filaments are part of a bundle, is called the bundle region 58 polymers can change their length and the myosin crosslinkers can perform directed walks along the , in the absence of ATP, the system is passivated and the actin polymers have a fixed length.Filaments with lengths in the range of the persistence length can be created with the help of they are biocompatible and network depends sensitively on the filament density, the crosslinker-filament ratio, or temperature stimuli.We have studied the system with the help of Monte Carlo simulations.The filaments are represented by long hard spherocylinders, while the crosslinkers are modeled by short spherocylinders 6a.The spherocylinders for filaments and is determined by the center of mass and o the cylinder axis.Steric interactions between spherocylinders and that gets infinite if the two spherocylinders overlap.has adhesive sites at the both ends, with which it can physically bind to a .The definition of the adhesion potential is based on the shortest distances (35) crosslinker axis a distance , so that the adhesion (36) denotes that crosslinkers can only adhere to filaments (which have a monodisperse length distribution).With a potential range the crosslinker can only bind to the filament if the surface of the spherocylinders are closer than and if the crosslinker axis is almost perpendicular to the filament surface.The total energy is given by: (37) The filament hydrogel is modeled with the help of a Metropolis algorithm.At each step, position ment or a crosslinker is changed randomly by a small amount, and the new , where is the In the simulations, we have studied the occurrence of filament bundles and the percolation behavior of the system.A cluster is defined as a set of filaments that are interconnected by crosslinker bonds.A filament bundle is a cluster in which the filaments are suitably well aligned.We have analyzed the ent hydrogel in phase diagrams as a function of the filament volume The phase diagram region, in bundle region.The percolation threshold is determined by studying the cluster size distribution largest cluster : At the percolation threshold, the fraction of the form with the Fisher e at the percolation threshold.In Fig adhesion strength and the filament volume fraction A percolated network forms if the filament volume fraction and the adhesion strength are large enough.This means, that a percolated network can be created and destroyed, reversibly, by changing the hydrogel volume or the temperature.Furthermore, the system fo filament volume fractions.This can be explained as follows: The system favors crosslinkers that are bound to filaments on both ends.Parallel filaments can be interconnected by a large number of crosslinks forming a ladder-like structure.At low filament concentrations, only aligned groups of filaments make it possible that large numbers of crosslinkers can bind on both ends.At high filament concentrations, crosslinkers can bind on both ends, even if the filament network is di phase diagram shows that percolated network may form with and without pronounced bundling and bundling can form with and without the formation of a percolated network.
The network and bundle formation can also be induced by increasing the cr Another interesting aspect is the dependence of the percolation threshold on the filament length [ Interestingly, we found that for a given crosslinker fraction, the binding strength at which percolation sets in, is rather independent of the filament length that has been varied in the range between reproduced qualitatively with a simple analytic model.Note given at a fixed ratio of filament to crosslink number density.In practice, the ratio to filament monomers can be more important.For suitably large filament lengths one has so that the percolation threshold is independent of the filament length if threshold is determined by studying the cluster size distribution and the fraction of filaments in the : At the percolation threshold, the fraction of clusters containing with the Fisher exponent .For a sufficiently large system, one has In Figure 6a, phase diagram of the system is shown as a function of the and the filament volume fraction for fixed crosslinker colated network forms if the filament volume fraction and the adhesion strength are large enough.This means, that a percolated network can be created and destroyed, reversibly, by changing the hydrogel volume or the temperature.Furthermore, the system forms bundles for suitably low filament volume fractions.This can be explained as follows: The system favors crosslinkers that are bound to filaments on both ends.Parallel filaments can be interconnected by a large number of e structure.At low filament concentrations, only aligned groups of filaments make it possible that large numbers of crosslinkers can bind on both ends.At high filament concentrations, crosslinkers can bind on both ends, even if the filament network is di phase diagram shows that percolated network may form with and without pronounced bundling and bundling can form with and without the formation of a percolated network.
The network and bundle formation can also be induced by increasing the cr Another interesting aspect is the dependence of the percolation threshold on the filament length [ Interestingly, we found that for a given crosslinker-filament ratio and a fixed filament volume at which percolation sets in, is rather independent of the filament length that has been varied in the range between and reproduced qualitatively with a simple analytic model.Note, that the independence of the length is of filament to crosslink number density.In practice, the ratio to filament monomers can be more important.For suitably large filament lengths one has tion threshold is independent of the filament length if is kept fixed.
Model network with stiff filaments (blue) and short crosslinkers (yellow) with adhesive ends (red), studied with Monte Carlo simulations.and the fraction of filaments in the of clusters containing filaments is For a sufficiently large system, one has phase diagram of the system is shown as a function of the for fixed crosslinker-filament ratio .colated network forms if the filament volume fraction and the adhesion strength are large enough.This means, that a percolated network can be created and destroyed, reversibly, by changing rms bundles for suitably low filament volume fractions.This can be explained as follows: The system favors crosslinkers that are bound to filaments on both ends.Parallel filaments can be interconnected by a large number of e structure.At low filament concentrations, only aligned groups of filaments make it possible that large numbers of crosslinkers can bind on both ends.At high filament concentrations, crosslinkers can bind on both ends, even if the filament network is disordered.The phase diagram shows that percolated network may form with and without pronounced bundling and The network and bundle formation can also be induced by increasing the crosslink-filament ratio.Another interesting aspect is the dependence of the percolation threshold on the filament length [25].and a fixed filament volume at which percolation sets in, is rather independent of the filament .The results could be that the independence of the length is of filament to crosslink number density.In practice, the ratio of crosslinkers to filament monomers can be more important.For suitably large filament lengths one has is kept fixed.
Model network with stiff filaments (blue) and short crosslinkers (yellow) with ) Snapshot of a network ) Phase diagram as a function of the adhesion strength over temperature and R. Chelakkot, R. Lipowsky and T. Gruhn, Soft permission of The Royal Society of Chemistry.

Summary and Conclusions
Various macroscopic models and molecular simulation techniques have been developed, which allow the study of various kinds of stimuli-sensitive polymer networks.With these methods, the response of polymer networks on external stimuli can be analyzed and optimized on various time and length scales.The behavior of a hydrogel can depend, in a complex way, on many physical and chemical properties of the environment, such as the pH value, ion concentrations, the concentration of specific molecules, external electric fields, temperature, or illumination.The development and optimization of tailored, stimuli-sensitive polymer networks requires a good understanding of the interplay of all these aspects.Here, numerical modeling and simulations are crucial.Molecular simulations require a comparably low amount of parameters.For atomistic simulations only reasonable interatomic potentials are needed, which, in many cases, are provided by simulation software.However, time scales and system sizes that can be studied with atomistic simulations are small so that only a small number of network meshes can be studied.More coarse-grained simulations require model interaction parameters, which are often not available for the polymer network of interest.The parameters can be extracted from atomistic simulations on a small length scale or by comparing material properties of the model system with experimental measurements.If the parameters are found, very detailed information can be found about system properties like the swelling behavior of the network, the breaking and healing of bonds, changes of the polymer configurations and the network structure or, for networks with suitably small mesh sizes, the elastic and inelastic responses of the network on applied strains.Results of molecular simulations help to interpret the experimental measurements but they can also support macroscopic models, which depend on material properties such as the elastic constants of the polymer network for a given composition of the solvent or the permeability of a hydrogel for a specific molecule.For systems with chain lengths that are too large for molecular simulations, approaches like the self-consistent field theory can be a good alternative, especially, if the system consists of copolymers and forms nanostructures.Our extended self-consistent field theory allows the modeling of polymer systems on an intermediate length scale and helps to fill the gap between macroscopic theories and molecular simulations.Thus far, only a very small number of multi-scale methods have been used for the modeling of stimuli-sensitive networks on various length scales.This is one of the major tasks for the future.
One particularly important aspect of smart polymer networks is the response time.Quantitative predictions of absolute time scales are a big challenge for numerical methods.Many macroscopic models are restricted to stationary or steady-state systems.Macroscopic numerical studies of network dynamics are based on parameters like diffusion constants that have to be taken from molecular simulations or experiments.This shows, once more, the necessity of numerical studies on various levels of detail and of multiscale studies that combine these methods.
Another challenge is the numerical study of polymer networks with spatial inhomogeneities on the microscale rather than on the nanometer scale.Experimentally, such inhomogeneities can be inferred, for example, by differential photo-crosslinking of polymers [110] or by applying an irradiation pattern on a thermosensitive gel [111].
phase, ].A multiphase approach, based on Macroscopic hydrogel theories are frequently based on conservation laws and constitutive ].If the ion, solvent, and polymer phases are incompressible, mass conservation (1) are the volume fractions of the polymer network and the solvent phase, while ion phases.The volume fractions are proportional to .In the absence of chemical reactions, the concentrations obey the (2) .In this article, we use the summation 41,43-45].If negatively , electric neutrality requires:

( 7 )
tial geometry of the hydrogel.Strong, ].If the deformations are (8) are components of the full pressure tensor in the hydrogel region and the unit are the first and second Lamé coefficients of the

Figure 3 .
Figure 3. Evolution of the gel fraction with time.In the shown time period, the boundary ion concentration is linearly decreased from to , and decreased again to 289, 2011, 513-521, Phase field stimulation, Li, D., Yang, H.L., Emmerich, H., kind permission from Springer Science+Business Media B.V.

Figure 3
Figure 3 demonstrates the swelling, shrinking and swelling dynamics of the hydrogel as the boundary ion concentration is decreased to step 163,500), and decreased again to the initial value.The gel fraction follows the boundary ion concentration with a short delay, which is too small to be visible in Fig investigated in the numerical study [ shrinking with increasing ion concentr experiments, in which the pH is kept constant and the ion concentration is varied [Our new phase field model is an alternative approach for studying the dynamics of hydrogel swelling induced by the ion concentration of the bath.With the chosen decrease or increase rate for the boundary ion concentration, the hydrogel swells and shrinks without significant hysteresis.The method goes without explicit boundary conditions at the hydrogel and can be easily adapted to other hydrogel systems.For example, one can study a hydrogel shell that includes an active agent.Then, a concentration field of drug molecules can be added to the given model in order to study the drug release kineti Evolution of the gel fraction with time.In the shown time period, the boundary ion concentration is linearly decreased from to , and decreased again to [Springer, Colloid Polymer Science, 521, Phase field model simulations of hydrogel dynamics under stimulation, Li, D., Yang, H.L., Emmerich, H., Figure 7, © Springer-Verlag 2011 kind permission from Springer Science+Business Media B.V. 3 demonstrates the swelling, shrinking and swelling dynamics of the hydrogel as the boundary (time step 54,000), increased to and decreased again to the initial value.The gel fraction follows the boundary ion concentration with a short delay, which is too small to be visible in Fig investigated in the numerical study

52
Evolution of the gel fraction with time.In the shown time period, the boundary , then increased [Springer, Colloid Polymer Science, model simulations of hydrogel dynamics under chemicalVerlag 2011] with 3 demonstrates the swelling, shrinking and swelling dynamics of the hydrogel as the boundary (time step 54,000), increased to (time and decreased again to the initial value.The gel fraction follows the changes of the boundary ion concentration with a short delay, which is too small to be visible in Figure2but has been

Figure 4 .
Figure 4. Nanostructures formed by a polymer blend of A and B homopolymers and reversibly crosslinked AB block copolymers.Shown are sketches of ( and (b) the lamellar phase.Blue and yellow lines represent copolymers, yellow circles denote crosslinkers.A and B homopolymer densities are indicated by the light blue and yellow background colors.Results of self theory calculations for the monomer density structure [22].Figure (a) and ( Emmerich: Phase behavior of polymer bl field theory study.Journal of Material Research; preparation, doi:10.1557/jmr.2013.315 Nanostructures formed by a polymer blend of A and B homopolymers and slinked AB block copolymers.Shown are sketches of (a) the hexagonal phase ) the lamellar phase.Blue and yellow lines represent the A and B parts of the yellow circles denote crosslinkers.A and B homopolymer densities are he light blue and yellow background colors.Results of self theory calculations for the monomer density in (c) a hexagonal, and ( ) and (b) reproduced with permission from Thomas Gruhn Phase behavior of polymer blends with reversible crosslinks-Journal of Material Research; Materials Research Society, 10.1557/jmr.2013.315.(b) 55 and

Figure 5 .
Figure 5. Phase diagram of the phase boundaries between the hexagonal phase on the inside (small phase on the outside.Phase boundaries are shown for (red), (green), Reproduced with permission from Thomas Gruhn, Heike Emmerich polymer blends with reversible crosslinks Material Research; Materials Research Society,

57 simulations,
in which the dynamics of interpenetrating spheres resemble the flow behavior of different ]).DPD simulations have been used to investigate the release of ].Polyelectrolytes in electric fields have been theories, molecular simulations provide more to much smaller time and length scales.In multiscale simulations, material based or other coarse-grained methods Atomistic simulations of polymers are very time consuming.Consequently, only a few atomistic 101]

Figure 6 .
Figure 6.Model network with stiff filaments (blue) and short crosslinkers (yellow) with adhesive ends (red), studied with Monte Carlo simulations.( structure.(b) Phase diagram as a function of the adhesion strength over temperature the filament volume fraction Matter, 2009, 5, 1504)-Reproduced (a) Snapshot of a network ) Phase diagram as a function of the adhesion strength over temperature the filament volume fraction .Figures from (R. Chelakkot, R. Lipowsky and T. Gruhn, Reproduced with permission of The Royal Society of (b)