2.1 Governing equations for the gel

The equations for the gel are formulated in terms of Eulerian coordinates $\boldsymbol{x} = r \boldsymbol{e}_r(\theta ) + z \boldsymbol{e}_z$ associated with the current state of the system, where $r$ , $\theta$ and $z$ denote the radial, angular and axial coordinates, respectively, and $\boldsymbol{e}_r$ , $\boldsymbol{e}_\theta$ and $\boldsymbol{e}_z$ are the corresponding cylindrical basis vectors. An Eulerian coordinate system enables the equations to be written in a physically intuitive way and it facilitates coupling the gel and bath models via boundary conditions. A detailed account of Eulerian-based hydrogel modelling is provided by Bertrand et al. [Reference Bertrand, Peixinho, Mukhopadhyay and MacMinn3]. We let $\boldsymbol{X} = R(r) \boldsymbol{E}_{R}(\theta ) + Z(z) \boldsymbol{E}_Z$ denote the Lagrangian coordinates associated with the stress-free reference configuration, which we assume is a dry gel. Here, $R$ , $\Theta = \theta$ and $Z$ represent the radial, angular, and axial Lagrangian coordinates, and $\boldsymbol{E}_R$ , $\boldsymbol{E}_{\Theta }$ and $\boldsymbol{E}_Z$ are the basis vectors.

The deformation gradient tensor $\boldsymbol{\mathsf{F}}$ describes the distortion of material elements relative to the dry state of the gel. For an axisymmetric geometry which remains cylindrical, the deformation gradient tensor can be written as $\boldsymbol{\mathsf{F}} = \lambda _r\, \boldsymbol{e}_r\otimes \boldsymbol{E}_R + \lambda _\theta \, \boldsymbol{e}_\theta \otimes \boldsymbol{E}_\Theta + \lambda _z \boldsymbol{e}_z \otimes \boldsymbol{E}_Z$ , where

(2.2) \begin{align} \lambda _r = \left (\frac{{\textrm{d}} R}{{\textrm{d}} r} \right )^{-1}, \quad \lambda _\theta = \frac{r}{R}, \quad \lambda _z = \left (\frac{{\textrm{d}} Z}{{\textrm{d}} z}\right )^{-1} \end{align}

denote the radial, orthoradial and axial stretches, respectively. The axial stretch $\lambda _z$ is imposed, whereas the radial and orthoradial stretches $\lambda _r$ and $\lambda _\theta$ are unknown and must be solved for. The determinant $J = \det \boldsymbol{\mathsf{F}} = \lambda _r \lambda _\theta \lambda _z$ characterises volumetric changes in material elements. Both the polymers and the imbibed salt solution are assumed to be incompressible. As a result, any volumetric change in a solid element must be due to a variation in the amount of fluid contained within that element. This leads to the so-called molecular incompressibility condition

(2.3) \begin{align} J = \left (1 - \sum _{m\in \mathbb{M}} \phi _m\right )^{-1}, \end{align}

where $\phi _k$ represent the volume fraction of species $k$ . Since $J$ describes the volume of swollen material elements relative to their dry volume, we also refer to it as the swelling ratio. It will be convenient to formulate the incompressibility condition as

(2.4) \begin{align} R \frac{{\textrm{d}} R}{{\textrm{d}} r} = \frac{\lambda _z r}{J}, \end{align}

where $J$ is given by (2.3).

The chemical potentials of the mobile species, $\mu _m$ , can be written as

(2.5a) \begin{align} \mu _s &= \Pi _s +{\mathcal{G}} p - \frac{\omega ^2}{r}\frac{{\textrm{d}} }{{\textrm{d}} r}\left (r \frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right ), \end{align}

(2.5b) \begin{align} \mu _{\pm } &= \Pi _{\pm } +{\mathcal{G}} p + z_{\pm } \Phi, \end{align}

where $z_{\pm }$ is the valence of the ions, $p$ is the mechanical pressure in the gel and $\Pi _m$ are osmotic pressures defined as

(2.6a) \begin{align} \Pi _s &= \log \phi _s + \chi \,J^{-1}(1 - \phi _s) + J^{-1}, \end{align}

(2.6b) \begin{align} \Pi _\pm &= \log \phi _\pm + J^{-1}(1 - \chi \phi _s). \end{align}

Here, $\chi$ is the Flory interaction parameter, which describes (unfavourable) enthalpic interactions between the solvent molecules and the polymers. Due to our assumption that the system is in equilibrium, the chemical potentials are spatially uniform. Hence, $\mu _m$ are constants that will be specified below.

The electric potential satisfies Poisson’s equation

(2.7) \begin{align} -\frac{\varepsilon ^2}{r}\frac{{\textrm{d}} }{{\textrm{d}} r}\left (r \frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right ) = z_+\phi _+ + z_{-} \phi _- + \frac{z_f \varphi _f}{J}, \end{align}

where $\varphi _f$ is the nominal volume fraction of fixed charges on the polymer network and $z_f$ denotes the valence of these charges. We will focus on cationic gels with positive fixed charges, $z_f \gt 0$ .

The conservation of linear momentum in the gel leads to

(2.8) \begin{align} \frac{{\textrm{d}} \mathsf{T}_{rr}}{{\textrm{d}} r} + \frac{\mathsf{T}_{rr} - \mathsf{T}_{\theta \theta }}{r} = 0, \end{align}

where $\mathsf{T}_{rr}$ and $\mathsf{T}_{\theta \theta }$ are the radial and orthoradial components of the Cauchy stress tensor. These stresses can be expressed as

(2.9a) \begin{align} \mathsf{T}_{rr} &= \mathsf{T}_{e,rr} + \mathsf{T}_{K,rr} + \mathsf{T}_{M,rr} - p, \end{align}

(2.9b) \begin{align} \mathsf{T}_{\theta \theta } &= \mathsf{T}_{e,\theta \theta } + \mathsf{T}_{K,\theta \theta } + \mathsf{T}_{M,\theta \theta } - p. \end{align}

The first contributions, $\mathsf{T}_{e,rr}$ and $\mathsf{T}_{e,\theta \theta }$ , represent elastic stresses, which are calculated by assuming the polymer network behaves as a neo-Hookean material. This leads to

(2.10a) \begin{align} \mathsf{T}_{e,rr} &= J^{-1}(\lambda _r^2 - 1), \end{align}

(2.10b) \begin{align} \mathsf{T}_{e,\theta \theta } &= J^{-1}(\lambda _\theta ^2 - 1). \end{align}

The second and third contributions to the Cauchy stresses in (2.9) correspond to Korteweg ( $\mathsf{T}_K$ ) and Maxwell ( $\mathsf{T}_M$ ) stresses, respectively, which capture the forces generated within the bulk of the gel due to internal interfaces and electric fields. The radial and orthoradial components of these stresses are given by

