You are here
Effect of constraint on swelling of hydrogels and formation of surface creases
Inspired by recent works by Wei Hong , Xuanhe Zhao, Zhigang Suo, and their coworkers, we started a project on hydrogels, with particular interest in various instability patterns observed in experiments. The attachment is our first manuscript on this subject. Through this work we hope to achieve the following:
(1) Learn the theoretical approach developed by Hong et al. as well as its theoretical background in rubber elasticity and polymer solutions. As a result, we re-formulate the theory in a variational form in the attached manuscript.
(2) Develop explicit formula for the true stress and tangent modulus for a specific material model of hydrogel. This is necessary for a finite element implementation using a UMAT subroutine in ABAQUS, which differs from the UHYPER implementation by Hong et al. The present implementation removes a restriction that the initial state must be isotropic in the previous UHYPER implementation. An anisotropic initial state becomes necessary for swelling of hydrogels under constraints.
(3) Study the effect of geometric constraint on swelling of surface-attached hydrogel lines with different aspect ratios. It is found that the effect of constraint can be tuned between two homogeneous limits. More interestingly, beyond a critical aspect ratio, crease-like surface instability occurs upon swelling. It is suggested that a critical condition for creasing on the surface of hydrogel can be expressed in terms of two material parameters of the hydrogel and the external chemical potential.
In addition, a few interesting points are noted through this study:
(1) While the Flory-Rehner free energy function has been widely used by polymer scientists, the volumetric term in the elastic free energy differs (by a factor of 2) from that commonly used by mechanicians. According to Treloar (1975), the debate over Flory's statistical mechanics model for rubber elasticity was inconclusive, and his book used exclusively a different free energy function. The difference results in different predictions for the osmotic pressure and the chemical potential in hydrogels, which may be significant and may be not.
(2) The phenomenon of surface creasing was observed much earlier than Tanaka (1987). In their experiments with rubber vulcanizates, Southern and Thomas (1965) observed surface instability (which they called wrinkles) when the degree of swelling is greater than a critical value.
(3) Due to the surface instability along with self-contact and possible singularity at the tip of a crease, it remains a numerical challenge to simulate formation and evolution (growth) of surface creases. 3D simulations are extremely demanding with the current computational facility and method. We hope to overcome this barrier in the future so that the numerical results can be directly compared to experimental observations.
We welcome constructive and instructive comments and suggestions!
RH
Update on February 14, 2010: The attached manuscript has been accepted for publication in Journal of Applied Mechanics. A follow-up study focusing on swell induced surface instabilty and creasing has also been submitted for review.
| Attachment | Size |
|---|---|
 hydrogel2009.pdf | 791.91 KB |
