Feedback-controlled solute transport through chemo-responsive polymer membranes

The precise and selective control of molecular transport through membranes is of fundamental importance for various applications in industry and medicine, such as water purification,^1–3 food processing,^4,5 nanocatalysis,^6–9 drug delivery,^10–13 and tissue engineering.^14,15 Modern membrane technology becomes increasingly inspired by responsive biomembranes with nonlinear potential-, pressure-, or flux-gated permeabilities, a bistable behavior, and memristive properties.^16–23 Such features allow the design of highly selective membrane devices that efficiently control the molecular transport, autonomously regulate the chemical milieu, and may act as logical operators, artificial synapses, or analogous filters for electrical or chemical signals. Moreover, the possible memristive properties create the foundation for information storage, adaptive responses to stimuli based upon past events, and neuromorphic systems.^24–26 

In general, such self-regulation requires a feedback mechanism controlling the transport properties in a nonlinear fashion.^27–30 In the scope of membrane applications, this may arise from various system-dependent effects, such as autocatalysis, substrate, or product inhibition;^31,32 the interplay of voltage and hydrodynamic pressure;^33,34 or, as highlighted in this work, the reciprocal impacts of molecular fluxes and membrane permeability.^35–42 In this regard, many polymeric compounds offer great potential as they are versatile in their response to various physico-chemical stimuli and environmental conditions, such as temperature, electric field, and solvent quality.^43–45 For example, the polymer can respond with a volume phase transition, either from a swollen to a collapsed state, or vice versa, in which the polymer volume fraction, ϕ, may change by orders of magnitude.^36,46–51 Such a drastic change in the polymer’s physical features, in turn, has substantial, nonlinear effects on the solute permeability of the membrane device.

Very illustrative examples are the so-called smart gating membranes,^52–59 which are (rather solid) porous membranes with polymer-coated channels that can reversibly open and close, triggered by external stimuli or, through autonomous feedback, by molecular recognition. Moreover, the literature on the solution–diffusion model^60–63 suggests that the use of more flexible, responsive polymeric membranes enables feedback-controlled solute transport with further valuable features, such as multiple steady states and hysteresis transitions.^35–38,64 However, more research is needed here to understand the role of the membrane feedback, especially in the presence of external driving, and how hysteresis transitions can occur.

For nonresponsive polymer membranes, we have previously shown that the Smoluchowski equation⁶⁵ well describes the solute flux and concentration profiles under stationary nonequilibrium conditions.⁶⁶ Therein, we reported that the membrane’s solute uptake, c_in, not only is a function of the polymer volume fraction, ϕ, the membrane permeability, $Pmem(ϕ)$⁠, and bulk concentration, c₀, but also is tunable in nonequilibrium by an external driving force, f. The latter leads to a nonlinear flux, j(f), with significant differences between the low- and high-force regimes. The nonlinear intermediate crossover was quantified using the newly introduced system’s differential permeability, $PsysΔ∝dj/df$⁠.

Motivated by the above features and open questions, in this work, we turn our attention to polymer membranes that are responsive to the penetrants and highlight the key differences compared to nonresponsive membranes. Specifically, we include a mean-field model for the polymer response in the Smoluchowski framework; i.e., ϕ → ϕ(c_in) is a sigmoidal function of the average solute uptake, which enters $Pmem(ϕ)$ and, in turn, controls c_in, leading to a membrane-intrinsic feedback mechanism (Fig. 1). Eventually, we use empirical expressions for $Pmem(ϕ)$ to study the feedback effect on j and $PsysΔ$ as a function of c₀ and f. Compared to nonresponsive membranes, we find a substantial enhancement in the nonlinear characteristics, such as an order of magnitude change in j due to a very small change in f, and report the emergence of multiple steady sates, bifurcations, and hysteresis in the force–flux relations.

In fact, there exist many extended versions of empirical formulas for D_in and $K$⁠, which take into account further microscopic details, such as the chemistry, shape and size of the solutes, the solvent and the membrane types, and the architecture of the polymer network.^62,69–76 However, we use Eqs. (5) and (6) to explain the feedback effects of responsive polymers at the simplest level of detail.