(2.11a) \begin{align} \mathsf{T}_{K,rr} &={\mathcal{G}}^{-1} \omega ^2 \left [\frac{\phi _s}{r}\frac{{\textrm{d}} }{{\textrm{d}} r}\left (r \frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right ) - \frac{1}{2}\left (\frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right )^2\right ], \end{align}

(2.11b) \begin{align} \mathsf{T}_{K,\theta \theta } &={\mathcal{G}}^{-1} \omega ^2 \left [\frac{\phi _s}{r}\frac{{\textrm{d}} }{{\textrm{d}} r}\left (r \frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right ) + \frac{1}{2}\left (\frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right )^2\right ], \end{align}

(2.11c) \begin{align} \mathsf{T}_{M,rr} &= \frac{1}{2}{\mathcal{G}}^{-1} \varepsilon ^2\left (\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right )^2, \end{align}

(2.11d) \begin{align} \mathsf{T}_{M,\theta \theta } &= -\frac{1}{2}{\mathcal{G}}^{-1} \varepsilon ^2\left (\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right )^2. \end{align}

The final contribution to the Cauchy stresses represents the stress induced by the fluid pressure.

2.2 Governing equations for the bath

The spatially uniform chemical potentials of the solvent and ions are

(2.12a) \begin{align} \mu _s &= \log \phi _s + \epsilon _r \varepsilon ^2 p, \end{align}

(2.12b) \begin{align} \mu _{\pm } &= \log \phi _{\pm } + \epsilon _r \varepsilon ^2 p + z_{\pm } \Phi, \end{align}

where $\epsilon _r = \epsilon ^{\text{bath}}/ \epsilon ^{\text{gel}}$ . The electric potential satisfies

(2.13) \begin{align} -\frac{\epsilon _r \varepsilon ^2}{r}\frac{{\textrm{d}} }{{\textrm{d}} r}\left (r \frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right ) = z_+\phi _+ + z_{-} \phi _-. \end{align}

Conservation of linear momentum in the bath implies that

(2.14) \begin{align} \frac{{\textrm{d}} \mathsf{T}_{rr}}{{\textrm{d}} r} + \frac{\mathsf{T}_{rr} - \mathsf{T}_{\theta \theta }}{r} = 0, \end{align}

where the components of the Cauchy stress tensor are

(2.15a) \begin{align} \mathsf{T}_{rr} &= \mathsf{T}_{M,rr} - p, \end{align}

(2.15b) \begin{align} \mathsf{T}_{\theta \theta } &= \mathsf{T}_{M,\theta \theta } - p. \end{align}

The Maxwell stresses are given by

(2.16a) \begin{align} \mathsf{T}_{M,rr} &= \frac{1}{2}\left (\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right )^2, \end{align}

(2.16b) \begin{align} \mathsf{T}_{M,\theta \theta } &= -\frac{1}{2}\left (\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right )^2. \end{align}

2.3 Boundary conditions

At the centre of the gel, we impose

(2.17) \begin{align} R|_{r=0} = 0, \quad \left .\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right |_{r=0} = 0, \quad \left .\frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right |_{r=0} = 0. \end{align}

The first condition ensures that the origin in Lagrangian coordinates is mapped to the origin in Eulerian coordinates. The second and third can be viewed as symmetry conditions. The boundary condition on $\phi _s$ is needed due to the presence of a second derivative in the expression for the chemical potential of solvent in the gel (2.5a).

Far from the bath, $r \to \infty$ , we set the electric potential to $\Phi ^{\text{bath}}$ and the pressure to zero. In addition, we assume that the volume fraction of cations has been fixed to $\phi ^{\text{bath}}_{+}$ . Assuming electroneutrality then requires that the volume fraction of anions is given by $\phi ^{\text{bath}}_{-} = |z_{+}/ z_{-}| \phi ^{\text{bath}}_{+}$ . The far-field solvent fraction is then fixed at $\phi ^{\text{bath}}_s = 1 - (1 + |z_{+}/z_{-}|) \phi ^{\text{bath}}_{+}$ . Consequently, the far-field chemical potentials are given by

(2.18a) \begin{align} \mu ^{\text{bath}}_{s} = \log \phi ^{\text{bath}}_{s}, \quad \mu ^{\text{bath}}_{\pm} = \log \phi ^{\text{bath}}_{\pm} + z_{\pm } \Phi ^{\text{bath}}. \end{align}

The radius of the deformed gel, and hence the position of the gel–bath interface, is denoted by $a$ . We use the notation $r \to a^{\pm }$ to describe approaching the interface from the bath ( $+$ ) and gel ( $-$ ). Due to the formulation of the model in terms of Eulerian coordinates, the deformed radius of the gel, $a$ , is unknown. Since the undeformed radius of the gel has been scaled to unity, the deformed radius is implicitly determined by the equation

(2.19) \begin{align} R|_{r= a^{-}} = 1. \end{align}

At the gel–bath interface, $r = a$ , thermodynamic equilibrium demands that the chemical potentials are continuous. Moreover, since the chemical potentials must also be spatially uniform, we have that

(2.20) \begin{align} \mu _m = \mu ^{\text{bath}}_m \end{align}

in both the gel and the bath. We also impose the condition

(2.21) \begin{align} \left .\frac{{\textrm{d}} \phi _s}{{\textrm{d}} r}\right |_{r = a^{-}} = 0 \end{align}

at the gel surface, which is needed due to the second derivative in (2.5a). From a physical point of view, (2.21) implies that the solvent does not preferentially wet or dewet the interface, both of which would lead to a localised gradient in the solvent composition. The balance of radial stresses at the interface leads to

(2.22) \begin{align} \left .{\mathcal{G}} \mathsf{T}_{rr} \right |_{r = a^{-}} = \left .\epsilon _r \varepsilon ^2 \mathsf{T}_{rr}\right |_{r = a^{+}}. \end{align}

We assume there are no surface charges on the interface and therefore impose continuity of the electric potential and electric displacement:

(2.23a) \begin{align} \left .\Phi \right |_{r = a^{-}} &= \left .\Phi \right |_{r = a^{+}}, \end{align}

(2.23b) \begin{align} \left .\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right |_{r= a^{-}} &= \epsilon _r\left .\frac{{\textrm{d}} \Phi }{{\textrm{d}} r}\right |_{r = a^{+}}. \end{align}

2.4 Parameter estimation

We assume that the molecular volume is $\nu \sim 10^{-28}$ m $^{3}$ [Reference Yu, Landis and Huang34], the system is held at a temperature of $T = 300$ K and dry radius of the gel is $a_0\sim 1$ cm. Horkay et al. [Reference Horkay, Tasaki and Basser14] measured the shear moduli of polyelectrolyte gels to be around $G \sim 10$ kPa, which leads to ${\mathcal{G}} \sim 10^{-4}$ . Yu et al. [Reference Yu, Landis and Huang34] reported values of ${\mathcal{G}} \sim 10^{-3}$ . The Flory interaction parameter $\chi$ is generally a function of the gel composition and temperature. However, we treat $\chi$ as a constant, which is a common simplification in the literature. Yu et al. [Reference Yu, Landis and Huang34] use constant values of $\chi$ that range from 0.1 to 1.6.

