Coil-stretch-like transition of elastic sheets in extensional flows
CCoil-stretch-like transition of elastic sheets in extensionalflows
Yijiang Yu a and Michael D. Graham a ∗ The conformation of a long linear polymer dissolved in fluid and exposed to an extensional flow is well-known to exhibit a “coil-stretch" transition, which for sufficiently long chains can lead to bistability.The present work reports computations indicating that an analogous “compact-stretched" transitionarises in the dynamics of a thin elastic sheet. Sheets of nominally circular, square or rectangular shapeare simulated in planar and biaxial flows using a finite element method for the sheet conformationsand a regularized Stokeslet method for the fluid flow. If a neo-Hookean constitutive model is usedfor the sheet elasticity, the sheets will stretch without bound once a critical extension rate, ascharacterized nondimensionally by a capillary number, is exceeded. Nonlinear elasticity, representedwith the Yeoh model, arrests the stretching, leading to a highly-stretched steady state once thecritical capillary number is exceeded. For all shapes and in both planar and biaxial extension, aparameter regime exists in which both weakly stretched (compact) and strongly stretched states canbe found, depending on initial conditions. I.e. this parameter regime displays bistability. As in thelong-chain polymer case, the bistable behavior arises from the hydrodynamic interaction betweendistant elements of the sheet, and vanishes if these interactions are artificially screened by use of aBrinkman model for the fluid motion. While the sheets can transiently display wrinkled shapes, allfinal shapes in planar and biaxial extension are planar.
In recent years, thin-structured materials have gained consider-able attention in many fields. For example, freely suspendedsheets immersed in a fluid environment arise the synthesis andprocessing of polymer networks . Polymer sheets are often usedas planar substrates in soft origami-inspired devices, which fur-ther assemble into complicated structures . Thin polymer filmsare also used as the basis of self-propelled swimming devices,with living cells deposited on the sheet to drive the motion .Highly extensible polymer films, hydrogels structures and theircomposites have shown great potential in biomedical and elec-tronic applications . Many of above examples involve com-plicated fluid-structure interactions, but little is known about thedynamics of soft sheet-like particles in fluid. The present workaims to shed light on some aspects of these dynamics throughsimulations of soft elastic sheets in two canonical flow fields thatarise in many applications such as coating operations: planar andbiaxial extensional flows.While soft elastic sheets have not been studied in these flows,long-chain polymer molecules have been. An important phe- a Department of Chemical & Biological Engineering, University of Wisconsin-Madison,1415 Engineering Drive, Madison, WI 53706, USA. ∗ [email protected] nomenon in this situation is the coil-stretch transition : poly-mers show abrupt change in conformation from a coil state to astrongly stretched state when reaching a critical flow strength inplanar extension. The discontinuity in the stretched length marksa bistable regime that corresponds to the existence of multiple sta-ble (or more precisely, metastable) conformations. This behaviorwas first predicted by De Gennes in 1974 , based on the ideathat, compared to a fully stretched state, in the coiled state theinner monomers are hydrodynamically screened and thereforeless exposed to the strong stretching flow. Fuller and Leal mod-eled this effect using a bead-spring dumbbell with conformation-dependent drag on the beads, showing that this simple approachcould yield coil-stretch hysteresis . Schroeder et al. experimen-tally confirmed this prediction, visualizing the hysteresis behaviorwith fluorescent dyed Escherichia coli DNA molecules in dilute so-lution under planar extension . They also probed the hysteresisphenomenon with a Brownian dynamics simulation of a bead-spring polymer model and found a double-welled effective poten-tial of polymer conformational energy around the transition flowstrength that leads to hysteresis in conformation. The coil-stretchtransition and hysteresis are also observed in simulations beyondthe dilute limit. Mohammad et al. applied molecular dynam-ics simulation to study atomistic entangled polyethylene melts inplanar extension . They found that the bistable region in thepolymer melts is wider than in dilute solution, with a higher tran- Journal Name, [year], [vol.] ,,
In recent years, thin-structured materials have gained consider-able attention in many fields. For example, freely suspendedsheets immersed in a fluid environment arise the synthesis andprocessing of polymer networks . Polymer sheets are often usedas planar substrates in soft origami-inspired devices, which fur-ther assemble into complicated structures . Thin polymer filmsare also used as the basis of self-propelled swimming devices,with living cells deposited on the sheet to drive the motion .Highly extensible polymer films, hydrogels structures and theircomposites have shown great potential in biomedical and elec-tronic applications . Many of above examples involve com-plicated fluid-structure interactions, but little is known about thedynamics of soft sheet-like particles in fluid. The present workaims to shed light on some aspects of these dynamics throughsimulations of soft elastic sheets in two canonical flow fields thatarise in many applications such as coating operations: planar andbiaxial extensional flows.While soft elastic sheets have not been studied in these flows,long-chain polymer molecules have been. An important phe- a Department of Chemical & Biological Engineering, University of Wisconsin-Madison,1415 Engineering Drive, Madison, WI 53706, USA. ∗ [email protected] nomenon in this situation is the coil-stretch transition : poly-mers show abrupt change in conformation from a coil state to astrongly stretched state when reaching a critical flow strength inplanar extension. The discontinuity in the stretched length marksa bistable regime that corresponds to the existence of multiple sta-ble (or more precisely, metastable) conformations. This behaviorwas first predicted by De Gennes in 1974 , based on the ideathat, compared to a fully stretched state, in the coiled state theinner monomers are hydrodynamically screened and thereforeless exposed to the strong stretching flow. Fuller and Leal mod-eled this effect using a bead-spring dumbbell with conformation-dependent drag on the beads, showing that this simple approachcould yield coil-stretch hysteresis . Schroeder et al. experimen-tally confirmed this prediction, visualizing the hysteresis behaviorwith fluorescent dyed Escherichia coli DNA molecules in dilute so-lution under planar extension . They also probed the hysteresisphenomenon with a Brownian dynamics simulation of a bead-spring polymer model and found a double-welled effective poten-tial of polymer conformational energy around the transition flowstrength that leads to hysteresis in conformation. The coil-stretchtransition and hysteresis are also observed in simulations beyondthe dilute limit. Mohammad et al. applied molecular dynam-ics simulation to study atomistic entangled polyethylene melts inplanar extension . They found that the bistable region in thepolymer melts is wider than in dilute solution, with a higher tran- Journal Name, [year], [vol.] ,, a r X i v : . [ c ond - m a t . s o f t ] N ov ition flow strength.Similar conformational transitions have been observed in sys-tems other than polymers as well. Kantsler et al. performed ex-periments to study the stretching of a single tubular vesicle sus-pended in planar extensional flow . They found a transition inconformation from a tubular to a dumbbell shape when reachinga critical extension rate. With further increasing flow strength,the shape evolves to a pearling state due to unstable higher or-der shape modes. This behavior is also found in a computa-tional study by Narsimhan et al. , who also considered otherextensional flows. Kumar et al. recently developed a microfluidicStokes trap to make a vesicle remain at the stagnation point andextended the work of Kantsler et al. by giving a more detailedanalysis of the vesicle shape transition in the parameter space ofreduced volume and bending modulus .In contrast to the substantial literature on polymers or vesi-cles, there is only a small amount of work on the dynamics ofdeformable freely suspended sheets in flow and no prior studyhas probed the possibility of hysteretic conformational transitionsof extensible sheets in extensional flow . (Though there existmany studies on flapping sheets with a clamped edge in flow(flapping flags) , this is not our interest here.) Motivated byapplications of relatively inextensible sheetlike materials such asgraphene and boron nitride, Green and Xu investigated the dy-namics of freely suspended nanosheets in both extensional andshear flow . They modeled a nominally square sheet witha bead-rod network inspired by coarse-grained models of poly-mers in solution, simulating the model with a Brownian dynam-ics method. Because of the inextensibility of the rods, the modelthey consider is highly resistant to stretching deformations, somuch of the dynamics in this system are determined by the rela-tive importance of flow, bending resistance and Brownian motion.In biaxial extensional flow, a sheet with large bending stiffnesstends to stay at a flat conformation. With small bending stiff-ness, Brownian motion leads to crumpling . In shear flow atlow bending stiffness, they found that sheets undergo a cycliccrumple-stretch-crumple motion . Botto and co-workers per-formed molecular dynamic simulations to study the liquid phaseexfoliation of graphite sheets in shear. They studied the criticalshear rate for exfoliation to occur in different solvent and devel-oped a theoretical model that can predict the critical shear ratein simulation . They also applied molecular dynamics simula-tions to explore the effects of hydrodynamic slip on the dynam-ics of a sheet in shear . Interestingly, in the presence of slip,they found that a rigid nanosheet can take on a stable steady ori-entation with a finite orientation angle with respect to the flowdirection. Finally, Dutta and Graham investigated the dynamicsof a piecewise rigid creased sheet – a model of an origami figureknown as a Miura pattern . Depending on geometry and initialconditions, this sheet can show steady, periodic or quasiperiodicorientational and conformational dynamics in shear.The present work describes the behavior of sheetlike particlesin a parameter regime that has not been explored in the studiesabove. We consider soft (highly extensible) elastic sheets, takingBrownian forces to be negligible compared to hydrodynamic andelastic ones, and focus on planar and biaxial extensional flows. We consider three different rest shapes: disc, square, and rect-angle. The rest of the paper is organized as follows: Section 2summarizes the model we use and the discretization and solutionmethods. Section 3 shows detailed results for the sheet dynamics,illustrating a bistability phenomenon analogous to the coil-stretchtransition. Concluding remarks are presented in Section 4. We consider a very thin, neutrally buoyant, elastic sheet sus-pended in an unbounded , incompressible Newtonian fluid withviscosity η and density ρ , and subjected to planar and biaxialextensional flow fields. Expressed in Cartesian coordinates, thefluid velocity v ∞ in the absence of the particle is given in pla-nar extension by v ∞ = ˙ ε [ x , , − z ] T and in biaxial extension by v ∞ = ˙ ε [ x / , y / , − z ] T , with ˙ ε the strain rate.The sheet is modeled as an elastic continuum. Fluid cannotpass through the sheet, and each point on the sheet is taken tomove with the local velocity of the fluid. That is, the usual no-penetration and no-slip boundary conditions are applied for thefluid velocity at the sheet surfaces. In this case, no specificationof the details of fluid-solid interfacial forces is necessary. As notedabove, for atomically smooth surfaces, the no-slip boundary con-dition may not apply, in which case interesting dynamics havebeen seen to arise in simulations of rigid sheets in shear flow .That case is not considered here. In this work, we focus onsheets with three rest shapes: (a) disc with radius a , (b) squarewith edge a and (c) rectangle with short edge a and long edge a . The rest thickness h is much less than a . The sheets are takento be sufficiently thin that the fluid exerts no traction along theiredges.The mechanical response of a very thin sheet is split into twoparts: in-plane shearing elasticity, with strain energy E s , and out-of-plane bending elasticity, with strain energy E b . The total strainenergy is thus E = E s + E b . (1)The in-plane shearing energy E s on the sheet surface Γ can bewritten as E s = (cid:90) Γ W dS , (2)where W is the areal shear strain energy density that is definedpointwise on the sheet surface . We consider the sheet as isotropicand incompressible, and describe the in-plane energy with theYeoh model, a model used to simulate rubber-like materials .Its strain elastic energy W is given by: W = G (cid:16) λ + λ + λ − (cid:17) + cG (cid:16) λ + λ + λ − (cid:17) . (3)Here, G is a two-dimensional shear elasticity modulus given by G h , where G is the usual shear modulus. The quantities λ and λ are the local principal stretch ratios along the (local) tan-gential directions of the sheet surface, while the normal stretchratio λ is determined by incompressibility of the sheet material: λ λ λ = . The parameter c determines the weight of the cu-bic terms in the energy formulation and determines whether the Journal Name, [year], [vol.] , aterial is strain-softening or strain-hardening. When c = , theYeoh model reduces to the strain-softening neo-Hookean (NH)model.For bending energy, we apply a simple linear bending modelthat penalizes deflection between neighboring discretized ele-ments. Details of the discretization will be covered in laterparagraphs. This energy has the form E b = ∑ ad jacent α , β k b (cid:2) − cos ( θ αβ − θ ) (cid:3) , (4)where k b is a bending constant, θ αβ is the angle between twoneighboring elements and θ is the angle at equilibrium. In thisstudy, we assume the sheet is flat at rest, so θ = for all adjacentelement pairs. The bending constant k b = √ ( K B + K B ) denotesan averaged bending modulus, which is derived from the Hel-frich bending energy of a spherical shell with a zero spontaneouscurvature . Here K B is the bending modulus and is related tothe Gaussian curvature modulus K B by K B / K B = ν − , where thePoisson ratio ν = / for an incompressible sheet . A nondimen-sional bending stiffness can be defined: ˆ K B = K B / a G ,The simulations performed and the results presented are innondimensional form. Lengths are scaled with a and time withthe inverse extension rate / ˙ ε . Forces are nondimensionalized byviscous drag force η ˙ ε a . With these choices, two nondimensionalparameters arise that characterize the deformability of the sheetby the flow: the in-plane deformability is determined by the capil-lary number Ca = η ˙ ε a / G , with the competition between flow andbending deformations characterized by Ca / ˆ K B . It is worth point-ing out here that in all cases studied here, the sheet becomes flatat long times, in which case the bending energy becomes irrele-vant. It does, however, have an effect on dynamics at short times,determining the transient wrinkling behavior that the sheet ex-hibits while orienting with the flow.Our results below indicate that complex behavior in extensionalflow occurs once Ca (cid:38) . . We present here some estimates to il-lustrate how these results will be relevant to microscale sheets ofspecific materials. Materials like polymer hydrogels usually haveshear modulus G D ≈ ) . If a disc-shaped polyethyleneglycol (PEG) hydrogel sheet with radius of 100 µ m and thicknessof 100 nm is suspended in water ( η ∼ − Pa · s ), a modest ex-tension rate of about s − will yield Ca = . .Other interesting two-dimensional materials, such as grapheneflakes, have much too large a shear modulus ( G D ≈
300 GPa ) toyield a capillary number in the range where they will be stronglystretched by flow, except in extreme conditions.In addition, the thermal fluctuation for polymer sheets at roomtemperature is considered negligible compared to strain energy.If we still consider x‘the above PEG hydrogel, its estimated shearstrain energy ( ∼ k B T ) and bending energy ( ∼ × k B T ) aremuch larger than thermal energy. Therefore, thermal fluctuationswill negligible influence on the dynamics. The numerical method for the elasticity problem is adapted fromCharrier et al. . We simulate the sheet by keeping track of material points as nodes on the sheet surface. The sheet sur-face is discretized into triangular elements with a node at eachcorner. From the above two energies, we obtain the elastic force F e exerted on each node from the first variation of the total en-ergy with respect to the nodal displacements. Each discretizedelement of the sheet is assumed to have homogenous deforma-tion, so the element edges always remain linear. The deformedelement is compared to its equilibrium shape under a local coor-dinate transformation via a rigid body rotation and the displace-ment for any point inside the element is obtained by linear inter-polation from nodes. The detailed implementation can be foundin references . The total nodal force F e is evaluated bysumming the elastic force F e , i due to deformation of each sur-rounding element shared by the node: F e = ∑ i F e , i , where thesum is over all elements meeting at the node. For the resultsshown, we discretize the disc with 1600 elements and 841 nodes,the square with 2048 elements and 1089 nodes, and the rectanglewith 1024 elements and 561 nodes. We have verified that changesin mesh resolution lead to only small quantitative changes in theresults and no qualitative changes. We consider the case of very small sheets, such that the particleReynolds number Re = ρ ˙ ε a / η is much less than unity and thefluid motion is governed by the Stokes equation. Neglecting theinertia of the sheet, the total elastic force F e exerted on each nodeby the surrounding solid is balanced by the force F h exerted by thefluid on the node: F e + F h = . (5)Thus, the sheet exerts a force on the fluid at each nodal posi-tion. To account for the fact that the forces are not completelylocalized to the nodal positions, we use the method of regular-ized Stokeslets : the force F i exerted by node i on the fluid cor-responds to a regularized force density f κ = F δ κ ( x − x i ) , where δ κ ( x ) is a regularized delta function with regularization parame-ter κ . Thus, the governing equations are the Stokes equation withregularized nodal forces and the continuity equation: − ∇ p + η ∇ v + ∑ i F i δ κ ( x ) = ∇ · v = . (6)Here the sum is over the nodal positions. The velocity fieldgenerated due to a regularized point force f = F δ κ ( x ) can be rep-resented using a regularized Stokeslet G κ : v κ ( x ) = G κ ( x ) · F . (7)As / κ → , G κ reduces to the usual Stokeslet operator G ( x ) = / πη r ( I + xx / r ) , where r = | x | . There are many ways to regu-larize a delta function δ κ ( x ) ; we choose a regularization functionfor which the deviation between G and G κ decays exponentiallyas κ r → ∞ : δ κ ( x ) = κ √ π exp ( − κ r ) (cid:20) − κ r (cid:21) . (8) Journal Name, [year], [vol.] ,,
300 GPa ) toyield a capillary number in the range where they will be stronglystretched by flow, except in extreme conditions.In addition, the thermal fluctuation for polymer sheets at roomtemperature is considered negligible compared to strain energy.If we still consider x‘the above PEG hydrogel, its estimated shearstrain energy ( ∼ k B T ) and bending energy ( ∼ × k B T ) aremuch larger than thermal energy. Therefore, thermal fluctuationswill negligible influence on the dynamics. The numerical method for the elasticity problem is adapted fromCharrier et al. . We simulate the sheet by keeping track of material points as nodes on the sheet surface. The sheet sur-face is discretized into triangular elements with a node at eachcorner. From the above two energies, we obtain the elastic force F e exerted on each node from the first variation of the total en-ergy with respect to the nodal displacements. Each discretizedelement of the sheet is assumed to have homogenous deforma-tion, so the element edges always remain linear. The deformedelement is compared to its equilibrium shape under a local coor-dinate transformation via a rigid body rotation and the displace-ment for any point inside the element is obtained by linear inter-polation from nodes. The detailed implementation can be foundin references . The total nodal force F e is evaluated bysumming the elastic force F e , i due to deformation of each sur-rounding element shared by the node: F e = ∑ i F e , i , where thesum is over all elements meeting at the node. For the resultsshown, we discretize the disc with 1600 elements and 841 nodes,the square with 2048 elements and 1089 nodes, and the rectanglewith 1024 elements and 561 nodes. We have verified that changesin mesh resolution lead to only small quantitative changes in theresults and no qualitative changes. We consider the case of very small sheets, such that the particleReynolds number Re = ρ ˙ ε a / η is much less than unity and thefluid motion is governed by the Stokes equation. Neglecting theinertia of the sheet, the total elastic force F e exerted on each nodeby the surrounding solid is balanced by the force F h exerted by thefluid on the node: F e + F h = . (5)Thus, the sheet exerts a force on the fluid at each nodal posi-tion. To account for the fact that the forces are not completelylocalized to the nodal positions, we use the method of regular-ized Stokeslets : the force F i exerted by node i on the fluid cor-responds to a regularized force density f κ = F δ κ ( x − x i ) , where δ κ ( x ) is a regularized delta function with regularization parame-ter κ . Thus, the governing equations are the Stokes equation withregularized nodal forces and the continuity equation: − ∇ p + η ∇ v + ∑ i F i δ κ ( x ) = ∇ · v = . (6)Here the sum is over the nodal positions. The velocity fieldgenerated due to a regularized point force f = F δ κ ( x ) can be rep-resented using a regularized Stokeslet G κ : v κ ( x ) = G κ ( x ) · F . (7)As / κ → , G κ reduces to the usual Stokeslet operator G ( x ) = / πη r ( I + xx / r ) , where r = | x | . There are many ways to regu-larize a delta function δ κ ( x ) ; we choose a regularization functionfor which the deviation between G and G κ decays exponentiallyas κ r → ∞ : δ κ ( x ) = κ √ π exp ( − κ r ) (cid:20) − κ r (cid:21) . (8) Journal Name, [year], [vol.] ,, ith this choice, G κ ( x ) = erf ( κ r ) πη r (cid:16) I + xx r (cid:17) + κ e − κ r π / η (cid:16) I − xx r (cid:17) . (9)In the simulations, κ must be chosen to scale with the mini-mum node-to-node distance l min . We take κ l min = . , whichis obtained using a validation case discussed below. Due to theno-slip and the no-penetration boundary conditions, the velocityof each node on the sheet surface equals the fluid velocity at thenodal position. This velocity is the sum of the free velocity v ∞ in the absence of the sheet and the perturbation velocity v p generated by the nodal forces: v ( x ) = v ∞ ( x ) + v p ( x ) = v ∞ ( x ) + ∑ i G κ ( x − x i ) · F i . (10)The position of each node is determined by integrating the ve-locity with a fourth-order Runge-Kutta method. The time stepapplied in simulation follows ∆ t = . l min , with an upper limitof × − strain unit. Here we present a validation case for the regularized Stokesletformulation: a rigid disc with radius m and no-slip boundary con-ditions lying in the stretching x − y plane of unbounded biaxialextensional flow. We set Ca = . for a disc (discretized by 1600triangular elements with l min = . ) to ensure negligible defor-mation on the sheet surface. For this case, an analytical solutioncan be found by solving the stream function form of the Stokesequation . In cylindrical coordinates with r and z the usualradial and axial directions, the solution ψ p for the perturbationaway from pure biaxial extension is ψ p ( r , z ) = π r z ˙ ε (cid:18) cot − ( λ ) − λλ + (cid:19) , (11)where λ = (cid:115) r + z − m + (cid:112) ( r + z − m ) + m z m . (12)Figure 1 compares the perturbation velocity norm along an ar-bitrarily chosen path away from disc center, showing excellentagreement between the analytical and regularized Stokeslet solu-tions. The perturbation velocity decays as r − because to leadingorder in / r , a neutrally-buoyant particle acts like a force dipole. We begin the discussion of sheets in planar extensional flow withan important general observation. Consider an initial conditionwhere the sheet is flat, with arbitrary orientation. We find thatthis initial condition will almost always evolve at long time to aflat shape aligned with the x − y plane, where the x - and y -axesare the extensional and neutral directions for the flow, respec-tively. The only exception is the case where the sheet is initiallyoriented perfectly in the x − z plane. By symmetry, this sheet’s Fig. 1
Perturbation velocity norm | v p | vs. distance r along a line awayfrom the origin for a disc in biaxial flow. The line chosen follows thevector [ , , ] . The black dashed line is proportional to r − . orientation will remain in this plane for all time. However, thissituation seems to be unstable for a sheet of any shape with anyfinite deformability, as exemplified in Figure 2 for a disc. In allcases studied, the transient dynamics, and in particular the de-gree of transient wrinkling as the sheet is compressed along the z -direction, depend on Ca and ˆ K B , but eventually a flat state alignedwith the x − y plane is reached. Accordingly, all subsequent resultsconsider sheets with this orientation and with ˆ K B = × − . Onlythe in-plane deformability ( Ca ) is varied in the rest of the paper,as ˆ K B has no influence on the final conformation.Another stability issue arises for sheets whose initial shapes arenot circular. Figure 3 shows steady states for the three shapeswe consider here, for low Ca and the neo-Hookean constitutivemodel (Eq.3 with c = ). For a disc, the steady state is roughlyellipsoidal, as shown in Figure 3a. For a non-radially symmetricgeometry like a square, we also need to consider its stable orien-tation. Even at very small Ca , an initial condition where the sidesof the square are aligned with the x and y coordinate axes willeventually rotate until the diagonal aligns with the flow direc-tion, reaching a steady state as indicated in Figure 3b. (At Ca = ,where the sheet is perfectly rigid, Stokes flow reversibility wouldprevent this reorientation.) By contrast, the rectangular shape westudy here remains symmetric at low Ca , as illustrated in Figure3c).Continuing with the neo-Hookean model, Figure 4 shows thesteady-state extension l s as a function of Ca for the three shapes.In each case, there is a critical value of Ca beyond which the sheetwill stretch without bound, not reaching a steady state shape.The appearance of the singularity is caused by the strain-softeningproperty of the neo-Hookean model – the increasing flow strengthovercomes the elastic response of the sheet – and is analogous tothe behavior found in a bead-spring dumbbell model of a poly-mer molecule when the spring force obeys Hooke’s law and theWeissenberg number exceeds a critical value. In the case of thedumbbell model, if the Hookean spring model is replaced by afinitely extensible spring, the singularity vanishes and the chain Journal Name, [year], [vol.] , eaches a steady state where the chain is highly stretched. Thisis the simplest version of the “coil-stretch” transition of flexiblepolymers. In the sheet case, we will refer to conformations at Ca values below the transition as “compact".To further explore the compact-stretched transition, we turnto the strain-hardening Yeoh model. At low Ca , a finite, smallvalue of c leads to negligible change in conformation. However,once the critical capillary number is exceeded, the the behavioris very different: the singularity is replaced by a sudden jump inthe stretched length within a small Ca interval around the critical Ca . Figure 5 shows steady states for our three rest shapes with c = × − at values of Ca just above the singularity in the neo-Hookean case, i.e., Ca = . , . and . for disc, square, andrectangle, respectively. All cases take on strongly stretched con-formations, and in addition, the rectangle takes on an asymmetricshape: reflection symmetry across the extensional axis is broken.The following paragraphs describe the parameter dependence infurther detail.Figure 6a shows steady state stretching length of a disc vs. Ca for c = × − . At sufficiently low Ca the only steady state iscompact (blue). Once Ca exceeds a critical value Ca c , the com-pact steady state branch loses existence and the length evolvesto a stretched state (red), as indicated by the upward black ar-row. This type of transition is known as a saddle-node bifurca-tion. The critical value Ca c splits the result into two branches inthe parameter space: a compact branch and the stretched branch.Meanwhile, the discontinuity in the stretched length suggests theexistence of a bistable region that corresponds to hysteresis as Ca is increased or decreased quasistatically. Indeed this is thecase – if we use an initial condition corresponding to a stretchedstate, then as long as Ca is greater than some lower limit, denoted Ca s , then the final state will be highly stretched. For Ca < Ca s ,the stretched steady state branch loses existence, and an initiallystretched sheet will relax to a compact steady state, as indicatedby the downward black arrow. In general, this type of bistable be-havior implies the existence of a branch of unstable steady statesintermediate between the upper and lower ones. This branch isindicated schematically by the black dashed lines – we have notexplicitly computed it. Figures 6c-e show the evolution of eitheran initially compact or an initially stretched state for Ca < Ca s , Ca s < Ca < Ca c , and Ca > Ca c , respectively.In polymer dynamics, the origin of coil-stretch hysteresis is un-derstood to be the conformation-dependent hydrodynamic inter-actions. The situation here is somewhat analogous. When the sheet remains at a compact state, the strain is relatively small andthere is not a substantial change in surface area for the fluid toexert stresses on. In the stretched state, however, the surface areais substantially larger, leading to larger forces on the sheet fromthe flow, preventing relaxation back to a compact state. Conse-quently, the sheet maintains a deformed state that is larger in sizethan the compact state under the same Ca . In Section 3.3, weverify this picture by examination of an artificial model in whichhydrodynamic interactions are artificially screened out.This bistability phenomenon only exists for small c . The steady-state stretch for c = × − is shown in Figure 6b; while thereis a sharp increase in l s as Ca increases, there no longer exists adiscontinuous transition or a multiplicity region. Figure 7a showsthe steady states over a range of c . All choices of c share a similarcompact branch, as the cubic energy term does not play a rolewhen the strain is small. The discontinuity between the compactand the stretched branch becomes more evident with decreasing c , where hysteresis arises. In the linear polymer case, bistable be-havior in extension is only observed for very long chain – highlyextensible – polymers, where the equilibrium coil size and thecontour length differ by orders of magnitude. In the present case,the extensibility is represented by the nonlinearity parameter c inthe strain energy. Thinking of the sheet as a crosslinked polymernetwork, small c roughly corresponds to a material with a largenumber of polymer segments between crosslinks: c = is the neo-Hookean limit where the material can stretch without limit. As c increases, large deformations are increasingly penalized, as in thecase of a linear chain with finite extensibility. Consistent with thelinear polymer case, for the sheets, bistability vanishes as exten-sibility decreases ( c increases).Figure 7b summarizes the values of Ca c and Ca s as a functionof c . Hysteresis exists in the gray region between the two curves.Note that Ca c increases only slightly with c , as the nonlinear elas-tic behavior only has a small influence for a compact state. The Ca s branch increases faster with c as a result of the increasingenergy penalty for a stretched state as c increases. Based on thechosen mesh resolution, we note that Ca c and Ca s are not sen-sitive to further increase in resolution; increasing the number ofelements by 10 % shifts Ca s and Ca c in the third decimal place.Since we found compact-stretched transition and hysteresis be-havior for a disc sheet in planar extensional flow, we expect asimilar transition to be found for a different sheet geometries. Wenow briefly present the behavior of a square in planar extensionalflow. As previously mentioned, the square shape preferentially Fig. 2
Time evolution of a deformable disc initially aligned with the x − z plane in planar extensional flow ( Ca = .
2; ˆ K B = × − ). Journal Name, [year], [vol.] , ig. 3 Stable steady-state conformation of a sheet in planarextensional flow with neo-Hookean elasticity at Ca = . : (a) disc. (b) asquare. (c) rectangle. The scale bar 1 dimensionless length unit. Fig. 4
Steady-state stretched length l s vs. Ca for deformable sheetswith neo-Hookean elasticity in planar extensional flow. The data pointsfilled with light color correspond to the shapes shown in Figure 3. Thecritical Ca for each geometry are: disc ( Ca = . ), square ( Ca = . ),rectangle ( Ca = . ). The subplot indicates how length l is measuredand the preferred stable orientation. Fig. 5
The steady-state conformation of a stretched sheet in planarextensional flow with Yeoh elasticity with c = × − : (a) a disc at Ca = . . (b) a square at Ca = . . (c) a rectangle at Ca = . . Thescale bar is 2 dimensionless length units. orients with diagonally opposed corners in the extension direction(Figure 3b and Figure 5b) for both the compact and the stretched state. We provide an example of compact-stretched transition andhysteresis behavior for a square in Figure 8a, for c = × − . Fig-ure 8b shows the region in Ca − c space where hysteresis can befound.In contrast, the rectangle is less symmetric than the disc or thesquare. From the previous steady-state analysis, we have alreadynoted the existence of an asymmetric stretched steady state. Thisresult indicates that there exists a symmetry breaking bifurcationthat interacts with the bifurcation behavior shown above. Figure9a illustrates the length evolution for a rectangle initially in itsrest state and oriented with the extension direction. This sheet isinitially stretched by the flow symmetrically relative to the x − z plane until it reaches a symmetric stretched state. The sheet staysat this unstable symmetric state for some time, but then evolvesto a tilted asymmetric stable conformation. This observation sug-gests that there is a stretched symmetric state that is unstablewith respect to symmetry-breaking perturbations, leading the fi-nal steady state to be asymmetric.Figures 9b and c shows bifurcation diagrams for c = × − and × − , respectively. (We can compute symmetric steadystates that unstable in the full space by enforcing reflection sym-metry.) In Figure 9b, multiplicity behavior exists for both the sym-metric and the (unstable) symmetric stretched steady states. Thebifurcation diagram indicates two saddle node bifurcations onthe symmetric branch and a hypothesized subcritical symmetry-breaking (pitchfork) bifurcation that also undergoes a saddle-node bifurcation, yielding the stable stretched asymmetric state.This is the simplest bifurcation scenario consistent with the com-putational results – the orange dashed curve has not been com-puted but is consistent with this scenario. For the larger value of c considered in 9c, there is no bistability of symmetric solutions,but the asymmetric branch still bifurcates subcritically and turnsaround to yield a stable asymmetric stretched state. If we furtherincrease c , bistability will vanish for the asymmetric solution andwe obtain a supercritical bifurcation. We present a case in Figure9d to show the existence of supercritical bifurcation.Figure 9e summarizes the behavior, showing the parameterspace where hysteresis arises for both the asymmetric and thesymmetric stretched states. The conformational hysteresis disap-pears faster as c increases for symmetric conformations than tiltedones. This section briefly turns to the case of biaxial extension, in whichwe shall see many of the features described above in planar exten-sion. To begin with, similar to planar extension, the strong com-pressive flow inward along the z direction suppresses any out-of-plane deformations, so at long times the deformations only occurin the x − y plane. Different from the planar extension case, dueto the radial symmetry of the biaxial flow in the x − y plane, thereis no issue of a preferred orientation for a square and a rectangle.Examples of stable steady states are illustrated in Figure 10.As in planar extension, if neo-Hookean elasticity is used, thereis again a singularity in stretching as Ca increases, as shown inFigure 11. In general, Ca c is larger compared to planar exten- Journal Name, [year], [vol.] , ig. 6 (a) Final stretched length l s vs. Ca for a disc with c = × − in planar extensional flow. Blue symbols and curves represent a compact finalconformation and red represent a stretched final conformation. The points filled with light color correspond to the transient evolution shown inFigure 6c-e. The black dashed line indicates the unstable steady state. The light blue dashed line shows the critical Ca with neo-Hookean elasticity. Ca c and Ca s are two critical Ca marked besides the arrows. (b) Final stretched length l s vs. Ca for a disc with c = × − in planar extensional flow.(c)-(e): Examples of stretched length evolution of a disc from either a stretched state or a compact state with c = × − in planar extensional flow:(c) Ca = . (d) Ca = . (e) Ca = . . Fig. 7 (a) Final stretched length l s vs. Ca for a disc in planar extensional flow with various c . All c have a similar (overlapped) compact branch. Theblue dashed line refers to Ca s with neo-Hookean elasticity. (b) Phase diagram of Ca c and Ca s vs. c in planar extensional flow. The gray area marksthe parameter regime where hysteresis can be observed. Journal Name, [year], [vol.] , ig. 8 (a) Bifurcation diagram of stretched length l s vs. Ca for a square in planar extensional flow, c = × − . (b) Phase diagram of Ca c and Ca s vs c for a square in planar extensional flow, with the gray area indicating the hysteresis region. Fig. 9 (a) Transient stretched length evolution of a rectangle in planar extensional flow, c = × − , Ca = . . (b) Bifurcation diagram of l s vs. Ca for a rectangle with c = × − .(c) Bifurcation diagram of l s vs. Ca for a rectangle sheet with c = × − . The data point filled with light colorcorresponds to the evolution shown in Figure 9a (d) Bifurcation diagram of l s vs. Ca for a rectangle with c = × − . (e) Phase diagram Ca c and Ca s vs c for a rectangle in planar extensional flow; here Ca s , a represents the stretched critical Ca for an asymmetric conformation, and Ca s , s representsthe stretched critical Ca for a symmetric conformation. The dark gray area indicates the hysteresis region for symmetric states and the light grayarea marks the hysteresis region without symmetry constraints, so the stretched state is tilted. Journal Name, [year], [vol.] , ig. 10 Steady-state conformation of a stable steady-state sheet withneo-Hookean elasticity in biaxial extensional flow at Ca = . : (a) disc.(b) square. (c) rectangle. Scale bar is 1 dimensionless length unit. sion. With the Yeoh model, bistability behavior similar to thatfound above is again observed. In biaxial extension, all threeshapes have a stretched state similar to its compact state, justmore stretched. No symmetry-breaking has been observed. Fig-ure 12a-c show bifurcation diagrams for the disc, square and rect-angle cases, respectively. Fig. 11
Steady-state stretched length l s vs. Ca for deformable sheetswith neo-Hookean elasticity in biaxial extensional flow. The data pointsfilled with light color correspond to the conformations shown in Figure10. The critical capillary numbers are: disc: Ca c = . , square: Ca c = . , rectangle: Ca c = . . The subplot indicates the lengthmeasured in each geometry. Figure 12d-f give the corresponding phase diagrams, showingwhere where hysteresis occurs. Compared with planar exten-sional flow results, the disc and rectangle in biaxial flow havelarger effective hysteresis regimes. Interestingly, the rectangle hasa substantially smaller hysteresis regime in c . This may arise fromthe smaller extent of the rectangle in one of the two characteristicstretching directions of biaxial extension. Recall that the presence of hysteresis in the coil-stretch transi-tion is due to the different hydrodynamic drag forces by the poly-mer in its coiled and stretched states . To assess the role ofhydrodynamic forces on bistability in the present case, we nowconsider an artificial model that weakens the hydrodynamic in-teraction between different parts of the sheet. In this model, we take the medium in which the sheet is moving to be a Brinkmanporous medium rather than a simple viscous fluid. In a porousmedium, hydrodynamic interactions decay quickly with distance,because the obstacles in the medium absorb momentum . Withthis model, the velocity field generated by point force F is de-termined not by the Stokeslet, but by the so-called Brinkmanlet: B ( x ) = πηα r (cid:20) (cid:0) − ( + α r ) e − α r (cid:1) xx r + (cid:16)(cid:16) + α r + α r (cid:17) e − α r − (cid:17)(cid:16) I − xx r (cid:17)(cid:21) . (13)Here a new length scale α − , the Brinkman screening length,arises. When α r (cid:28) , the Brinkmanlet reduces to the Stokeslet.When α r (cid:29) , the velocity decays as r − , as in a Darcy’s lawporous medium . In our model, we apply the Brinkmanlet,rather than the regularized Stokeslet, to the interactions betweendifferent nodes, and examine the effect of α . In dimensionlessform, α is scaled with sheet size a , so once α > , hydrodynamicwill be important only between nearby positions on the sheet.We take a square ( c = × − ) in both planar and biaxialextension as an example to illustrate how the HI influences thehysteresis behavior. Figure 13 shows the steady state stretchinglength vs. Ca for multiple choices of α . Compared with the case α = (Figure 8a for planar extension and Figure 12b for biaxialextension at c = × − ), the bistable region becomes narrowerand shifts towards smaller Ca as α increases. When α = , thehysteresis region vanishes and we obtain a smooth monostablecurve. We conclude that, as with linear polymers, bistable behav-iors arise from the hydrodynamic interactions between differentparts of the sheet surface. In this study, we systematically explored the dynamics of freelysuspended elastic sheets with three different rest shapes in planarand biaxial extensional flows. In both flow fields, the sheet alwaystakes on a flat steady conformation in the stretching plane. If neo-Hookean elasticity is applied, we observed a weakly stretched“compact" sheet conformation at low Ca , while above a critical Ca c the degree of stretching diverges. With a nonlinear Yeoh elas-ticity model, this singularity vanishes, but for sufficiently smallnonlinearity parameter c , we observe a “compact-stretched" tran-sition, where the sheet suddenly becomes highly stretched oncea critical capillary number is exceeded. The discontinuity instretched length marks a bistable or hysteretic regime, where thesheet has either a compact or a stretched conformation based onits deformation history. This behavior occurs for all three restshapes considered. The origin of such bistable behavior is a re-sult of hydrodynamic interaction between different parts of thesheet. The dynamics become monostable if hydrodynamic inter-actions are artificially screened out. The highly nonlinear behav-ior predicted here may play a substantial role in the behavior andperformance of flow processes involving soft sheet-like particles. Journal Name, [year], [vol.] ,,
Steady-state stretched length l s vs. Ca for deformable sheetswith neo-Hookean elasticity in biaxial extensional flow. The data pointsfilled with light color correspond to the conformations shown in Figure10. The critical capillary numbers are: disc: Ca c = . , square: Ca c = . , rectangle: Ca c = . . The subplot indicates the lengthmeasured in each geometry. Figure 12d-f give the corresponding phase diagrams, showingwhere where hysteresis occurs. Compared with planar exten-sional flow results, the disc and rectangle in biaxial flow havelarger effective hysteresis regimes. Interestingly, the rectangle hasa substantially smaller hysteresis regime in c . This may arise fromthe smaller extent of the rectangle in one of the two characteristicstretching directions of biaxial extension. Recall that the presence of hysteresis in the coil-stretch transi-tion is due to the different hydrodynamic drag forces by the poly-mer in its coiled and stretched states . To assess the role ofhydrodynamic forces on bistability in the present case, we nowconsider an artificial model that weakens the hydrodynamic in-teraction between different parts of the sheet. In this model, we take the medium in which the sheet is moving to be a Brinkmanporous medium rather than a simple viscous fluid. In a porousmedium, hydrodynamic interactions decay quickly with distance,because the obstacles in the medium absorb momentum . Withthis model, the velocity field generated by point force F is de-termined not by the Stokeslet, but by the so-called Brinkmanlet: B ( x ) = πηα r (cid:20) (cid:0) − ( + α r ) e − α r (cid:1) xx r + (cid:16)(cid:16) + α r + α r (cid:17) e − α r − (cid:17)(cid:16) I − xx r (cid:17)(cid:21) . (13)Here a new length scale α − , the Brinkman screening length,arises. When α r (cid:28) , the Brinkmanlet reduces to the Stokeslet.When α r (cid:29) , the velocity decays as r − , as in a Darcy’s lawporous medium . In our model, we apply the Brinkmanlet,rather than the regularized Stokeslet, to the interactions betweendifferent nodes, and examine the effect of α . In dimensionlessform, α is scaled with sheet size a , so once α > , hydrodynamicwill be important only between nearby positions on the sheet.We take a square ( c = × − ) in both planar and biaxialextension as an example to illustrate how the HI influences thehysteresis behavior. Figure 13 shows the steady state stretchinglength vs. Ca for multiple choices of α . Compared with the case α = (Figure 8a for planar extension and Figure 12b for biaxialextension at c = × − ), the bistable region becomes narrowerand shifts towards smaller Ca as α increases. When α = , thehysteresis region vanishes and we obtain a smooth monostablecurve. We conclude that, as with linear polymers, bistable behav-iors arise from the hydrodynamic interactions between differentparts of the sheet surface. In this study, we systematically explored the dynamics of freelysuspended elastic sheets with three different rest shapes in planarand biaxial extensional flows. In both flow fields, the sheet alwaystakes on a flat steady conformation in the stretching plane. If neo-Hookean elasticity is applied, we observed a weakly stretched“compact" sheet conformation at low Ca , while above a critical Ca c the degree of stretching diverges. With a nonlinear Yeoh elas-ticity model, this singularity vanishes, but for sufficiently smallnonlinearity parameter c , we observe a “compact-stretched" tran-sition, where the sheet suddenly becomes highly stretched oncea critical capillary number is exceeded. The discontinuity instretched length marks a bistable or hysteretic regime, where thesheet has either a compact or a stretched conformation based onits deformation history. This behavior occurs for all three restshapes considered. The origin of such bistable behavior is a re-sult of hydrodynamic interaction between different parts of thesheet. The dynamics become monostable if hydrodynamic inter-actions are artificially screened out. The highly nonlinear behav-ior predicted here may play a substantial role in the behavior andperformance of flow processes involving soft sheet-like particles. Journal Name, [year], [vol.] ,, ig. 12 Bifurcation diagrams of l s vs. Ca for sheets of different shapes in biaxial extensional flow: (a) disc: c = × − . (b) square: c = × − ,(c) rectangle: c = × − . The light blue dashed line refers to Ca c with neo-Hookean elasticity. Phase diagram of the hysteresis regime in biaxialextensional flow: (d) disc, (e) square, (f) rectangle. Fig. 13
Bifurcation diagrams of l s vs. Ca for a square sheet ( c = × − ) with increasing α : (a) planar extensional flow. (b) biaxial extensionalflow. Conflicts of interest
There are no conflicts to declare.
Acknowledgements
This material is based on work supported by the National Sci-ence Foundation under grant No. CBET-1604767. This work
10 | 1–11
Journal Name, [year], [vol.] , sed the Extreme Science and Engineering Discovery Environ-ment (XSEDE) SDSC Dell Cluster with Intel Haswell Pro-cessors (Comet) through allocations TG-CTS190001 and TG-MCB190100.
Notes and references
Advanced materials , 2013, , 1248–1253.2 B. Ni, M. Huang, Z. Chen, Y. Chen, C.-H. Hsu, Y. Li, D. Pochan,W.-B. Zhang, S. Z. Cheng and X.-H. Dong, Journal of the Amer-ican Chemical Society , 2015, , 1392–1395.3 S. So and R. C. Hayward,
ACS Applied Materials & Interfaces ,2017, , 15785–15790.4 M. Jamal, A. M. Zarafshar and D. H. Gracias, Nature commu-nications , 2011, , 1–6.5 A. M. Hubbard, R. W. Mailen, M. A. Zikry, M. D. Dickey andJ. Genzer, Soft Matter , 2017, , 2299–2308.6 J.-H. Na, A. A. Evans, J. Bae, M. C. Chiappelli, C. D. San-tangelo, R. J. Lang, T. C. Hull and R. C. Hayward, AdvancedMaterials , 2015, , 79–85.7 A. W. Feinberg, A. Feigel, S. S. Shevkoplyas, S. Sheehy, G. M.Whitesides and K. K. Parker, Science , 2007, , 1366–1370.8 K.-Y. Chun, S. H. Kim, M. K. Shin, Y. T. Kim, G. M. Spinks,A. E. Aliev, R. H. Baughman and S. J. Kim,
Nanotechnology ,2013, , 165401.9 Y. Zhao, J. Guan and S. J. Wu, Macromolecular Rapid Commu-nications , 2019, , 1900389.10 A. K. Gaharwar, S. A. Dammu, J. M. Canter, C.-J. Wu andG. Schmidt, Biomacromolecules , 2011, , 1641–1650.11 J. Liu, C. Chen, C. He, J. Zhao, X. Yang and H. Wang, ACSnano , 2012, , 8194–8202.12 S. M. Mehta, T. Jin, I. Stanciulescu and K. J. Grande-Allen, Acta Biomaterialia , 2018, , 52–62.13 P. De Gennes, The Journal of Chemical Physics , 1974, ,5030–5042.14 C. M. Schroeder, H. P. Babcock, E. S. Shaqfeh and S. Chu, Science , 2003, , 1515–1519.15 G. Fuller and L. Leal,
Journal of Non-Newtonian Fluid Mechan-ics , 1981, , 271–310.16 M. H. Nafar Sefiddashti, B. J. Edwards and B. Khomami, TheJournal of Chemical Physics , 2018, , 141103.17 V. Kantsler, E. Segre and V. Steinberg,
Physical review letters ,2008, , 048101.18 V. Narsimhan, A. P. Spann and E. S. Shaqfeh,
Journal of FluidMechanics , 2015, , 1–26. 19 D. Kumar, C. M. Richter and C. M. Schroeder,
Soft Matter ,2020, , 337–347.20 M. J. Shelley and J. Zhang, Annual Review of Fluid Mechanics ,2011, , 449–465.21 S. Alben and M. J. Shelley, Physical review letters , 2008, ,074301.22 Y. Xu and M. J. Green,
The Journal of chemical physics , 2014, , 024905.23 Y. Xu and M. J. Green,
Journal of Polymer Science Part B: Poly-mer Physics , 2015, , 1247–1253.24 S. Gravelle, C. Kamal and L. Botto, The Journal of ChemicalPhysics , 2020, , 104701.25 C. Kamal, S. Gravelle and L. Botto,
Nature communications ,2020, , 1–10.26 S. Dutta and M. D. Graham, Soft Matter , 2017, , 2620–2633.27 O. H. Yeoh, Rubber Chemistry and technology , 1993, , 754–771.28 D. A. Fedosov, B. Caswell and G. E. Karniadakis, ComputerMethods in Applied Mechanics and Engineering , 2010, ,1937–1948.29 D. A. Fedosov,
Multiscale modeling of blood flow and soft mat-ter , Citeseer, 2010.30 S. P. Timoshenko and S. Woinowsky-Krieger,
Theory of platesand shells , McGraw-hill, 1959.31 X. Liu, T. H. Metcalf, J. T. Robinson, B. H. Houston andF. Scarpa,
Nano letters , 2012, , 1013–1017.32 J. Charrier, S. Shrivastava and R. Wu, The Journal of StrainAnalysis for Engineering Design , 1989, , 55–74.33 V. Pappu and P. Bagchi, Comput Biol Med , 2008, , 738–753.34 A. Kumar and M. D. Graham, Journal of ComputationalPhysics , 2012, , 6682–6713.35 R. Cortez, L. Fauci and A. Medovikov,
Physics of Fluids , 2005, , 031504.36 J. P. Hernández-Ortiz, J. J. de Pablo and M. D. Graham, Phys-ical review letters , 2007, , 140602.37 M. D. Graham, Microhydrodynamics, Brownian motion, andcomplex fluids , Cambridge University Press, 2018, vol. 58.38 J. Happel and H. Brenner,
Low Reynolds number hydrodynam-ics: with special applications to particulate media , Springer Sci-ence & Business Media, 2012, vol. 1.39 J. Towns, T. Cockerill, M. Dahan, I. Foster, K. Gaither,A. Grimshaw, V. Hazlewood, S. Lathrop, D. Lifka, G. D. Pe-terson et al. , Computing in science & engineering , 2014, ,62–74. Journal Name, [year], [vol.] ,,