Furthermore, the presented framework assumes that the equilibrium findings for $K$ and D_in are also valid under moderate nonequilibrium conditions; i.e., they are independent of the flux or the force. The validity of this assumption was demonstrated in our preceding work with nonequilibrium coarse-grained simulations of membrane–bulk systems.⁶⁶ 

The detailed behavior of Eq. (10) is described in Sec. III with an appropriate choice of the model parameters.

Note that Eq. (11) assumes continuous transitions although hysteresis has been reported in experimental studies.^49,78,80,81 However, this work will demonstrate that hysteresis and bistability can result from the mutual dependencies between ϕ and c_in.

Furthermore, Eq. (11) is a mean-field approach as it does not resolve spatial inhomogeneities of ϕ and c_in [cf. the example concentration profile in Fig. 2(a)]. We assume that the system is small and that the membrane width does not change in the direction of the solute flux. Despite the multiple assumptions, our simplified model enables the investigation of the effect of responsive membrane permeability on the transport.

All length scales are expressed in units of σ, the effective diameter of one monomer. We set the system size to L = 100σ and fix the membrane width to d = 90σ; i.e., the boundary layers between the membrane and the two bulk reservoirs with concentration c₀ have a width of 5σ. The concentrations c₀, c_in, and the transition width Δc are rescaled by the crossover concentration c_c of the volume phase transition [Eq. (11)]. We choose three different transition widths, $Δc/cc∈0.1,1.0,10$⁠. The force, βfσ, is rescaled by the thermal energy and the solute size. The solute diffusivity inside the membrane, D_in, and the permeability, $P$⁠, are expressed in units of the solute bulk diffusivity D₀. The parameters A and B, which enter D_in(ϕ) [Eq. (5)] and $K(ϕ)$ [Eq. (6)], respectively, and the limits of the polymer volume fraction, ϕ_min = 0.05 and ϕ_max = 0.35, are based on our group’s coarse-grained simulations.^69,71 We fix A = 5, assuming that the diffusion is dominated by steric exclusion. The interaction parameter B is chosen in a way to yield three different values for the equilibrium membrane permeability at ϕ_c = 0.2, and, hence, we denote the membranes as repulsive $(Pmem(ϕc)/D0=0.1)$⁠, neutral $(Pmem(ϕc)/D0=1.0)$⁠, and attractive $(Pmem(ϕc)/D0=10.0)$⁠. In fact, due to typical canceling effects of $K(ϕ)$ and D_in(ϕ),^69,71 we can safely assume that the permeability of our neutral membrane does not significantly deviate from unity throughout the range of ϕ. In a system with Lennard-Jones (LJ) interactions between the solutes and the membrane monomers of equal size, the characteristic LJ interactions strengths would take the approximate values βɛ ≈ 0.03 (repulsive), βɛ ≈ 0.55 (weakly attractive), and βɛ ≈ 0.9 (attractive), respectively. All parameter values and relevant quantities are summarized in Table I.

From Eq. (10), the low- and high-force limits for the mean solute concentration inside the membrane, c_in, can be deduced, which has been discussed and substantiated with the concentration profiles from theory and coarse-grained simulations in our previous work.⁶⁶ Here, we recapture the main findings and discuss the results for the parameters used in this work.

In Fig. 3, we depict c_in [Eq. (10)] as a function of ϕ for different values of $βfσ∈0.01,10$ (color-coded) and for three different values of $Pmem(ϕc)/D0∈0.1,1.0,10.0$ (different panels). In the zero force limit, c_in reduces to the expected equilibrium value, $limf→0cin=c0K(ϕ)$⁠, which monotonously decreases for the repulsive membrane [panel (a)] and increases otherwise [panels (b) and (c)]. The same limiting result is obtained, if $Pmem(ϕ)=D0∀ϕ$⁠, which applies approximately for the “neutral” membrane [panel (b)]. In the high-force limit, the concentration profiles become piecewise constant with $limf→∞cin=c0D0Din−1(ϕ)$ and one further finds that lim_f→∞j = c₀βfD₀ = c_inβfD_in for all membrane types, since it is independent of $K$⁠.