We assume that the electrical permittivity of the gel and the bath are approximately the same as water due to the ions being dilute. Thus, we set $\epsilon ^{\text{gel}} \simeq \epsilon ^{\text{bath}} \simeq 80\,\epsilon _0$ , where $\epsilon _0$ is the permittivity of free space. Hence, $\epsilon _r = \epsilon ^{\text{gel}}/ \epsilon ^{\text{bath}} \simeq 1$ . The non-dimensional width of the EDL is then $\varepsilon \sim 10^{-8}$ , corresponding to a dimensional value of $0.1$ nm. However, we will show in Sec. 5 that the value of $\varepsilon$ underestimates the width of the EDL because it is based on an imprecise estimate of the ionic volume fractions.

The dimensionless parameter $\omega$ is difficult to estimate due to uncertainties in the values of the Kuhn length. Hua et al. [Reference Hua, Mitra and Muthukumar15] set $L_K = 0.9$ nm in their modelling study. Similarly, Wu et al. [Reference Wu, Jha and de la Cruz31] take $L_K = 1$ nm. Both values lead to an estimate of $\omega \sim 10^{-7}$ .

3.1 The outer problem

We now consider the limit $\varepsilon \to 0$ with $r = O(1)$ . The outer variables are expanded as $f(r) = f^{(0)}(r) + O(\varepsilon )$ , where $f$ is an arbitrary quantity.

3.1.1 Electroneutral equations for the bath

Taking $\varepsilon \to 0$ in (2.13) leads to the electroneutrality condition for the bath

(3.1) \begin{align} z_+ \phi _+^{(0)} + z_- \phi _-^{(0)} = 0. \end{align}

The chemical potentials (2.12) reduce to

(3.2a) \begin{align} \mu ^{\text{bath}}_s &= \log \phi _s^{(0)}, \end{align}

(3.2b) \begin{align} \mu ^{\text{bath}}_{\pm } &= \log \phi _{\pm }^{(0)} + z_{\pm } \Phi ^{(0)}. \end{align}

Solving these four equations shows that the outer solution in the bath corresponds to a homogeneous mixture with the same composition and voltage as the far field:

(3.3) \begin{align} \phi ^{(0)}_{m}(r) = \phi ^{\text{bath}}_m, \quad \Phi ^{(0)}(r) = \Phi ^{\text{bath}}. \end{align}

The radial stress balance (2.14) reduces to ${\textrm{d}} p^{(0)}/{\textrm{d}} \xi = 0$ . Solving and matching to the far field implies that

(3.4) \begin{align} p^{(0)}(r) = 0, \quad \mathsf{T}_{rr}^{(0)}(r) = 0. \end{align}

Thus, the bath is stress free to leading order.

3.1.2 Electroneutral equations for the gel

Motivated by the constant outer solutions for the bath, as well as past works in the literature [Reference Drozdov and deClaville Christiansen9, Reference Drozdov, deClaville Christiansen and Sanporean10, Reference Horkay, Tasaki and Basser14, Reference Hua, Mitra and Muthukumar15, Reference Ohmine and Tanaka23, Reference Yu, Landis and Huang34, Reference Zhang, Dehghany and Hu35], we assume that the outer solution for the gel represents a homogeneous state. We therefore write the leading-order contributions to the outer solution as $f^{(0)}(r) = f^{\text{gel}}$ , where $f^{\text{gel}}$ is a constant, for all variables except the Lagrangian radius $R^{(0)}$ , which retains a dependence on $r$ .

The $O(1)$ contributions to (2.7) lead to the electroneutrality condition for the gel,

(3.5) \begin{align} z_+ \phi ^{\text{gel}}_+ + z_- \phi ^{\text{gel}}_- = - z_f \varphi _f/ J^{\text{gel}}. \end{align}

The solvent and ionic chemical potentials (2.5) are given by

(3.6a) \begin{align} \Pi ^{\text{gel}}_s +{\mathcal{G}} p^{\text{gel}} = \mu ^{\text{bath}}_s, \end{align}

(3.6b) \begin{align} \Pi ^{\text{gel}}_{\pm } +{\mathcal{G}} p^{\text{gel}} + z_{\pm } \Phi ^{\text{gel}} = \mu ^{\text{bath}}_{\pm }. \end{align}

For a homogeneous gel that is free to swell in the radial and orthoradial directions, the radial and orthoradial stretches are equal; thus, $\lambda ^{\text{gel}}_r = \lambda ^{\text{gel}}_{\theta } = (J^{\text{gel}}/\lambda _z)^{1/2}$ . The leading-order contribution to the Lagrangian coordinate can then be obtained from (2.4) as $R^{(0)}(r) = (\lambda _z/ J^{\text{gel}})^{1/2} r$ . We define

(3.7) \begin{align} R^{\text{gel}} \equiv R^{(0)}(a^{-}) = (\lambda _z/ J^{\text{gel}})^{1/2} a, \end{align}

where $R^{\text{gel}}$ must be determined by matching to the inner solution. The final quantity to determine is the mechanical pressure in the gel. The radial stress balance (2.8) implies that the radial Cauchy stress $\mathsf{T}_{rr}^{\text{gel}}$ must be a constant. Using asymptotic matching, we will show that $\mathsf{T}_{rr}^{\text{gel}} = 0$ . Hence, from (2.9a), we have that $p^{\text{gel}} = 1/ \lambda _z - 1/ J^{\text{gel}}$ .

3.2 The inner problem

The inner problem is formulated by introducing the change of variable $r = a + \varepsilon \xi$ , where $\xi$ is a radial coordinate that is localised to the free surface of the gel. By definition, $\xi \gt 0$ corresponds to the regions in the bath whereas $\xi \lt 0$ corresponds to regions in the gel. Tildes are used to denote dependent variables in the inner region, which are generally expanded as $\tilde{f} = \tilde{f}^{(0)} + \varepsilon \tilde{f}^{(1)} + O(\varepsilon ^2)$ , where $\tilde{f}$ is an arbitrary quantity. However, additional rescaling is required in some cases; this will be made explicit in the proceeding discussion. Near the interface, the outer solutions for the bath and gel are expanded in terms of inner variables, respectively, as

(3.8a) \begin{align} &f(a + \varepsilon \xi ) = f^{(0)}(a^{+}) + O(\varepsilon ) = f^{\text{bath}} + O(\varepsilon ), \end{align}