- Rui Huang's blog
- Log in or register to post comments
- 31842 reads
Comments
Interesting work
Interesting work, Rui! I just had time to read it briefly and here are my 3 questions (comments):
1) In equilibrium, the field equation of chemical potential, (22-2), is a set of homogeneous first order ODEs, the only solution is mu=mu^, no boundary condition is needed. The boundary conditions are valid for dynamic cases though.
2) The anisotrpoic initial state, Eq. (37) can be obtained by prescribing a uniform deformation to an isotropic free-swelling state, as long as the gel itself (or the free-energy function) is isotropic. Maybe it could be useful for intrinsically anisotropic gel, i.e. the gels that swells differently in different direction.
In section 4.2, the example you calculated (lambda1=lambda2=1) is different from the one you compared to(lambda1=lambda2=1.5). A comparable calculation can be performed using the UHYPER program by starting from the lambda=1.5 isotropic initial state, and biaxially compress it to the lambda1=lambda2=1 state. I believe the result will be identical to yours as long as both programs are correct. It might be better to compare to that case for verification.
3) What is the critical swelling ratio for the formation of creases in your case? I remember in similar caculations we did (for a gel in an equilibrium state), the critical value is actually higher than the one we reported recently (for instantaneous response). It would be nice to see an dynamic analysis.
Please correct me if I made any mistake.
Wei
chemical boundary condition etc
Dear Wei,
Thanks for your comments and questions. My answers follow.
1) The boundary conditions come out naturally from the variational principle, in which the concentration C is treated as the state variable (like displacement). If one chooses to solve the coupled field equations for the displacement and concentration simultaneously, both the mechanical and chemical boundary conditions would be needed. In the present method (following your work), both the equilibrium equation and the boundary condition for the chemical potential are automatically satisfied.
2) Yes you can deform an isotropic free-swelling state into an anisotropic state. But with UHYPER, you have to start with the isotropic initial state in your numerical simulation. For the example in Section 4.2, you could start with an isotropic state (lambda1 = lambda2 = lambda3 = 1.5) and gradually increase the pressure (or displacement) in the 1 and 3 directions to deform the gel until it reaches the laterally constrained state (lambda1 = lambda3 = 1), then followed by increasing the chemical potential. So, you need to run two numerical steps in ABAQUS and the boundary conditions at the bottom face will have to be different in the two steps. If all done correctly, the results should be identical to the analytical solution with lambda1 = lambda3 = 1 as well as the numerical solution with UMAT. Wo do mention in the text that the UHYPER results plotted in Fig. 3 correspond to a different analytical solution (with lambda1 = lambda3 = 1.5).
3) Since the swelling deformation is inhomogeneous for the surface-attached hydrogel lines, we do not think the swelling ratio can be used as a critical condition for creasing. As you see in Fig. 8, the volume swelling ratio is high for low aspect ratios (W/H), where no creasing occurs. Our analysis is not a dynamic analysis (no time is involved), but it is not the instantaneous response as you called it (by comparing to an incompressible elastomer). My student has run many simulations with different aspect ratios and different combinations of the material parameters, with which we are trying to sort out the critical condition for creasing in hydrogels. I would be surprised if it turns out identical to that for a neo-Hookean elastomer.
RH
Chemical potential B.C.s
Thank you for the explaination, Rui!
I still don't quite agree with 1). Even when C is the independent variable, it still satisfies an algebraic equation: mu(C)=mu^, no additional boundary condition is needed.
Your Eq. (23-2) is not a valid boundary condition, for:
a) If mu=mu^ is a boundary condition, then I can specify mu=mu_1 on one side, and mu=mu_2 on the other side, then no solution could be found. This is all because the field equation is a first order differential equation, which can never consist a boundary value problem. You are only allowed to specify the condition at 1 point (for 1 continuous domain), and integrate from that point.
b) The condition N_K deltaI_K =0 does not mean zero flux. Actually any constant flux would satisfy this "boundary condition".
This set of boundary condition is suitable for a steady-state case, rather then the equilibrium state described by Eq. (22-2).
Wei
Re: Chemical potential B.C.s
Wei,
I do not disagree with you, but think about the following:
a) The mechanical equilibrium equation (22-1) is also a first order different equation in terms of stress. The only difference is that stress is a second-order tensor (while chemical potential is a scalar) and we cannot fully specify the stress state at the boundary.
b) As we mentioned in the paper, no equilibrium state exists if the gel is in contact with an environment of different chemical potentials, just like heat transfer when the temperature is not uniform at the boundary of a body. This does not mean that the boundary condition is invalid.
c) Zero flux is a special case of prescribed constant flux, just like zero displacement as to prescribed displacement (essential boundary condition). Suppose the gel is only partly immersed in the solvent. There would be no molecular transport (flux) through the boundary where it is not in contact with the solvent.
d) All the boundary conditions are derived from the general variational statement in Eq. (18).
Thanks.
RH
Variational Statement
Rui,
I think I found the origin of the apparent disagreement:
In Eq. (18), delta_I_K is not really an arbitrary vector. It must satisfy mass conservation which is not present in this equation.
delta_I_K is like delta_strain in mechanics, the latter must satisfy the compatibility equations. We don't usually do variation on strain, but convert it into delta_displacement with Green't theorem.
I think we need to convert delta_I_K into delta_C or something alike here as well.
Wei
Re: Variational statement
Wei,
Eqs. (13) and (14) relate the variation in I to the variation in C as a result of mass conservation. At the end of Section 2.1, we emphasize that the variational statements include two parts, one for the free energy and one for mass conservation.
Best regards,
RH
Re: Variational statement
Rui,
I was refering to this statement in your paper: (after Eq. 18)
"For Eq. (18) to hold for any arbitrary variations, it necessarily requires that ... "(the strong form + B.C.s)
Here, delta_I_K can not be arbitrary, and we can not obtain (19)&(20) from (18).
Wei
Re: Variational statement
Wei,
I disagree. Eq. (18) is reached with the condition of mass conservation. From there on, it is the standard variational analysis. I do not see any other alternatives. Remember, at this point, the assumption of molecular incompressiblity has not been imposed. Consequently, the concentration C is considered independent from the deformation.
RH
Re: Variational statement
Rui,
Please rethink about it. I don't think the derivation after Eq. (18) is standard variational analysis.
If a variational analysis is OK, we can do a similar analysis like int(stress*delta_strain) and arrive at a result that stress=0 for any case when no body force is present.
Again, I_K is something like strain, not displacement.
Wei
mass transport in variational analysis
Wei,
Looks like we are both stubborn.
I hope someone else (Zhigang would be the best) could jump in to help us. I believe similar approaches have been used for other type of mass transport, such as surface diffusion. Perhaps this paper by Zhigang and Wei Lu could help: Z. Suo and W. Lu, "Forces that drive self-assembly on solid surfaces", Journal of Nanoparticles Research 2, 333-344 (2000). Pay particular attention to Eq. (9) and the paragraph after it for the equilibrium condition.
Let's hope our discussion will lead to a better understanding of the variational approach, one way or another.
Regards,
RH
variational statement for a gel in a state of equilibrium
Dear Rui:
I'm not sure I can really be of much help, but only wish that we can all be in the same room, in front of a whiteboard. In any event, I usually lost to Wei when he and I argued, and he and I argued all the time!
Perhaps this line of thoughts may help us re-orient the discussion. In this IJSS paper, we have tried to recap the variational statement for a gel in equilibrium with a solvant and a set of mechanical loads. A longer derivation is contained in my notes on poroelasticity, which was updated in the Spring of 2009 when I taught the class for the second time.
Are you OK with these derivations? Or do you see aspects that make you uncomfortable and wish to address in your re-formulation?
variational approach to derive equilibrium equations
Dear Zhigang,
Thank you for joining the discussion. I am comfortable with your derivation in the IJSS paper as well as your note. In my re-formulation, I try to derive the equilibrium equations and corresponding boundary conditions all together from the variational statements, without knowing or assuming anything a prior (other than delta(G) = 0). For example, before I reach the equilibrium equations, I do not know the chemical potential has to be a constant in the gel. All I know is the external chemical potential and the internal free energy function U. I do not assume that the external chemical potential is a constant either. Only later (after I have reached the equilibrium equations and the boundary conditions) do I conclude that an equilibrium state exists only when the external chemical potential is a constant. Consequently I leave the chemical potential inside the surface integral for the chemical work and do not have the volume integral in Eq. (3) of your IJSS paper.Eventually, the same equilibrium conditions are reached. It may be just a matter of style that I like to clarify when an equilibrium state exists and when it does not and why, all from the variational analysis.
The debate with Wei Hong is perhaps more fundamental than it appears. Our disagreement is whether the variation in the flux components (I_K) can be considered independent like the displacement components.
RH
Alternative derivation
Rui,
An alternative derivation is here. (It seems that a comment cannot have attachement on imechanica)
http://www.public.iastate.edu/~whong/variational.pdf
The compoenents of the flux vector are not indepedent variables. One can not simply say that they are arbitrary so that the coefficients must be 0
Physical meaning of the variational approach
Rui,
After discussion, I gave the problem another thought. And here is my another argument why delta_I_K can not be used in the same way as the displacement in the variational approach:
Physically, the variational principle corresponds to an energy mimimum, so that delta_G/delta_x = 0. Here this x should be a state variable, so that the change of x would represent another state. e.g. change in the displacement field would result in a different deformation state; change in the concentration field would also result in a different state in terms of solvent content.
However, the flux-vector field is not a state variable. Change in the flux may or may not represent another state. (I can have a same piece of gel in exactly the same deformation state and solvent content, but have the same amount of molecules traveling in from one side and out from the other side. The flux would be different from the case when no flux is present, while the gel is in a same state discribed by F and C or mu.) In fact the total free energy is not a functional of I_K.
You might argue that the change in the flux might still introduce the change in the total free energy, as in the approach we used for the kinetic model. While the corresponding change in total free energy, if any, is actually the energy dissipation. From the second law of thermodynamics, we know the dissipation can not be negative, (for arbitrary delta_I_K). However, this statement can not be extended to the equilibrium case -- delta_G = 0 for arbitrary delta_I_K simply means the system has no dissipation and/or the process is an reversible process -- it is not the equilibrium condition.
Hope this argument makes more physical sense to you rather than saying that delta_I_K is not independent. I believe I convinced myself, let me know if I have convinced you this time.
Wei
are you sure?
"change in the displacement field would result in a different deformation state"? How about rigid-body motion?
RH
I am sure flux is not a state variable
Rigid-body motion does actually change the state of the system. (sorry about the wrong word deformation state.)
If we have a displacement boundary condition, rigid-body motion is not allowed in delta_u_i.
If we have a traction boundary condition, rigid-body motion changes the potential of the external load.
If the state of the system does not change with respect to a rigid-body motion, we usually say that the problem is ill constrained. In this case, delta_G/delta_u_i=0 will give you a group of solution rather than a single one.
These are not the case for our flux, even in a steady state. Flux describes a process, rather than a state.
Wei
Rui, Very interesting
Rui,
Very interesting work. I went through it quickly and hope to study it in details later. I am especially impressed by your simulation on crease formation. You mentioned, "surface creases form automatically in the numerical simulations for hydrogel lines beyond the critical aspect ratio, without any perturbation." Could you please give more details about this point? For example, how to deal with the "snap-through" instability during crease formation? Thanks.
XH
formation of creases
Dear Xuanhe,
Thank you for your comments. In our simulations, when the line aspect ratio is small (W/H < 10), ABAQUS runs smoothly till the end without any problem, and no creases form. At W/H = 12 and beyond, the ABAQUS calculation becomes unstable and would abort (due to divergence, I think) before the chemical potential reaches zero. We take this as an indication of the onset of crease instability. To simulate the growth of creases beyond the critical point, we had to introduce some stabilization scheme (provided within ABAQUS). This remains a numerical issue to be addressed, as we do not fully understand how the stabilization method works in ABAQUS and it appears artificial to me (with some viscosity introduced to dissipate energy). We are still working on this matter and hope to be able to find a reasonable solution soon.
RH
instability of a crease
The instability of the crease may be understood in the scaling relation, eq. (2), in our short paper on crease.
This equation is valid when the size of the crease L is much smaller than any other size in the probelm, so that L is the only relevant length. According to this scaling relation, when the critical condition is passed, the reduction in the free energy will be larger for a larger crease.
When the size of the crease is significant compared to the thickness of the layer, Eq. (2) is no longer valid, and the crease may stabilize again. Thus, crease is a snap-through instability.
Re: instability of a crease
Dear Zhigang,
I think I understand your point. I actually like the scaling equation in your paper very much. But I have concerns in applying it to hydrogels. For example, what is the homogeneous swelling state for a hydrogel half-space? As for whether creasing is a snapthough instability, it may also depend on how it is loaded (load control vs displacement control). I still remember the video shown by Ryan Hayward more than a year ago at the ACS meeting: creasing is a dynamic (or kinetic) process, much like growth of channel cracks in thin films.
RH
kinetics of creasing in a gel as it swells
I agree with your concern about the kinetic process of creasing in a gel as it swells. Ryan came to Harvard several times to discuss with the group about his experimental observations. Several students and I have also talked about it on a number of occasions. We have not reached a clean picture yet.
The scaling relation (2) was for elastomers. The relation is also correct for a swollen gel in equilibrium with a solvent. The function f should have another list of independent variables. How to apply the scaling relation to a swelling gel needs more thought.
A more complete formulation of the problem should involve two time-dependent processes: diffusion and inertia. One then starts with a initial defect, and evolves the initial-boundary value problem. But one wishes to gain qualitative understanding from this detailed model.
critical swelling ratio for hydrogel thin films
Inspired by the discussions here with Wei Hong, Xuanhe, and Zhigang, we have performed a linear perturbation analysis for hydrogel thin films swelling under lateral constraints. The approach is similar to that by Biot (1963) for instability of rubber under compression. It is found that the homogeneous swelling of a hydrogel thin film (swelling in the thickness direction only) becomes unstable at a critical swelling ratio, and the critical swelling ratio varies between 2.0 and 3.5, depending on the two dimensionless material parameters in the model. The range of the critical swelling ratio agrees well with the scattering experimental data by several groups. A manuscript is in preparation. (update 12 February 2010: the manuscript has been submitted for review.)
RH