The solute uptake of the repulsive membrane at fixed ϕ increases with f, while it decreases in the attractive membrane [see Figs. 3(a) and 3(c)]. In the “neutral” membrane [panel (b)], c_in shows no significant force dependence.

With Eqs. (10) and (11), the feedback loop depicted in Fig. 1 is closed. We numerically obtain the steady-state solutions, $(cin*,ϕ*)$⁠, by finding the intersection points of c_in(ϕ, f) and ϕ(c_in) in the phase plane. Here, we show the results with only the attractive membrane and demonstrate the general procedure. (For the repulsive membrane, we show a representative phase plane in Fig. 7 in  Appendix C.)

In Fig. 4, the black lines depict the polymer’s volume phase transition of ϕ(c_in), Eq. (11), for three different values of $Δc/cc∈0.1,1,10$⁠. The color-coded lines in Fig. 4 depict c_in(c₀, f, ϕ), Eq. (10), for three different bulk concentrations, c₀, from high [panel (a)] to low values [panel (c)]. Obviously, changing c₀ performs a shift of c_in(c₀, f, ϕ) along the horizontal axis. In panel (a), we choose c₀ = c_cD_in(ϕ_c)/D₀ ≈ 0.29 such that the high-force limit of c_in intercepts with the crossover point (c_c, ϕ_c). In panel (c), we impose that the low-force limit intercepts with the crossover point, i.e., $c0=cc/(K(ϕc)≈0.03$⁠. In panel (b), the geometric mean, $c0=ccDin(ϕc)/(D0K(ϕc))≈0.09$⁠, is used as the intermediate probe concentration.

For fixed f, c₀, and Δc, we find one or three interception points of c_in(c₀f, ϕ) and ϕ(c_in), yielding the steady-state solutions, $cin*,ϕ*$⁠. In Fig. 4(b), for instance, the low-force limit (blue dotted line) intercepts with the black solid line (Δc/c_c = 1) only once at ϕ* ≈ ϕ_∞, while it has three interceptions with the black dotted line (Δc/c_c = 0.1) at $ϕ1*≈ϕmin$⁠, $ϕ2*≈0.15$⁠, and $ϕ3*≈ϕmax$⁠. In the case of the triple solutions, the intermediate one is an unstable solution, while the other two are stable solutions. Precisely, the latter correspond to asymptotically stable solutions of the time-dependent Smoluchowski equation, $ċ(z,t)=−∂j(z,t)/∂z$⁠; i.e., the steady-state solution, c*(z), is restored after a small perturbation. A more detailed discussion on the stability of multiple solutions and the consequences for the bistable domains is provided in a separate section below.

We summarize the steady-state solutions for the attractive membrane by plotting ϕ*(f) for different values of c₀ and Δc in Fig. 5. We observe a swelling (a decrease in ϕ) with hysteresis due to an increase in f [see panels (a), (b), and (e)]. In more detail, a higher c₀ can shift the force-induced ϕ-transition to higher force values [e.g., compare panels (a) and (b)] and it determines whether transitions may occur at all. For instance, as in the case of a low transition sharpness, Δc/c_c = 10, and a low solute concentration [panel (c)], there is no significant effect on ϕ*. Similarly, for a moderate sharpness, Δc/c_c = 1, no transition is induced if c₀ is too high [panel (d)].

Furthermore, Δc plays an important role as it tunes the width of the bistable domains; e.g., while bistability is observed only in small force ranges for weakly responsive membranes (Δc/c_c = 10) [see Figs. 5(a) and 5(b)], it can exist in the entire force range for sufficiently sharp transitions (Δc/c_c = 0.1) [see Figs. 5(g) and 5(h)].

We already conclude that the membrane’s feedback response can lead to large bistable domains in ϕ tuned by f and c₀, which is characterized by a drastic switching of the membrane properties, such as the permeability, due to the bifurcations at the critical values.

We make use of $PsysΔ(f)$ to highlight the novel nonlinear effects on j(f) caused by the membrane’s feedback response. We limit ourselves to the very sharp membrane response (Δc/c_c = 0.1) and point out the significant difference between the fluxes in attractive and repulsive membranes.

In Fig. 6, we present the heat maps of $PsysΔ$ in the f–c₀ plane for the attractive [panel (a)] and the repulsive membrane [panel (h)]. The white lines labeled with Roman numerals depict the selected values of c₀ and correspond to the panels below, showing j and $PsysΔ$⁠. The heat maps share the same color-scale (see the color bars), allowing a direct comparison between the results of the attractive and repulsive membranes.

The attractive membrane exhibits, in general, larger $PsysΔ$ values than the repulsive one, particularly in the low-force and collapsed regime (ϕ = ϕ_max), in which the influence of $Pmem$ is the greatest. For f → 0, we find $PsysΔ≈7D0$ and $PsysΔ≈0.01D0$ for the attractive and repulsive membranes in the collapsed state, respectively. In the high-force limit, the system permeability converges to D₀ irrespective of the volume phase. If the membrane is swollen (ϕ = ϕ_min), the permeability of the repulsive and attractive membranes is of the same order of magnitude, i.e., $Pmem(ϕmin)≈0.6D0$ and $Pmem(ϕmin)≈2D0$⁠, respectively, and $PsysΔ$ does not deviate significantly from bulk diffusivity D₀, even for low forces (compare the blue lines in the lower panels of Fig. 6).

Due to the sharp response with Δc = 0.1c_c, the membrane is either fully swollen (ϕ_min) or fully collapsed (ϕ_max). Moreover, this also leads to large bistable domains in the c₀–f plane, visualized as striped patterns in Figs. 6(a) and 6(h). Crossing the boundary of these domains can lead to a discontinuous volume phase transition accompanied with an order of magnitude change in the solute flux (examples I, III, V, and VI). In the case of the repulsive membrane, the flux can be switched even by two orders of magnitude, particularly for small f.

In examples I–VI [Figs. 6(b)–6(g) and 6(i)–6(n)], we also depict the results for nonresponsive membranes in the fully swollen and collapsed case for comparison. In the nonresponsive case, $PsysΔ$ can also be tuned in the same range by only controlling f. Nonetheless, the membrane’s responsiveness brings about more dramatic effects, such as bistability (I, III, V, and VI) and hysteresis (V), yielding new control mechanisms to switch between two flux states. While nonresponsive membranes require large forces to exhibit bulk-like properties (D₀), a transition to this neutral state can be achieved in responsive membranes in a much sharper fashion and for lower force values; e.g., see I, V, and VI in Fig. 6. In V, for example, the crossover from the low- to high-permeability state occurs abruptly around βfσ ≈ 0.2 (the red solid line), whereas for the collapsed, nonresponsive membrane, a gradual change is observed in the range βfσ ≈ 0.1 − 1.0.

Furthermore, even if the polymer volume phase transition occurs without bifurcation, nonlinearities in the force–flux relations can be significantly enhanced. For instance, in panel (l) (example IV), we find a tenfold maximization of $PsysΔ≈12D0$ at roughly βfσ ≈ 0.5, which even exceeds the maximum differential system permeability measured for the attractive membrane.

Our model results in well-defined force–flux relations in the domains with unique steady-state solutions. In the bistable domains, however, the question arises whether the states coexist or whether only one survives under real conditions. So far, the solutions have been simply deduced from a deterministic interpretation of the macroscopic model equations, i.e., by evaluating the self-consistency equation c_in = R(c_in), with R(c_in) = d⁻¹ ∫_inc(z, ϕ(c_in))dz. The steady state is asymptotically stable if $dR(cin)/dcin|cin*<1$⁠.⁸² However, our approach neglects larger fluctuations in ϕ and c_in and does not further analyze nonequilibrium extremum principles.⁸³ In the following, we first discuss the consequences of the deterministic perspective and, then, briefly review alternative interpretations.

In the case of negligible fluctuations, the (deterministic) transition between states can be either reversible or irreversible.³⁵ One reversible transition is example V [Figs. 6(j) and 6(m)]. Here, the membrane is in the collapsed state (red line) for a small f. With increasing f, the membrane is driven into the bistable regime, yet remains in the collapsed state. Only if f exceeds the second bifurcation line, the membrane swells. In the same example V, if f is decreased from high force values, the membrane stays swollen in the bistable domain and returns to the collapsed state after the first bifurcation line is passed. Hence, we find a reversible transition with hysteresis between the two states in V.

In contrast, consider example I or VI and assume that the membrane is in the collapsed state at f = 0, an increase in f leads to a swelling when the bifurcation line is crossed. Decreasing the force again, however, does not induce a collapse, and only the swollen case survives. Hence, this transition can be termed irreversible in the deterministic interpretation. Analogously, in example III, a collapsed-to-swollen transition cannot be induced by increasing f.

Although two stable states may coexist in the deterministic model, one of them could be metastable and practically unoccupied under experimental conditions. In the literature, one finds nonequilibrium principles, e.g., based on the maximization of entropy, the minimization of entropy production (least dissipation), or the minimization of power, providing various possible routes.^31,83–91

Such extremum principles may lead to unique solutions in the bistable regime and to different values for f and c₀, where the switching between the high and low flux states occurs. For example, it should be the minimum flux, if the least-dissipation principle applies. This has direct consequences on the flux–force relation and the critical transition values of f and c₀, where the phase transition in example V would occur always at the first bifurcation line, i.e., without bistability and hysteresis.

Furthermore, the presented diffusion process can also be modeled with the stochastic Smoluchowski equation⁹² and, possibly, further coarse-grained into a stochastic differential equation for $ċin$⁠.^93,94 Hence, given the fluctuations are large enough, a stochastic switching between the two steady states may be observed in the bistable domains, and the effective force–flux relations are determined by the averaged values of ϕ and c_in. Consequently, changing f results in a continuous transition between the two states, implying that example V does not exhibit a hysteresis behavior but is rather similar to the transition in example IV.

Ultimately, the appropriate stability interpretation remains to be verified and is likely specific to the membrane material and the experimental nonequilibrium conditions. In all scenarios, a strong amplification of nonlinear characteristics and a critical switching in the force–flux relations can be expected due to the membrane’s responsiveness.

The present work focuses on the steady-state results; however, in principle, the membrane response also affects the dynamical properties as described by the time-dependent Smoluchowski equation. Precisely, the instantaneous change in the membrane permeability determines how the steady state is reached.

Numerical time-dependent solutions (not presented) and the general stability analysis for the nonlinear Fokker–Planck equations⁸² suggest that the relaxation to the stable steady states is usually overdamped. However, the fine-tuning of the model parameters may lead to underdamped oscillations. For example, consider an initially unpopulated repulsive membrane in the collapsed state. The driven solutes first accumulate in the boundary layer due to the low $Pmem$ and are only slowly driven into the membrane. When the polymer swells, $Pmem$ increases and the accumulated solutes are suddenly pushed into the membrane, leading to an overshoot of the inside concentration $(cin(t)>cin*)$⁠. If the driving is high enough, i.e., if the concentration profile cannot relax fast enough, the solutes are pushed out of the membrane, resulting in an instantaneous inside concentration that is lower than the steady-state solution $(cin(t)<cin*)$⁠. If this decrease in concentration, in turn, causes the membrane to shrink, the cycle restarts. However, since some of the particles remain in the membrane during this process, the solute accumulation in the boundary layer is reduced with each cycle and the oscillations eventually become extinct.

In order to generate sustained oscillations within our model framework, the coupling to further time-dependent variables is required. This could possibly be achieved by resolving the swelling and shrinking relaxations, $ϕ̇(t)$⁠, of the polymer or by engaging further solute species and chemical reactions. In fact, employing the permeability hysteresis in the control of (otherwise, non-oscillatory) reaction–diffusion systems is known to produce excitability and autonomous oscillations, as first proposed by theory^31,32,35,64 and eventually validated by experiments.^36–38 

We have investigated the driven steady-state solute transport through polymeric membranes with a sigmoidal volume phase response to the penetrant uptake. The change in the polymer volume fraction is decisive for the membrane permeability, which we modeled with exponential functions. This, in turn, impacts the solute uptake, leading to novel feedback-induced effects in force–flux relations that cannot be achieved by nonresponsive membranes. We quantified our findings in terms of system’s differential permeability.

The feedback effects of responsive membranes are most pronounced in the low-force regime, where the bulk concentration largely tunes the membrane density between the swollen and the collapsed state. Increasing the force can lead to a membrane swelling accompanied with a strong amplification of nonlinear characteristics and critical switching in the force–flux relations. For instance, the swelling of membranes with repulsive polymer–solute interactions can be caused by a small change in the driving force, for which we report an increase in the flux by two orders of magnitude and a pronounced maximization of the differential permeability, i.e., a tenfold increase compared to the case of nonresponsive membranes.

Moreover, of particular note is the feedback-induced coexistence of two stable steady states, while the size of the bistable domains increases with the sharpness of the sigmoidal polymer response. The bifurcations from mono- to bistability occur at critical values of the driving force and the solute bulk concentration, leading to discontinuous changes in the flux of up to two orders of magnitude.

Thus, the force-dependent switching between high and low flux states provides a valuable control mechanism that can be fine-tuned to control the appearance of hysteresis.

Hysteresis transitions found in the literature^{36–38,49,78,80,81} are usually rationalized by a bistability in the polymer’s conformational free energy^80,95,96 and attributed to the complex microscopic interactions or the competition between entropic and energetic contributions.^80,95,96 Nonetheless, many polymers exhibit a hysteresis-free response, for which the presented feedback mechanism provides a novel explanation of how hysteresis transitions can be generated and tuned in polymer membranes.

Moreover, our model is promising for the further investigation of dynamical consequences of the membrane response. In particular, the observed hysteresis in the permeability facilitates biomimetic features, such as membrane excitability and autonomous oscillations, when coupled to (non-oscillatory) chemical reactions.^{31,32,35–38,64}

Adaptations employing more complex functions for partitioning, diffusivity, and polymer response, to differing feed and permeate bulk concentrations as well as extensions to different spatial arrangements and geometries, or the investigation of the dynamics under full nonequilibrium conditions could be interesting for future studies. In fact, controlling the dynamical behavior, the critical flux switching, and the hysteresis range is crucial for the development of novel soft materials for self-regulation,⁹⁷ signal processing,^98,99 information storage, adaptation and learning,¹⁰⁰ or neuromorphic computing.²⁴ In this regard, our model approach provides accessible insights for addressing these future challenges.

The authors thank Matej Kanduč for fruitful discussions. This work was supported by the Deutsche Forschungsgemeinschaft (DFG) via the Research Unit FOR 5099 “Reducing complexity of nonequilibrium systems.” W.K.K. acknowledges the support by the KIAS Individual Grant Nos. CG076001 and CG076002 at Korea Institute for Advanced Study.

The authors have no conflicts to disclose.

Sebastian Milster: Conceptualization (equal); Data curation (lead); Formal analysis (equal); Investigation (equal); Visualization (lead); Writing – review & editing (equal). Won Kyu Kim: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Writing – review & editing (equal). Joachim Dzubiella: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

For moderate to large forces, we find S(f) ≈ exp(−βf(L − d)/2). With increasing f, the denominator in Eq. (9) converges to unity, governed by (L − d). So, the larger the d with respect to L, the higher the f has to be in order to reach the high-force limit. Obviously, the onset of the high-force limit also depends on the membrane permeability; precisely, large values of $Pmem$ will shift the crossover to smaller force values.

In Fig. 7, c_in(c₀, f, ϕ) [Eq. (10)] and ϕ(c_in) [Eq. (11)] are presented in the c₀–ϕ plane. The interception points of c_in and ϕ correspond to the force-dependent steady-state solution. The results were used to calculate ϕ*(c₀, f), which enters the flux and the differential permeability as depicted in Figs. 6(h)–6(n). The bulk concentration, c₀, and the external force, f, yield, in general, an increase in the mean inside concentration, c_in, resulting in a swelling of the membrane. For very sharp transitions, e.g., Δc = 0.1c_c, one can find multiple solutions for a fixed c₀ and f, giving rise to bistability and hysteresis.