(3.8b) \begin{align} &f(a + \varepsilon \xi ) = f^{(0)}(a^{-}) + O(\varepsilon ) = f^{\text{gel}} + O(\varepsilon ), \end{align}

which will be used for asymptotic matching.

3.2.1 Inner problem for the bath

3.2.1.1 Mechanics

After changing variables, we anticipate that the Maxwell stresses (2.16) scale as $\mathsf{T}_M = O(\varepsilon ^{-2})$ because the electric potential $\Phi$ should remain $O(1)$ in size across the EDL. Moreover, we expect that the pressure will scale like the Maxwell stress so that $p = O(\varepsilon ^{-2})$ ; the rationale behind this choice is discussed in the introduction of Sec. 2. Therefore, we write $\mathsf{T}_M = \varepsilon ^{-2} \tilde{\mathsf{T}}_M$ and $p = \varepsilon ^{-2} \tilde{p}$ . Consequently, the Cauchy stresses must also be scaled as $\mathsf{T} = \varepsilon ^{-2} \tilde{\mathsf{T}}$ . Expanding $\tilde{\mathsf{T}}_{rr}$ and $\tilde{\mathsf{T}}_{M,rr}$ in powers of $\varepsilon$ and matching to the far field leads to the stress-free condition $\tilde{\mathsf{T}}_{rr}^{(0)}\to 0$ as $\xi \to \infty$ . The leading-order contribution to the radial stress balance in the bath (2.14) leads to ${\textrm{d}} \tilde{\mathsf{T}}^{(0)}_{rr}/{\textrm{d}} \xi = 0$ . Integrating and imposing the far-field condition leads to the conclusion that $\tilde{\mathsf{T}}_{rr}^{(0)} = 0$ . In terms of the original scaling, this means that, in the EDL, the total stress in the bath must be $O(\varepsilon ^{-1})$ in size. The pressure can be obtained from (2.15a) and is found to be

(3.9) \begin{align} \tilde{p}^{(0)} = \frac{1}{2}\left (\frac{\partial \tilde{\Phi }^{(0)}}{\partial \xi }\right )^2 = \tilde{\mathsf{T}}_{M,rr}^{(0)}. \end{align}

Thus, as expected, the pressure in the bath balances the Maxwell stresses.

3.2.1.2 Chemical equilibrium

Expanding the chemical potentials (2.12) gives, to leading order,

(3.10a) \begin{align} \log \tilde{\phi }_s^{(0)} + \epsilon _r \tilde{p}^{(0)} &= \mu ^{\text{bath}}_s, \end{align}

(3.10b) \begin{align} \log \tilde{\phi }_{\pm }^{(0)} + \epsilon _r \tilde{p}^{(0)} + z_{\pm } \tilde{\Phi }^{(0)} &= \mu ^{\text{bath}}_{\pm }. \end{align}

Equating (3.10) with (2.18a) and using (3.9) to eliminate the pressure leads to an expression for the ion fractions in the EDL,

(3.11) \begin{align} \tilde{\phi }_{\pm }^{(0)} = \phi ^{\text{bath}}_{\pm } \exp \left [z_{\pm }(\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)}) - \frac{\epsilon _r}{2}\left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2\right ]. \end{align}

This is a modification of the Boltzmann distribution for the ions, which arises from accounting for the pressure dependence of the ionic chemical potentials.

3.2.1.3 Electrostatics

The leading-order electrical problem is obtained by combining (2.13) with the ionic volume fractions (3.11) to obtain a modified Poisson–Boltzmann equation given by

(3.12) \begin{align} -\epsilon _r \frac{{\textrm{d}}^2 \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi ^2} = \exp \left [-\frac{\epsilon _r}{2} \left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2\right ]\sum _{i\in \mathbb{I}}z_i \phi ^{\text{bath}}_i\exp \left (z_{i}(\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)})\right ). \end{align}

Equation (3.12) can be integrated once and the conditions ${\textrm{d}} \tilde{\Phi }^{(0)}/{\textrm{d}} \xi \to 0$ and $\tilde{\Phi }^{(0)} \to \Phi ^{\text{bath}}$ as $\xi \to \infty$ used to obtain

(3.13) \begin{align} \frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi } = \mp \sqrt{\frac{2}{\epsilon _r}\log \left \{ 1 + \sum _{i\in \mathbb{I}}\phi ^{\text{bath}}_i\left [\exp \left (z_i(\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)})\right ) -1 \right ]\right \}}. \end{align}

The minus sign is taken if $\Phi ^{\text{gel}} - \Phi ^{\text{bath}} \gt 0$ , which will generally be the case if the fixed charges on the polymer chains are positive, as assumed here.

3.2.2 Inner problem for the gel

3.2.1.4 Chemical equilibrium

The chemical potential of solvent (2.5a) can be expanded as

(3.14) \begin{align} \mu ^{\text{bath}}_s = \tilde{\Pi }_s^{(0)} +{\mathcal{G}} \tilde{p}^{(0)} - \varepsilon ^{-2}\omega ^2\Bigg [&\frac{{\textrm{d}}^2 }{{\textrm{d}} \xi ^2}\left ( \tilde{\phi }_s^{(0)} + \varepsilon \tilde{\phi }_s^{(1)} + \varepsilon ^2 \tilde{\phi }_s^{(2)}\right ) + \frac{\varepsilon }{a}\frac{{\textrm{d}} }{{\textrm{d}} \xi }\left (\tilde{\phi }_s^{(0)} + \varepsilon \tilde{\phi }_s^{(1)}\right ) - \varepsilon ^2 \frac{\xi }{a^2} \frac{{\textrm{d}} \tilde{\phi }_s^{(0)}}{{\textrm{d}} \xi }\Bigg ] + O(\varepsilon ). \end{align}

Similarly, the boundary condition at the gel–bath interface (2.21) can be expanded to give ${\textrm{d}} \tilde{\phi }_s^{(n)}/{\textrm{d}} \xi = 0$ at $\xi = 0$ for $n = 0, 1, 2$ . The $O(\varepsilon ^{-2})$ and $O(\varepsilon ^{-1})$ contributions to (3.14) along with the boundary and matching conditions show that the solvent concentration is uniform to leading and next order,

(3.15) \begin{align} \tilde{\phi }_s^{(0)}(\xi ) = \phi _s^{\text{gel}}, \quad \tilde{\phi }_s^{(1)}(\xi ) = \phi _s^{(1)}, \end{align}

which is a distinguishing feature of the asymptotic limit in which $\varepsilon \to 0$ with $\omega$ fixed. Physically, this result is a consequence of gradients in the solvent concentration having a high energy cost when the Kuhn length is large. The ionic chemical potentials (2.5b) can be expanded as

(3.16) \begin{align} \mu ^{\text{bath}}_{\pm } = \log \tilde{\phi }_{\pm }^{(0)} + \frac{1}{\tilde{J}^{(0)}}(1 - \chi \phi ^{\text{gel}}_s) +{\mathcal{G}} \tilde{p}^{(0)} + z_{\pm } \tilde{\Phi }^{(0)}. \end{align}

By combining (3.16) and (2.18a) and using (3.15), we find that

(3.17) \begin{align} \tilde{\phi }_{\pm }^{(0)} = \phi ^{\text{bath}}_{\pm }\exp \left [z_\pm (\Phi ^{\text{bath}}-\tilde{\Phi }^{(0)}) -{\mathcal{G}}\tilde{p}^{(0)} - \frac{1}{\tilde{J}^{(0)}}(1 - \chi \phi ^{\text{gel}}_s)\right ]. \end{align}

Although (3.17) appears to be a closed-form expression for the volume fraction of ions, it is important to recall that the swelling ratio $\tilde{J}^{(0)}$ also depends on the these quantities; see (2.3).

3.2.1.5 Kinematics and incompressibility

The $O(\varepsilon ^{-1})$ part of the incompressibility condition (2.4) implies that ${\textrm{d}} \tilde{R}^{(0)}/{\textrm{d}} \xi = 0$ . Imposing the boundary condition $\tilde{R}^{(0)}(0) = 1$ leads to $\tilde{R}^{(0)} = 1$ . Matching the inner and outer solutions leads to the condition $R^{\text{gel}} = 1$ , which can be used in combination with (3.7) to find that the radius of the deformed gel is given by

(3.18) \begin{align} a = (J^{\text{gel}}/ \lambda _z)^{1/2}. \end{align}

The $O(1)$ part of (2.4) gives ${\textrm{d}} \tilde{R}^{(1)}/{\textrm{d}} \xi = \lambda _z a/ \tilde{J}^{(0)}$ . The leading-order radial stretch can then be found by expanding (2.2) to find that $\tilde{\lambda }_r^{(0)} = (\lambda _z J^{\text{gel}})^{-1/2}\tilde{J}^{(0)}$ . Similarly, the leading-order orthoradial stretch is given by $\tilde{\lambda }_{\theta }^{(0)} = (J^{\text{gel}}/\lambda _z)^{1/2}$ .

3.2.1.6 Mechanics

The leading-order part of the stress balance in the gel (2.8) is ${\textrm{d}} \tilde{\mathsf{T}}^{(0)}_{rr}/{\textrm{d}} \xi = 0$ . Hence, the total stress in the gel is constant across the EDL. Imposing stress continuity at the interface (2.22) and using the fact that the stress in the bath is $O(\varepsilon ^{-1})$ shows that the gel is stress free to leading order, $\tilde{\mathsf{T}}_{rr}^{(0)} = 0$ . Matching to the outer solution leads to $\mathsf{T}_{rr}^{\text{gel}} = 0$ , as previously claimed.

The mechanical pressure can be obtained from (2.9a) as

(3.19) \begin{align} \tilde{p}^{(0)} = \tilde{\mathsf{T}}_{e,rr}^{(0)} + \tilde{\mathsf{T}}_{K,rr}^{(0)} + \tilde{\mathsf{T}}_{M,rr}^{(0)}. \end{align}

The leading-order radial components of the elastic, Korteweg, and Maxwell stresses can be obtained from (2.10a), (2.11a), and (2.11c) as

(3.20a) \begin{align} \tilde{\mathsf{T}}_{e,rr}^{(0)} &= \frac{1}{\lambda _z} \frac{\tilde{J}^{(0)}}{J^{\text{gel}}} - \frac{1}{\tilde{J}^{(0)}}, \end{align}

(3.20b) \begin{align} \tilde{\mathsf{T}}_{K,rr}^{(0)} &= \frac{\omega ^2}{{\mathcal{G}}} \phi ^{\text{gel}}_s\frac{{\textrm{d}}^2 \tilde{\phi }_s^{(2)}}{{\textrm{d}} \xi ^2}, \end{align}

(3.20c) \begin{align} \tilde{T}_{M,rr}^{(0)} &= \frac{1}{2{\mathcal{G}}}\left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2. \end{align}

In order to evaluate the Korteweg stresses without explicitly solving for $\tilde{\phi }_s^{(2)}$ , the $O(1)$ contributions to the solvent chemical potential (3.14) can be used in (3.20b) to obtain

(3.21) \begin{align} \tilde{\mathsf{T}}_{K,rr}^{(0)} = \frac{\phi ^{\text{gel}}_s}{{\mathcal{G}}}\left (\tilde{\Pi }_s^{(0)} +{\mathcal{G}} \tilde{p}^{(0)} - \mu ^{\text{bath}}_s\right ). \end{align}

Note that setting $\omega = 0$ leads to $\tilde{\Pi }_s^{(0)} +{\mathcal{G}} \tilde{p}^{(0)} - \mu ^{\text{bath}}_s = 0$ from (3.14) and hence $\tilde{\mathsf{T}}_{K,rr}^{(0)} = 0$ . Substitution of (3.21) into (3.19) gives an algebraic relation for the pressure $\tilde{p}^{(0)}$ .

3.2.1.7 Electrostatics

The leading-order electrical problem in the gel is given by

(3.22) \begin{align} -\frac{{\textrm{d}}^2 \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi ^2} = z_+ \tilde{\phi }_{+}^{(0)} + z_-\tilde{\phi }_-^{(0)} + \frac{z_f\varphi _f}{\tilde{J}^{(0)}}, \end{align}

which is coupled to the algebraic equations for the volume fractions of ions (3.17) and the mechanical pressure (3.19). The electrical problems for the bath and gel can be decoupled by combining the first integral for the electric potential in the bath (3.13) with the electrostatic boundary conditions (2.23) to obtain

(3.23a) \begin{align} \left .\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right |_{\xi = 0^{-}} = \left .\mp \sqrt{2\epsilon _r \log \left \{ 1 + \sum _{i\in \mathbb{I}}\phi ^{\text{bath}}_i\left [\exp \left (z_i(\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)})\right ) -1 \right ]\right \}}\right |_{\xi = 0^{-}}, \end{align}

which acts as a boundary condition for (3.22). The electrical problem in the gel is closed by imposing the matching condition

(3.23b) \begin{align} \tilde{\Phi }^{(0)} \to \Phi ^{\text{gel}}, \quad \xi \to -\infty. \end{align}

3.3 Summary

The asymptotic analysis has produced a closed system of algebraic equations that determines the outer solution in the gel. The inner problems for the gel and the bath decouple. In the case of the gel, the inner problem can be condensed into a nonlinear system of differential-algebraic equations; this system will be presented in Sec. 5. The inner problem for the bath has been solved in terms of the electric potential; this can be obtained by integrating (3.13) once the electric potential at the gel surface has been determined by solving the inner problem for the gel.

5.1 Solution of the outer problem for the gel

In the outer region of the gel, the volume fraction of solvent and ions, as well as the electric potential, are determined from (3.5) to (3.6). This nonlinear algebraic system can be formulated as

(5.1a) \begin{align} \log \phi ^{\text{gel}}_s + \frac{1}{J^{\text{gel}}} + \frac{\chi (1 - \phi ^{\text{gel}}_s)}{J^{\text{gel}}} +{\mathcal{G}}\left (\frac{1}{\lambda _z} - \frac{1}{J^{\text{gel}}}\right ) = \log (1 - 2 \phi ^{\text{bath}}_+), \end{align}

(5.1b) \begin{align} \phi ^{\text{gel}}_{\pm } = \phi ^{\text{bath}}_+ \exp \left [\pm (\Phi ^{\text{bath}} - \Phi ^{\text{gel}}) -{\mathcal{G}}\left (\frac{1}{\lambda _z}-\frac{1}{J^{\text{gel}}}\right ) - \frac{1}{J^{\text{gel}}}\left (1 - \chi \phi ^{\text{gel}}_s\right )\right ], \end{align}

(5.1c) \begin{align} 2 \phi ^{\text{bath}}_+ \sinh (\Phi ^{\text{bath}} - \Phi ^{\text{gel}}) = -\frac{z_f \varphi _f}{J^{\text{gel}}}\exp \left [{\mathcal{G}}\left (\frac{1}{\lambda _z}-\frac{1}{J^{\text{gel}}}\right ) + \frac{1}{J^{\text{gel}}}\left (1 - \chi \phi ^{\text{gel}}_s\right )\right ], \end{align}

where $J^{\text{gel}}$ is given by (2.3). When the cation fraction in the bath is small, $\phi ^{\text{bath}}_{+} \ll 1$ , the nonlinear system (5.1) can be reduced to a single equation, as described in Appendix B.

We numerically solve (5.1) over a range of values of $\phi ^{\text{bath}}_{+}$ using pseudo-arclength continuation. Three values of $\lambda _z \leq 1$ are considered, corresponding to gels in axial compression. The equilibrium swelling ratios $J^{\text{gel}}$ computed from the solutions of (5.1) are shown as solid curves in Figure 2. The dashed black lines represent numerical solutions to the reduced model derived in Appendix B. The figure shows that for each value of $\lambda _z$ there are three branches of solutions, two of which intersect and then vanish as the salt fraction in the bath $\phi ^{\text{bath}}_{+}$ increases beyond a critical value. The loss of equilibrium solutions at this critical point indicates that a volume phase transition can occur, as the gel volume will undergo a discontinuous decrease as the salt content of the bath is increased. For a given value of $\lambda _z$ , the branch of solutions corresponding to the largest and smallest values of $J^{\text{gel}}$ are referred to as the swollen and collapsed branches, respectively. The branch of solutions corresponding to intermediate values of $J^{\text{gel}}$ is unstable [Reference Celora, Hennessy, Münch, Wagner and Waters4] and will not be considered further.

Figure 2. (a) Equilibrium swelling ratio $J^{\text{gel}}$ as a function of cation fraction in the bath $\phi ^{\text{bath}}_+$ showing swollen and collapsed branches. (b) The swelling ratio along the collapsed branch. Solid lines correspond to solutions of (5.1). Dashed lines represent solutions to the reduced equation (B.4) for a dilute concentration of cations. The parameter values are ${\mathcal{G}} = 5 \times 10^{-4}$ , $\chi = 1.2$ , $\varphi _f = 0.05$ $z_\pm = \pm 1$ , $z_f = 1$ .

For a fixed value of the salt fraction, increasing the axial compression reduces the degree of swelling that occurs. Moreover, increasing the axial compression also decreases the critical salt fraction at which the volume phase transition occurs. Both of these findings are in agreement with experimental observations [Reference Horkay, Tasaki and Basser14]. Due to the incompressibility of the gel, imposing an axial compression results in a radial stretch. The elastic energy cost of inserting a molecule into a pre-stretched gel is greater than for an unstretched gel. Hence, the balance between the mixing and elastic energies is established at smaller concentrations, resulting in the equilibrium swelling ratio $J^{\text{gel}}$ decreasing with the axial stretch $\lambda _z$ .

5.2 Formulation of the inner problems

The inner problem for the gel can now be constructed using the results from the previous sections. In particular, if $\omega = 0$ , then the governing equations for the gel can be condensed into

(5.2a) \begin{align} \tilde{\Pi }_s^{(0)}+{\mathcal{G}}\tilde{p}^{(0)} = \log (1 - 2 \phi ^{\text{bath}}_+), \end{align}

(5.2b) \begin{align} \tilde{\phi }^{(0)}_{\pm } = \phi ^{\text{bath}}_+ \exp \left [\pm (\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)}) -{\mathcal{G}}\tilde{p}^{(0)} - \frac{1}{\tilde{J}^{(0)}}\left (1 - \chi \tilde{\phi }^{(0)}_s\right )\right ], \end{align}

(5.2c) \begin{align} -\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi } = \tilde{\phi }_{+}^{(0)} - \tilde{\phi }_-^{(0)} + \frac{z_f \varphi _f}{\tilde{J}^{(0)}}, \end{align}

(5.2d) \begin{align} \tilde{p}^{(0)} = \frac{1}{\lambda _z} \frac{\tilde{J}^{(0)}}{J^{\text{gel}}} - \frac{1}{\tilde{J}^{(0)}} + \frac{1}{2{\mathcal{G}}}\left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2, \end{align}

(5.2e) \begin{align} \tilde{J}^{(0)} = (1 - \tilde{\phi }^{(0)}_s - \tilde{\phi }_{+}^{(0)} - \tilde{\phi }_{-}^{(0)})^{-1}, \end{align}

(5.2f) \begin{align} \tilde{\Pi }_s^{(0)} = \log \tilde{\phi }^{(0)}_s + \frac{\chi (1 - \tilde{\phi }_s^{(0)})}{\tilde{J}^{(0)}} + \frac{1}{\tilde{J}^{(0)}}. \end{align}

In the case $\omega \gg \varepsilon$ , equation (5.2a) is replaced with $\tilde{\phi }_s^{(0)} = \phi ^{\text{gel}}_s$ , resulting in the system

(5.3a) \begin{align} \tilde{\phi }^{(0)}_{\pm } = \phi ^{\text{bath}}_+ \exp \left [\pm (\Phi ^{\text{bath}} - \tilde{\Phi }^{(0)}) -{\mathcal{G}}\tilde{p}^{(0)} - \frac{1}{\tilde{J}^{(0)}}\left (1 - \chi \phi ^{\text{gel}}_s\right )\right ], \end{align}

(5.3b) \begin{align} -\frac{{\textrm{d}}^2 \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi ^2} = \tilde{\phi }_{+}^{(0)} - \tilde{\phi }_-^{(0)} + \frac{z_f \varphi _f}{\tilde{J}^{(0)}}, \end{align}

(5.3c) \begin{align}{\mathcal{G}} (1 - \phi ^{\text{gel}}_s) \tilde{p}^{(0)} ={\mathcal{G}}\left (\frac{1}{\lambda _z} \frac{\tilde{J}^{(0)}}{J^{\text{gel}}} - \frac{1}{\tilde{J}^{(0)}}\right ) + \phi ^{\text{gel}}_s \left (\tilde{\Pi }_s^{(0)} - \mu ^{\text{bath}}_s\right ) + \frac{1}{2}\left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2, \end{align}

(5.3d) \begin{align} \tilde{J}^{(0)} = (1 - \phi ^{\text{gel}}_s - \tilde{\phi }_{+}^{(0)} - \tilde{\phi }_{-}^{(0)})^{-1}, \end{align}

(5.3e) \begin{align} \tilde{\Pi }_s^{(0)} = \log \phi ^{\text{gel}}_s + \frac{\chi (1 - \phi ^{\text{gel}}_s)}{\tilde{J}^{(0)}} + \frac{1}{\tilde{J}^{(0)}}. \end{align}

In both cases, the boundary conditions for the electrical potential are given by (3.23a). The expression for the hoop stress in the gel is the same in both cases as well:

(5.4) \begin{align} \tilde{\mathsf{T}}^{(0)}_{\theta \theta } = \frac{1}{\lambda _z}\left (\frac{J^{\text{gel}}}{\tilde{J}^{(0)}} - \frac{\tilde{J}^{(0)}}{J^{\text{gel}}}\right ) -{\mathcal{G}}^{-1}\left (\frac{{\textrm{d}} \tilde{\Phi }^{(0)}}{{\textrm{d}} \xi }\right )^2. \end{align}

The first term represents the elastic contribution to the total hoop stress, which can be compressive or tensile. The second term captures the contribution from the Maxwell stresses, which is always compressive.

Figure 3. Numerical solutions of the inner problems with far-field conditions corresponding to the collapsed state. The regions of gel and bath are defined by $\xi \lt 0$ and $\xi \gt 0$ , respectively. Dashed lines represent the solution to (5.2) when $\omega = 0$ . Solid lines represent the solution to (5.3) when $\omega \gg \varepsilon$ . The parameter values are $\chi = 1.2$ , ${\mathcal{G}} = 5 \times 10^{-4}$ , $\varphi _f = 0.05$ , $\phi ^{\text{bath}}_+ = 10^{-5}$ , $\lambda _z = 1$ , $\epsilon _r = 1$ , $z_{\pm } = \pm 1$ , and $z_f = 1$ .

Figure 4. Numerical solution of the inner problem with far-field conditions corresponding to the swollen state when $\omega \gg \varepsilon$ . The regions of gel and bath are defined by $\xi \lt 0$ and $\xi \gt 0$ , respectively. The curves are obtained by solving (5.2). Parameter values are the same as in Figure 3: $\chi = 1.2$ , ${\mathcal{G}} = 5 \times 10^{-4}$ , $\varphi _f = 0.05$ , $\phi ^{\text{bath}}_+ = 10^{-5}$ , $\lambda _z = 1$ , $\epsilon _r = 1$ , $z_{\pm } = \pm 1$ , and $z_f = 1$ .

Figure 5. Phase separation in the inner region. (a)–(c) The swelling ratio and (d)–(f) total electric charge density for different salt fractions in the bath $\phi ^{\text{bath}}_{+}$ . The curves are obtained by numerically solving the inner problem in the intermediate asymptotic limit $\omega = O(\varepsilon )$ as $\varepsilon \to 0$ ; see (5.5). We have taken $\omega = \Omega \varepsilon$ with $\Omega = 10^{-1}$ . The remaining parameters are $\chi = 0.7$ , ${\mathcal{G}} = 4\times 10^{-3}$ , $\varphi _f = 0.04$ , $z_{\pm } = \pm 1$ , $z_f = 1$ , $\epsilon _r = 1$ , and $\lambda _z = 1$ .

5.3 The structure of the electric double layer

The systems (5.2) and (5.3) are discretised using finite differences and solved using Newton’s method. The asymptotic solution is validated against numerical solutions of the full problem in Appendix C.

We consider the case where the axial stretch and salt content in the bath are set to $\lambda _z = 1$ and $\phi ^{\text{bath}}_{+} = 10^{-5}$ , with the remaining parameters being the same as those in Figure 2. There are three possible solutions to the outer problem. We are only concerned with two of these, which correspond to the collapsed state ( $J^{\text{gel}} \simeq 1.447$ ) and the swollen state ( $J^{\text{gel}} \simeq 82$ ).

In Figure 3, we plot the inner solution when the outer solution corresponds to the collapsed state ( $J^{\text{gel}} \simeq 1.447$ ). The solid and dashed lines correspond to the cases $\omega \gg \varepsilon$ and $\omega = 0$ , respectively. In both cases, we see that our non-dimensionalisation underestimates the width of the double layer, which is about $10 \varepsilon$ in the gel (or 1 nm) and $1000 \varepsilon$ in the bath (or 100 nm). For this parameter set, the value of $\omega$ does not lead to noticeable changes in the electric potential and ion fractions; see Figure 3 (a)–(b). However, substantial differences arise in the gel pressure and the solvent fraction; see Figure 3 (c)–(d). When $\omega = 0$ , the gel pressure balances a large Maxwell stress. This large pressure causes a local decrease in the solvent fraction and a slight volumetric contraction of the gel (Figure 3 (d)), which can be rationalised in terms of equation (4.1a). At equilibrium, the osmotic pressure $\tilde{\Pi }_s$ must balance the mechanical pressure $\tilde{p}$ . To compensate for the increase in mechanical pressure that arises from the Maxwell stress, the osmotic pressure must decrease, which drives solvent out of the gel and causes it to shrink. When $\omega \gg \varepsilon$ , gradients in the solvent fraction are energetically penalised; thus, the solvent fraction remains uniform across the EDL. From a mechanical perspective, this penalisation occurs through the development of a large Korteweg stress, which counters the effect of the Maxwell stress in order to maintain a uniform solvent fraction. The mechanical contribution from the Korteweg stress manifests as an increase in the gel pressure compared to the $\omega = 0$ case, as seen in Figure 3 (c). Although the solvent fraction is constant across the EDL when $\omega \gg \varepsilon$ , the swelling ratio still decreases relative to the bulk value (Figure 3 (d)) due to the variation in ionic content (Figure 3 (b)).

The inset of Figure 3 (c) shows the total hoop stress in the gel, which is the same in both models owing to the strong similarities in the electric potential. Due to the large Maxwell stress, the gel experiences a substantial compressive hoop stress, which leads to the intriguing possibility of mechanical instabilities in the EDL.

In Figure 4, we show the numerical solution of the inner problem with $\omega \gg \varepsilon$ when the outer solution corresponds to the swollen state ( $J^{\text{gel}} \simeq 82$ ). The qualitative features of the inner solution are similar to those obtained when the outer solution corresponds to the collapsed state (Figure 3). However, an important difference is that the volume fraction of anions has decreased by more than a factor of ten. This decrease is driven by the reduction in the volume fraction of fixed charges on the polymers that occurs when the gel is highly swollen; see Figure 4 (b). Consequently, the EDL in the gel has increased in thickness by a roughly factor of ten to approximately $250 \varepsilon$ (or 25 nm). The gradient in the electric potential in the gel is therefore ten times weaker, resulting in a 100-fold reduction in the Maxwell stress and the total hoop stress; see Figure 4 (c). Despite these decreases, the pressure in the gel remains large because of the Korteweg stress. Due to convergence issues, it was not possible to compute the corresponding inner solution when $\omega = 0$ .

To understand the origin of these numerical difficulties, we consider an intermediate asymptotic limit where $\omega = \Omega \varepsilon$ , with $\Omega = O(1)$ as $\varepsilon \to 0$ . By assuming that the gel remains in a homogeneous, swollen state away from the EDL, the outer problem in this limit is still given by (5.1). The corresponding inner problem can be formulated by changing (3.14) or (4.1a) to

(5.5a) \begin{align} \tilde{\Pi }_s^{(0)} +{\mathcal{G}} \tilde{p}^{(0)} + \Omega ^2 \frac{{\textrm{d}}^2 \tilde{\phi }_s^{(0)}}{{\textrm{d}} \xi ^2} = \mu ^{\text{bath}}_s. \end{align}

The pressure (3.19) can be evaluated using a Korteweg stress given by

(5.5b) \begin{align} \tilde{\mathsf{T}}_{K,rr}^{(0)} ={\mathcal{G}}^{-1} \Omega ^2 \left [\tilde{\phi }_s^{(0)}\frac{\textrm{d} ^2 \tilde{\phi }_s^{(0)}}{\textrm{d} \xi ^2} - \frac{1}{2}\left (\frac{\textrm{d} \tilde{\phi }_s^{(0)}}{\textrm{d} \xi }\right )^2\right ]. \end{align}

We solve the intermediate asymptotic model based on (5.5) by imposing the boundary conditions ${\textrm{d}} \tilde{\phi }_s^{(0)}/{\textrm{d}} \xi = 0$ as $\xi \to -\infty$ and $\xi \to 0^{-}$ . A second parameter set is used to reduce the degree of swelling that occurs in the gel. This parameter set leads to the outer problem having single branch of equilibrium solutions that does not fold back on itself. Thus, the gel monotonically and continuously decreases in volume as the salt fraction in the bath $\phi ^{\text{bath}}_+$ increases.

Figure 6. Phase separation drives the breakdown of charge neutrality in the gel when the Debye length is comparable to the Kuhn length. (a) The swelling ratio and (b) the total electric charge density computed from the full steady problem in cylindrical coordinates. The gel self-organises into a highly swollen, negatively charged core (dark blue); a moderately swollen interior with alternating electric charge (blue), and a weakly swollen, positively charged shell (light blue). We set $\varepsilon = 10^{-2}$ , $\omega = 10^{-3}$ . The remaining parameter values correspond to Figure 5 (b) and (e) and are $\phi ^{\text{bath}}_+ = 6.6 \times 10^{-4}$ , $\chi = 0.7$ , ${\mathcal{G}} = 4\times 10^{-3}$ , $\varphi _f = 0.04$ , $z_{\pm } = \pm 1$ , $z_f = 1$ , $\epsilon _r = 1$ , and $\lambda _z = 1$ .

It is important to point out that the intermediate asymptotic model based on (5.5) is only fourth order in space. Therefore, it can be solved without explicitly imposing that its solutions tend to the homogeneous and electrically neutral outer solutions determined by (5.1). However, if the inner solutions tend to constants in the far field, then these constants must satisfy (5.1) and hence a match with the outer solution will be obtained.

The inner problem in the intermediate asymptotic limit is solved at three specific values of $\phi ^{\text{bath}}_{+}$ with $\Omega = 0.1$ . The swelling ratio $\tilde{J}^{(0)}$ and total charge density $\tilde{Q}^{(0)} = \tilde{\phi }^{(0)}_{+} - \tilde{\phi }_{-}^{(0)} + z_f \varphi _f/ \tilde{J}^{(0)}$ are computed and plotted in Figure 5. In this case, decreasing the salt fraction in the bath from $\phi ^{\text{bath}}_+ = 10^{-3}$ triggers the onset of phase separation, which gives rise to a periodic array of electrically charged phases that spans the entire inner region. Charge neutrality is not recovered in the far field, even if the domain used to numerically solve the inner problem is increased. Moreover, enforcing the boundary condition $\tilde{\phi }_s^{(0)} \to \phi ^{\text{gel}}_s$ as $\xi \to -\infty$ leads to convergence issues. Hence, the inner solution cannot be matched to the homogeneous and electroneutral outer solutions computed from (5.1). We therefore posit that homogeneous outer solutions do not always exist in the thin-EDL limit $\varepsilon \to 0$ if $\omega = O(\varepsilon )$ or $\omega = 0$ .

To explore the hypothesis that the bulk of the gel may not be homogeneous and electrically neutral at equilibrium, we solved the full steady problem in the regime $\omega = O(\varepsilon )$ by setting $\varepsilon = 10^{-2}$ and $\omega = 10^{-3}$ . The salt fraction in the bath was set to $\phi ^{\text{bath}}_+ = 6.6\times 10^{-4}$ , corresponding to the parameters in Figure 5 (b) and (e). The swelling ratio $J$ and the total charge density $Q$ , which are shown in Figure 6, reveal that phase separation occurs throughout the entire gel and gives rise to a periodic arrangement of electrically charged domains. Using numerical integration, we find that the total amount of electric charge contained within a pair of adjacent domains is on the order of $10^{-7}$ . Thus, the gel effectively separates into three distinct regions consisting of an electrically negative, highly swollen core ( $0 \lt r \lt 0.073$ ); a moderately swollen interior that is electrically neutral on average ( $0.073 \lt r \lt 2.0$ ); and a positively charged, collapsed shell ( $2.0 \lt r \lt 2.1$ ). Overall, the gel carries a net positive charge which exactly balances the net negative charge in the bath to ensure that charge neutrality holds on a global scale. The pointwise breakdown of charge neutrality across the gel indicates that it is not always appropriate to decompose the problem into inner and outer regions that are characterised by the local charge density of the gel.