Publication Cover
Applicable Analysis
An International Journal
Volume 95, 2016 - Issue 12
1,867
Views
23
CrossRef citations to date
0
Altmetric
Original Articles

A discontinuous Poisson–Boltzmann equation with interfacial jump: homogenisation and residual error estimate

&
Pages 2661-2682 | Received 04 May 2015, Accepted 06 Oct 2015, Published online: 04 Nov 2015

Abstract

A nonlinear Poisson–Boltzmann equation with inhomogeneous Robin type boundary conditions at the interface between two materials is investigated. The model describes the electrostatic potential generated by a vector of ion concentrations in a periodic multiphase medium with dilute solid particles. The key issue stems from interfacial jumps, which necessitate discontinuous solutions to the problem. Based on variational techniques, we derive the homogenisation of the discontinuous problem and establish a rigorous residual error estimate up to the first-order correction.

AMS Subject Classifications:

1. Introduction

In this paper, we consider the steady-state problem of a nonlinear Poisson–Nernst–Planck (PNP) system, which describes multiple concentrations of charged particles (e.g. ions) subject to a self-consistent electrostatic potential calculated from Poisson’s equation. In particular, we shall investigate the PNP model on a multiphase medium. The prototypical multiphase medium in mind consists of an electrolyte medium (pores), which surrounds disjoint solid particles. Such models have numerous applications describing electro-kinetic phenomena in bio-molecular or electro-chemical models, photo-voltaic systems and semiconductors, see e.g. [Citation1Citation6] and references therein. Our specific interests are motivated by models of Li-Ion batteries, see the relevant references [Citation7Citation10].

We shall deal with the nonlinearity of the model within an analytic framework, in which the PNP system can be transformed into an equivalent scalar Poisson–Boltzmann (PB) equation, see e.g. [Citation11Citation13]. This is possible, when reaction terms in the charged particle fluxes are omitted and the equations for the concentrations decouple since the charged particle concentrations are explicitly determined by the corresponding Boltzmann statistics. At the interface boundaries, this implies homogeneous Neumann conditions, which nevertheless allow for jumps in the concentrations of the species. For references applying linearisation of the PNP equations near the Boltzmann distribution see [Citation14,Citation15].

In order to be able to homogenise the equations over the entire pores-and-particles medium, a physically consistent continuation of the governing relations in the particles is needed. Here, we assume the Gouy–Chapman-Stern model for electric double layers (EDLs).[Citation16] By this, charged species underlie transport through diffusion and electrostatic drift in the pores and pure diffusion in the solid phase. This model proposes a jump of the electrostatic field across the interface (a voltage drop) and a current prescribed at the interior boundary of the solid particles.

In our general setting, we assume that the current may have a jump across the interface as well. Allowing for jumps extends the modelling to describe e.g. a separator in super-capacitors.

In the following, we will derive a discontinuous formulation of the PB equation (valid both on the volume occupied by the solid particles and on the surrounding porous space). The key feature addressed in this manuscript is the imposed inhomogeneous jump conditions, which complement the PB equation (see (8) below).

We emphasise that, from the point of view of partial differential equations, the considered PB equations are characterised by a strongly nonlinear term, which is unbounded and features exponential growth, see e.g. [Citation17].

A first aim of this paper is to establish a proper variational setting of the problem, while a second part deals with its rigorous homogenisation. With respect to the later, the averaged effective coefficients of the limit problem represent the macroscopic behaviour of the EDL, which is of primary practical importance.

For reference concerning the classic homogenisation theories, we refer to [Citation18Citation23]. The applied methods range from two-scale convergence (see e.g. [Citation24]) over Gamma-convergence (see e.g. [Citation25]) to unfolding (see [Citation26]) and others. While formal methods of averaging are widely used in the literature, their verification in terms of residual error estimates is a hard task.

From the point of view of homogenisation, the principal difficulty of discontinuous problems concerns the nonstandard boundary conditions with jumps: firstly, related jump conditions are inherent for cracks. For models and methods used in crack problems, we refer to [Citation27Citation29] and references therein. From a geometric viewpoint, cracks are open manifolds in the reference domain. Hence, classic Poincare–Friedrichs–Korn inequalities are valid in such situations. In contrast to cracks, the interfaces here are assumed to be closed manifolds disconnecting the reference domain. This fact requires discontinuous versions of Poincare–Friedrichs–Korn inequalities, which are then applied for semi-norm estimates.

Secondly, the interface boundary conditions are of Robin type. The homogenisation results known for linear problems with Robin (also called Fourier) conditions are crucially sensitive to the asymptotic rates of the involved homogenisation parameters. This issue concerns the coefficients in the boundary condition (cf. Lemma 1 below) and the volume fraction of solid particles in periodic cells (cf. Lemma2 below), see e.g. [Citation24,Citation30Citation32].

In the literature, on the one hand, the fully coupled PNP as well as decoupled PB equations are treated mostly under homogeneous boundary conditions of the Neumann type, see [Citation11,Citation12,Citation14]. On the other hand, even inhomogeneous Neumann and also Robin boundary conditions are homogenised for the most part of linear equations, see e.g. [Citation26,Citation31]. In e.g. [Citation4,Citation6] only the pore model is described and extended in the solid phase inexplicitly.

Homogenisation of transmission problems with interface jumps can be found in works concerning models of diffusion in a partially fractured porous medium, models of heat conduction in composite materials or in models of absorption of a dissolved chemical in a fluid flowing through a porous medium. We refer to [Citation33Citation36] for linear transmission problems as well as for linear problems under nonlinear [Citation37] and contact [Citation38] boundary conditions, and to [Citation39] for a semi-linear case with bounded nonlinearity.

In comparison to the works mentioned above, the principal challenge of the underlying, discontinuous PB equation is the strongly nonlinear term together with nonstandard jump conditions of Robin type, which allow not only discontinuous fields but also discontinuous fluxes across a separator. In the present contribution, we homogenise it and derive the averaged limit problem. A further major result is a rigorous estimate of the residual error up to the first-order correction.

For these purposes, we develop a variational technique based on orthogonal Helmholtz decomposition following the lines of [Citation21,Citation23]. In a periodic cell, we decompose oscillating coefficients (describing the electric permittivity) by using the nontrivial kernel in the space of vector valued periodic functions, which is represented by sums of constant and divergence free (and thus, skew symmetric) vector fields (cf. Lemma 3). Employing solutions of appropriately defined discontinuous cell problems, we obtain a regular decomposition of the homogenisation problem (see Theorem 2).

A second result establishes the critical rates of the asymptotic behaviour with respect to a homogenisation parameter for coefficients in the inhomogeneous transmission condition: we find on the one side that the critical rate for the coefficient by interfacial jumps is . This factor occurs in the discontinuous Poincare inequality (for the norm squared, cf. (Equation21) below) and is thus relevant for a coercivity estimate, which in return contributes to the solvability of the discontinuous problem and the subsequent estimate of the homogenisation error.

On the other side, the critical rate for the flux prescribed at the interior boundary of solid particles is . At this rate, the interior boundary flux induces an additional potential, which distributes over the macroscopic domain in the homogenisation limit . If the asymptotic rate is lower than the critical one, then this flux vanishes in the limit. Otherwise, if the asymptotic rate is bigger, then the flux term diverges.

From the above description, we summarise the key points of this paper as follows:

  • the study of inhomogeneous interfacial conditions describing EDL;

  • the combination of the strongly nonlinear term of exponential type, jumps and Robin conditions;

  • a variational framework of the discontinuous problem;

  • the performing of the homogenisation procedure with rigorous error estimates; and

  • the identification of the critical asymptotic rates of the boundary coefficients.

Outline: In Sections 2.12.3, we first present the problem geometry, the physical and the mathematical model. Section 2.3 establishes moreover the equivalence of the steady state of the PNP system with the scalar Poisson–Boltzmann (PB) equation and the existence of a unique solution to the PB equation (see Theorem 1).

In Section 3, we consider the homogenisation problem and the residual error estimate. At first, we state three auxiliary Lemmata before stating the main homogenisation Theorem 2.

Finally, Section 4 provides a brief discussion of the obtained results.

2. Statment of the Problem

We start with the description of the geometry.

2.1. Geometry

Let denote the domain (which is an open set) occupied by solid particles of general shape (either single or multiple particles), which are located inside the unit cell , . We assume that all particles are disjunctively located as well as bounded away from the boundary , i.e. .

We assume that the boundary is Lipschitz continuous with outer normal vector pointing away from the domain . Moreover, we distinguish the positive (outward orientated) surface and the negative (inward orientated) surface as the faces of the boundary , when approaching the boundary from outside, i.e. from or from the inside, i.e. from , respectively. For a two-dimensional example configuration see the illustration in Figure (a).

Figure 1. Two-dimensional example geometry with one star-shaped particle: (a) the unit cell and (b) the periodic disjoint domains .

Figure 1. Two-dimensional example geometry with one star-shaped particle: (a) the unit cell and (b) the periodic disjoint domains .

In the following, we consider a fixed, small homogenisation parameter and pave with periodic cells indexed by . The periodic cells are constructed from in the following way: The position of every spatial point can be decomposed as

into the integer-valued floor function coordinates and the fractional coordinates . We shall then enumerate all possible integer vectors by means of a natural ordering with the index . According to this index, we associate with the th cell and shall denote the local coordinates in all cells which correspond to .

We will denote by the respective solid particles obtained by means of the paving with for . We note that the rescaling does not change the unit outer normal vector .

Evidently, the periodic mapping , , is surjective. This construction can be generalised to an arbitrary orthotope , see [Citation26].

Let be the reference domain in with Lipschitz boundary and denote again the outer normal vector by . By reordering the index , it is then possible to account for all solid particles with the index set , see [Citation26,Citation40]. We remark that as .

By omitting solid particles which are ‘too close’ to the external boundary , we shall ensure a constant gap with the distance between and all particles . Thus, is divided into the multiple domains corresponding to all the solid particles located periodically in the reference domain and the remaining porous space .

In the following, we shall denote by the union of boundaries and introduce the multiply-connected domains

Moreover, for functions , which are discontinuous over the interface , we will denote the jump across the interface by

Here, summarises the positive faces (orientated towards the interior of the pore space ), and accounts for the negative faces (orientated towards the interior of the solid phase ).

2.2. Physical model

In the heterogeneous domain , which consist of the particle volumes and the porous space , we consider the electrostatic potential and components of concentrations of charged particles , . The physical consistency requires positive concentrations .

At the external boundary , we shall impose Dirichlet boundary conditions and corresponding to a surrounding bath and given by constant values and . We can then consider the normalised electrostatic potential and concentrations (i.e. for all ) and prescribe the following normalised Dirichlet conditions:(1)

In the following, all further relations will be formulated for the normalised potential and concentrations such that (Equation1) holds.

Let denote the electric charge of the th species with concentration for . For the -components of charges particles, we shall assume the following charge-neutrality(2)

because of the normalisation supposed in (Equation1). A necessary condition for (Equation2) is .

The charge-neutrality is used to justify the Boltzmann statistics when extending the PNP equations to the solid phase with Ohm’s law, see equations (5) and Proposition 2 below.

Moreover, the charge-neutrality assumption (Equation2) implies also the following strong monotonicity property(3)

for a constant , which follows directly from Taylor expansion with respect to , see [Citation11, Lemma 15].

We consider the following PNP steady-state system consisting of nonlinear, homogeneous equations:(4a) (4b) (5a) (5b) In both equations (4), , , denote symmetric and positive definite diffusion matrices, which are in general discontinuous over . In (Equation4b), is the Boltzmann constant, and is the temperature. We remark that the form of (Equation4b) is based on assuming the Einstein relations for the mobilities. Moreover, equation (Equation4a) models the effect of charges particles being included into the solid particles, which is well known, for instance, for -ions, see e.g. [Citation16].

In (5), denotes the symmetric and positive definite matrix of the electric permittivity, which oscillates periodically over cells according to and satisfies(6)

The entries of the permittivity matrix are discontinuous functions in the cell across the interface . A typical example considers piecewise constant in and in , with material parameters and , where denotes here the identity matrix in . In the following, we denote by , , the matrix entries of .

From a physical point of view, (Equation5a) represents Ohm’s law in the solid phase. Moreover, we remark that the equations on , i.e. (Equation4a) for and (Equation5a) for are linear while the equations (Equation4b) and (Equation5b) on form a coupled nonlinear problem on the porous space.

The modelling of boundary conditions at the interfaces is a delicate issue. For the charge carries fluxes in (4), we assume homogeneous Neumann conditions(7a) (7b) In fact, imposing inhomogeneous conditions in (7) will lead to quasi-Fermi statistics instead of the Boltzmann statistics, see (Equation16) below. In particular inter-species reaction terms would pose significant difficulties, which are out of the scope of the present work.

For the electrostatic potential in (5), we suppose the following inhomogeneous interfacial boundary conditions (see the related physics in [Citation16]):(8a) (8b) Here and are material parameters given at the interface. We note that by summing (Equation8a) and (Equation8b), we derive the relation(9)

implying that not only the electric potential but also fluxes are discontinuous functions across the interface for .

The asymptotic weights in front of and at the right-hand side of (8), which were already mentioned in the introduction, shall be discussed in detail during the below asymptotic analysis as .

We emphasise that the interface conditions (8) couple the porous phase with the solid phase by means of the jump in . In fact, the jump conditions (8) can be compared with the following two cases of simplified boundary conditions: First, if were continuous across , i.e. , then (Equation8a) and (Equation8b) would be decoupled into two usual Neumann boundary condition which do not represent the EDL. Second, if were known on the solid phase boundary , then the model would be reduced to a model on the porous space with the following inhomogeneous Robin (called also Fourier) boundary condition (see [Citation13])

However, the subsequent homogenisation of this alternative model on the porous space would nevertheless require a suitable continuation of onto .

2.3. Mathematical model

In the following, we shall amend the state variables with the superscript in order to highlight the dependency on the cell size.

The physical model will be described by the following weak variational formulation of the boundary value problem (Equation1), (4), (5), (7), (8): find an electrostatic potential and components of charge carrier concentrations such that the concentrations are positive and satisfy(10) (11) (12)

Here denotes the characteristic function of the set , and is the usual Sobolev -space defined on the multiply-connected domain and endowed with the standard norm:

Proposition 1:

For strong solutions , the variational system (Equation10)–(Equation12) and the boundary value problem (Equation1), (4), (5), (7), (8) are equivalent.

Proof:

The assertion can be verified by usual variational arguments, which we briefly sketch for the sake of the reader. The variational equations (Equation11) and (Equation12) are derived by multiplying the equations (4), (5) with test-functions and subsequent integration by parts over and due to boundary conditions (Equation1) and (7), (8).

In return, given strong solutions , the boundary value problem (4), (5), (7), (8) is obtained by varying the test-functions in (Equation11), (Equation12) and with the help of the following Green’s formulas: by recalling the denotes both the outer normal on and , we have for all (13a) (13b) which are valid on and , respectively. Hence, by suming (Equation13a) and (Equation13b), we obtain the Green’s formula representation(14)

which holds on the disjoint domain for all and , see e.g. [Citation29].

The following Proposition 2 states the crucial observation that introducing Boltzmann statistics allows to decouple the system of the homogeneous equations (Equation11) and derive an equivalent scalar PB equation.

Proposition 2:

The system (Equation10)–(Equation12) it is equivalent to the following nonlinear Poisson–Boltzmann equation: find such that(15a) (15b) together with the Boltzmann statistics determining from , i.e.(16)

Proof:

Starting with (Equation10)–(Equation12), we shall first prove the Boltzmann statistics (Equation16) by introducing the entropy variables (the chemical potentials)(17)

Then, Equation (Equation11) can be rewritten in terms of (Equation17) in divergence form as(18)

Due to the boundary condition (Equation10), we have on and the test-function can be inserted into (Equation18). Hence, by recalling that are symmetric and positive definite matrices and , we derive the identity a.e. in . Using again the boundary condition (Equation10), we conclude(19)

and is an arbitrary constant in . This fact together with (Equation17) implies (Equation16). By substituting the expressions (Equation16) into equation (Equation12) and by using the charge-neutrality (Equation2) on , equation (Equation15b) follows directly.

Conversely, the equations (Equation10)–(Equation12) follow evidently from (15) and (Equation16). This completes the proof.

We remark that the concentrations in (Equation16) are unique up to fixing the constant positive values within the solid particles .

By exploiting Proposition 2, we construct a solution for the variational problem (Equation10)–(Equation12) from the scalar problem (15) for the potential . The concentrations are afterwards explicitly determined by (Equation16).

Theorem 1:

There exists the unique solution to the semilinear problem (15) satisfying the following residual estimate(20)

which is uniform with respect to .

Proof:

We first emphasise that for the first two terms on the left-hand side of (Equation20) the following discontinuous version of Poincare’s inequality for homogeneous Dirichlet condition (Equation15a) holds on the multiple domains without interfaces (see e.g. [Citation34,Citation36]):(21)

Therefore, the lower estimate (Equation21) together with (Equation3) ensures the coercivity of the operator of the problem (Equation15b).

The main difficulty of the existence proof arises from the unbounded, exponential growth of the nonlinear term in (Equation15b). While classic existence theorems on quasilinear equations are thus not applicable here, the solution can nevertheless be constructed by a thresholding, see e.g. [Citation17] and references therein for the details.

To derive the estimate (Equation20), it suffices to insert as the test-function in the variational equation (Equation15b) and apply (Equation3) in order to estimate below the nonlinear term at the left-hand side of (Equation15b). Finally the right-hand side of (Equation15b) can be estimated by means of the following trace theorem(22)

see [Citation30] for the details. This completes the proof.

We remark that in the following Section 3, we will refine the residual error estimate (Equation20) by means of asymptotic analysis as and homogenisation.

3. Homogenisation and residual error estimate

We start the homogenisation procedure with three auxiliary cell problems. The first two cell problems serve to expand the inhomogeneous boundary traction and the volume potential of the variational problem (15) from the porous space onto the whole domain .

The third cell problem is needed to decompose the matrix of oscillating coefficients in the cells with respect to small . This procedure will result in a regular asymptotic decomposition of the perturbation problem with a subsequent error estimate of the corrector term.

For a generic cell , we introduce the Sobolev space of functions which can be extended periodically to . This requires matching traces on the opposite faces of . Moreover, we shall denote by those periodic functions, which are discontinuous, i.e. allow jumps across the interface .

3.1. Auxiliary results

We state the first auxiliary cell problem as follows: find such that(23)

In view of the homogenisation result stated in Theorem 2 in Section 3.2 below, the auxiliary problem (Equation23) serves to expand the inhomogeneity of the boundary condition (Equation8a) given by the material parameter in terms of the weak formulation stated in (Equation15b).

The existence of a unique solution in (Equation23) follows via standard elliptic theory from the assumed properties (Equation6) of . With its help, we are able to prove the following result.

Lemma 1:

(The cell boundary-traction problem) For all test-functions : on holds the following expansion(24)

where is a linear form satisfying(25)

Proof:

We apply the auxiliary cell problem (Equation23). By inserting a constant test-function , we calculate the average value(26)

Here, and denote the Hausdorff measures of the solid particle boundary in and of the cell in , respectively.

Subtracting from (Equation23), we rewrite it equivalently as(27)

where we have added to the residuum the trivial term

In the following, we shall apply the discontinuous Poincare inequality(28)

and the Trace Theorem(29)

with , which combine to the estimate(30)

By recalling that and by applying Cauchy’s inequality to the right-hand side of (Equation27) and subsequently applying estimate (Equation30) to and , we obtain the following estimate(31)

with from (Equation6) and from (Equation30).

For a proper test-function with , we insert into (Equation27) and apply the periodic coordinate transformation , , by paving such that (recall Section 2.1). After observing that , , , we also multiply (Equation27) with the constant and use (Equation26) in order to derive

which is (Equation24) with the following right-hand side term:(32)

where we denote and .

Similarly, (Equation30) transforms into the uniform estimate(33)

with from (Equation30). We note that the first line of (Equation33) expresses the -norm by the standard homogeneity argument, see e.g. [Citation22, Appendix, Lemma 1, p.370].

Therefore, the estimate (Equation31) of yields the following estimate of (34)

Here we used (Equation6) and inequalities (Equation30) for and (Equation33) for . Then, (Equation34) follows (Equation25) with the constant , which completes the proof.

Remark 1:

We remark that Lemma 1 justifies not only the a-priori estimate (Equation22), but also refines it by specifying the limiting asymptotic term as , which consists of the constant potential distributed uniformly over .

Another idea for the proof of Lemma 1 can be referred in [Citation30,Citation33].

The next auxiliary cell problem studies the asymptotic expansion of a volume force , which is given on the porous space surrounding the solid particle . It will be applied in particular to the nonlinear term in (Equation15b), i.e. we shall consider the specific volume force in Theorem 2 below.

With (recall Section 2.1), the following unfolding operator

is well defined, see [Citation26]. For its modification near the boundaries of nonrectangular domains , see [Citation40].

For , there exists a function piecewisely composed of solutions of the following -dependent cell problems (compare with (Equation23)): find such that(35)

Lemma 2:

(Unfolding of the cell volume-force problem) For all : on holds the following expansion(36)

where is a linear form satisfying(37)

Proof:

By inserting a constant test-function into the auxiliary cell problem (Equation35), we obtain the locally averaged value of (38)

Moreover, by using the average , we can expand(39)

See [Citation41] for the analysis of expansion (Equation39) in terms of Fourier series. For fixed the residual has zero average and estimates as(40)

due to the discontinuous Poincare inequality (Equation30). By inserting (Equation39) into (Equation38), we calculate

and thus derive by using again (Equation39), i.e. (41)

After multiplying the identity (Equation41) with and integrating it over , we subtract it from (Equation35) and rewrite (Equation35) equivalently as(42)

where we have added the trivial term and the residuum shortly denotes the right-hand side terms of (Equation42).

For fixed , Cauchy’s inequality yields for the first term on the right-hand side of (Equation42)(43)

Thus, by applying the estimates (Equation40) and (Equation43) to and the discontinuous Poincare inequality (Equation30) to and , we estimate at the right-hand side of (Equation42) as(44)

where and by recalling from (Equation6) and from (Equation30).

Next, we substitute as the test-function in (Equation42) and use the property of the unfolding operator. After applying the periodic coordinate transformation , to (Equation42) similar to the proof of Lemma 1, we arrive with and at (Equation36) with(45)

where and . Similarly to (Equation44), we estimate with

hence,(46)

where we have used (Equation30) for and (Equation33) for and . Thus, (Equation46) implies the estimate (Equation37) of the residual term given in (Equation45) with

This completes the proof.

Remark 2:

We remark that the factor in (Equation36) reflects the porosity of the cell due to the presence of the solid particles . In our particular geometric setting, we have and , respectively.

The third cell problem considers the solutions of the following system of linear equations: find a vector of periodic functions with componentwise zero average such that(47)

Here, denotes the space of periodic -functions and for stands for the row-wise gradient matrix of the vector , that is

Moreover in (Equation47), yields the identity matrix. The solvability of (Equation47) follows from the symmetry and positive definiteness assumption (Equation6). The uniqueness of the solution is provided due to the constraint . Indeed, since with an arbitrary constant solves also (Equation47), the zero average condition is sufficient (and necessary) to ensure the uniqueness of the solution, see e.g. [Citation21]. Finally, the solution is smooth locally in .

Remark 3:

We remark in particular, that if would hold, then the discontinuous cell problem (Equation47) would reduce to a standard, continuous cell problem.

The system (Equation47) is essential to determine the efficient coefficient matrix of the macroscopic model averaged over . In fact, following the lines of [Citation21,Citation23], we shall establish an orthogonal decomposition of Helmholtz type for the oscillating coefficients .

The Helmholtz type decomposition is based on the left-hand side of (Equation47) defining an inner product in . Due to , the variational equation (Equation47) reads as for all , which implies that belongs to the kernel of this topological vector space. Thus, the fundamental theorem of vector calculus (the Helmholtz theorem, see e.g. [Citation23]) permits the following representation as sum of a constant matrix and divergence free fields in :(48)

where has zero average, i.e.

Thus, we obtain the following lemma:

Lemma 3:

(The cell oscillating-coefficient problem) The constant matrix of effective coefficients is determined by averaging(49)

Moreover, is a symmetric and positive definite matrix with the entries:(50)

For the transformed solution vector , which depends only on since the coefficient also depends only on , the following decomposition holds:(51)

The transformed function is deduced from the symmetric matrix with zero average . Its entries , express divergence free fields (called solenoidal in 3d) obtained by combining the derivatives , of a third-order skew-symmetric tensor in the following way(52)

It follows in particular from (Equation52) that(53)

At the interface the following jump relations hold:(54)

Proof:

The constant values of stated in (Equation49) follow from averaging (Equation48) with over and by using . The formula (Equation50) can be checked directly from the th component of equation (Equation47) tested with the function , divided by , and from (Equation49) for :

The symmetry and positive definiteness of follow straightforward from the assumption in (Equation6) of being symmetric and positive definite. The formulas (Equation52) and (Equation53) describe the fact that the columns of are divergence free. Inserting the representation (Equation48) into (Equation47) and integrating by parts yields

due to the second equality in (Equation53). Then, by choosing test-functions satisfying either or , it follows(55)

Finally, we apply the periodic coordinate transformation , , with to (Equation48) and (Equation55). With , we have for the row-wise gradient matrix and . Thus, we arrive at (Equation51) and (Equation54). The proof is completed.

3.2. The main theorem

Based on the Lemmata 1–3, we formulate the main homogenisation result:

Theorem 2:

The homogenisation of the discontinuous nonlinear PB problem under the interfacial jump conditions (15) yields the following averaged (macroscopic) nonlinear PB problem: find such that(56)

For smooth solutions and of (Equation56) and (Equation47), in the limit , the solution of (15) and the corrector term to satisfy the residual error estimate (improving (Equation20)):(57)

Proof:

First, we remark that the left-hand side of (Equation57) defines a norm in due to the lower estimate (Equation21).

Secondly, the unique solution of (Equation56) can be established by following the arguments given in the proof of Theorem 1. Moreover, the solution is smooth inside by standard arguments of local regularity of weak solutions, see e.g. [Citation22, Appendix] and references therein.

Next, we prove the residual error estimate (Equation57). Integrating (Equation56) by parts on yields the strong formulation(58)

By applying the Green formulas (Equation13a) and (Equation13b) in and , respectively, we have for all : on

and for all :

By summing these two expressions and by using the continuity of across the interface , we insert the strong formulation (Equation58) into the above right-hand sides and rewrite problem (Equation56) in the disjoint domain as follows(59)

In the following, we expand the terms in (Equation59) based on the Lemmata 1–3. By applying the decomposition (Equation51) of Lemma 3 to the integrand of the first term in the left-hand side of (Equation59), we can represent it as the following sum(60)

where we have used that

Next, the integral of the last function on the right-hand side of (Equation60) can be integrated by parts by using (Equation52) and (Equation53) to calculate(61)

with . Substituting (Equation60) and (Equation61) in (Equation59), we rewrite it(62) where the bilinear continuous forms are given by(63a) (63b)

Next, we apply Lemma 2 with and obtain the following representation of the nonlinear term in (Equation62)(64)

The boundary integral in (Equation62) can be expanded by using (Equation24) in Lemma 1, i.e.

Next, we subtract equation (Equation62) for from the perturbed equation (Equation15b) for and use the notation . Moreover, for , we remark that at . Hence in view of (Equation54). Thus, after subtracting (Equation62) from (Equation15b), we calculate using the above relations(65)

One difficulty is that cannot be substituted as test-function into (Equation65) since at the boundary . For its lifting, we take a cut-off function supported in a -neighbourhood of such that at . Hence, and as . Due to the assumed -gap between and , we remark that does not intersect .

After substitution of with into (Equation65) and by using , we obtain the equality(66)

where we introduce the form due to the cut-off function as(67)

and the short notation stands for the following terms(68)

where the nonlinear form in (Equation68) is given by(69)

From (Equation69), it can be estimated uniformly as(70)

due to the Taylor series for small . To this end we remind that we have used the assumption of the higher regularity of the weak solutions and . Therefore, and are bounded uniformly, since does not depend on , and by the homogeneity argument.

The left-hand side of (Equation66) can be estimated from below by applying the coercivity of the matrix as assumed in (Equation6) and by observing that the third term on the left-hand side is nonnegative due to the strict monotonicity of the exponential function. Altogether with (Equation21), this implies that(71)

with after recalling from (Equation6) and from (Equation21).

At this point, we remark that the right-hand side of (Equation71) is a homogeneous function of degree one with respect to the norm as the following estimates will prove. Thus, the inequality (Equation71) implies directly that the norm is bounded, which reconfirms estimate (Equation20).

However, the following argument allows to refine the asymptotic residual estimate to obtain (Equation57) as . In particular, we shall estimate above the three terms at the right-hand side of (Equation71) and then apply Young’s inequality to obtain sums of sufficiently small terms of order and constant terms, which will constitute the refined residual estimate.

At first, from the estimates (Equation25), (Equation37), (Equation70) and due to the boundedness of the bilinear form (Equation63a) for , it follows that(72)

Since and recalling the properties of the cut-off function implying , we estimate that(73)

Therefore, specifically for , and by using Young’s inequality, it follows from (Equation72) and (Equation73) that(74)

For , by using again Young’s inequality we estimate (Equation67) with an arbitrary by(75)

and the form in (Equation63b) by(76)

with an arbitrary . Therefore, by applying the estimates (Equation74), (Equation75) and (Equation76) with to (Equation71) and for suitable and such that

we conclude

for all , which yields estimate (Equation57). This finishes the proof.

4. Discussion

In the following, we shall summarise the main observations concerning the presented results.

  • We observe that the first two terms on the right-hand side of (Equation71) express the residual error near and at . These terms are asymptotically of order (as can be see by setting in (Equation75) and (Equation76)) and thus constitute the leading order in the residual error estimate (Equation57). Therefore, by constructing corrector terms in form of the respective boundary layers, which is possible in practice for polyhedral domains with rational slopes of their flat boundaries with respect to the unit cell, the -estimate (Equation57) could be improved to the order .

  • The factor appears at the jump across interface in the left-hand side of microscopic equation (Equation15b). It is controlled by the coercivity condition (Equation21). In the homogenisation limit, this term disappears in the macroscopic equation (Equation56). Instead it rather contributes to the macroscopic coefficients in (Equation50).

  • The factor in front of the inhomogeneous material parameter , which is prescribed at the solid phase boundary , presents the critical order. After averaging this factor guarantees the presence of the potential distributed over the homogeneous domain in (Equation56).

  • For variable functions distributed periodically over the interface , the decomposition yields in the limit that the constant value replaces in the averaged problem (Equation56), see e.g. [Citation42].

  • The nonlinear term appearing in (Equation56) scales with the porosity coefficient .

Acknowledgements

The research results were obtained with the support of the Austrian Science Fund (FWF) in the framework of the project P26147-N26 ‘PION: Object identification problems: numerical analysis’ and SFB F32 ‘Mathematical Optimization and Applications in Biomedical Sciences’. The authors gratefully acknowledge partial support by NAWI Graz.

Notes

No potential conflict of interest was reported by the authors.

References

  • Arnold A, Carrillo JA, Desvillettes L, et al. Entropies and equilibria of many-particle systems: an essay on recent research. Monatsh. Math. 2004;142:35–43.
  • Brinkman D, Fellner K, Markowich P, et al. A drift-diffusion-reaction model for excitonic photovoltaic bilayers: asymptotic analysis and a 2-D HDG finite element scheme. Math. Models Methods Appl. Sci. 2013;23:839–872.
  • Glitzky A, Mielke A. A gradient structure for systems coupling reaction-diffusion effects in bulk and interfaces. Z. Angew. Math. Phys. 2013;64:29–52.
  • Ray N, Muntean A, Knabner P. Rigorous homogenization of a Stokes--Nernst--Planck--Poisson system. J. Math. Anal. Appl. 2012;390:374–393.
  • Schmuck M. Modeling and deriving porous media Stokes--Poisson--Nernst--Planck equations by a multi-scale approach. Commun. Math. Sci. 2011;9:685–710.
  • Schmuck M, Bazant MZ. Homogenization of the Poisson--Nernst--Planck equations for ion transport in charged porous media. SIAM J. Appl. Math. 2015;75:1369–1401.
  • Becker-Steinberger K, Funken S, Landstorfer M, et al. A mathematical model for all solid-state Lithium-ion batteries. ECS Trans. 2010;25:285–296.
  • Dreyer W, Guhlke C, Müller R. Overcoming the shortcomings of the Nernst--Planck model. Phys. Chem. Chem. Phys. 2013;15:7075–7086.
  • Kovtunenko VA. Electro-kinetic structure model with interfacial reactions. In: Ara\’uo AL, Mota Soares CA, Moleiro F, Mota Soares CM, Suleman A, editors. Proceeding 7th ECCOMAS thematic conference on smart structures and materials SMART 2015; Lissabon: IDMEC; 2015. 20p.
  • Latz A, Zausch J, Iliev O. Modeling of species and charge transport in Li-ion batteries based on non-equilibrium thermodynamics. In: Dimov I, Dimova S, Kolkovska N,editors, Numerical methods and application. Vol. 6046, Lecture notes computer science. Heidelberg: Springer; 2011. p. 329–337.
  • Allaire G, Dufreche J-F, Mikelic A, et al. Asymptotic analysis of the Poisson-Boltzmann equation describing electrokinetics in porous media. Nonlinearity. 2013;26:881–910.
  • Burger M, Schlake B, Wolfram M-T. Nonlinear Poisson--Nernst--Planck equations for ion flux through confined geometries. Nonlinearity. 2012;25:961–990.
  • Fellner K, Kovtunenko VA. A singularly perturbed nonlinear Poisson-Boltzmann equation: uniform and super-asymptotic expansions. Math. Methods Appl. Sci. Accepted. doi: 10.1002/mma.3593
  • Allaire G, Mikelic A, Piatnitski A. Homogenization of the linearized ionic transport equations in rigid periodic porous media. J. Math. Phys. 2010;51:123103. 18 pp.
  • Looker JR, Carnie SL. Homogenization of the ionic transport equations in periodic porous media. Transp. Porous Media. 2006;65:107–131.
  • Pichler F. Anwendung der Finite-Elemente Methode auf ein Lithium-Ionen Batterie Modell [Application of the finite element method to a lithium-ion battery model]. Graz: Hochschulschrift; 2011.
  • Hess P. A strongly nonlinear elliptic boundary value problem. J. Math. Anal. Appl. 1973;43:241–249.
  • Argatov II. Introduction to asymptotic modelling in mechanics. St.-Petersburg: Polytechnics; 2004. Russian.
  • Bakhvalov NS, Panasenko GP. Homogenisation: averaging processes in periodic media. Dordrecht: Kluwer; 1989.
  • Bensoussan A, Lions J-L, Papanicolaou G. Asymptotic analysis for periodic structures. Providence (RI): AMS; 2011.
  • Oleinik OA, Shamaev AS, Yosifian GA. Mathematical problems in elasticity and homogenization. Amsterdam: North-Holland; 1992.
  • Sanchez-Palencia E. Nonhomogeneous media and vibration theory. Berlin: Springer-Verlag; 1980.
  • Zhikov VV, Kozlov SM, Oleinik OA. Homogenization of differential operators and integral functionals. Berlin: Springer-Verlag; 1994.
  • Allaire G, Damlamian A, Hornung U. Two-scale convergence on periodic surfaces and applications. In: Bourgeat AP, Carasso C, Luchhaus S, Mikelic A, editors. Proceeding of international conference on mathematical modelling of flow through porous media. Singapore: World Scientific; 1996. p. 15–25.
  • Dal Maso G. Gamma-convergence and homogenization. In: Francoise J-P, Naber G, Tsun TS, editors. Encyclopedia of mathematical physics. Oxford: Elsevier; 2006. p. 449–457.
  • Cioranescu D, Damlamian A, Donato P, et al. The periodic unfolding method in domains with holes. SIAM J. Math. Anal. 2012;44:718–760.
  • Fremiot G, Horn W, Laurain A, et al. On the analysis of boundary value problems in nonsmooth domains. Vol. 462, Dissertationes mathematicae. Warsaw: Institute of Mathematics Polish Academy of Sciences; 2009.
  • Hintermüller M, Kovtunenko VA, Kunisch K. A Papkovich--Neuber-based numerical approach to cracks with contact in 3D. IMA J. Appl. Math. 2009;74:325–343.
  • Khludnev AM, Kovtunenko VA. Analysis of cracks in solids. Southampton: WIT-Press; 2000.
  • Belyaev AG, Piatnitski AL, Chechkin GA. Averaging in a perforated domain with an oscillating third boundary condition. Sb. Math. 2001;192:933–949.
  • Cioranescu D, Donato P, Zaki R. Asymptotic behavior of elliptic problems in perforated domains with nonlinear boundary conditions. Asymptotic Anal. 2007;53:209–235.
  • Oleinik OA, Shaposhnikova TA. On the homogenization of the Poisson equation in partially perforated domains with arbitrary density of cavities and mixed type conditions on their boundary. Atti Accad. Naz. Lincei Cl. Sci. Fis. Mat. Natur. Rend. Lincei (9) Mat. Appl. 1996;7:129–146.
  • Donato P, Monsurro S. Homogenization of two heat conductors with an interfacial contact resistance. Anal. Appl. 2004;2:247–273.
  • Hummel H-K. Homogenization for heat transfer in polycrystals with interfacial resistances. Appl. Anal. 2000;75:403–424.
  • Lipton R. Heat conduction in fine scale mixtures with interfacial contact resistance. SIAM J. Appl. Math. 1998;58:55–72.
  • Orlik J. Two-scale homogenization in transmission problems of elasticity with interface jumps. Appl. Anal. 2012;91:1299–1319.
  • Amar M, Andreucci D, Gianni R, et al. Evolution and memory effects in the homogenization limit for electrical conduction in biological tissues. Math. Models Methods Appl. Sci. 2004;14:1261–1295.
  • Cioranescu D, Damlamian A, Orlik J. Homogenization via unfolding in periodic elasticity with contact on closed and open cracks. Asymptotic Anal. 2013;82:201–232.
  • Neuss-Radu M, Jäger W. Effective transmission conditions for reaction-diffusion processes in domains separated by an interface. SIAM J. Math. Anal. 2007;39:687–720.
  • Francu J. Modification of unfolding approach to two-scale convergence. Math. Bohem. 2010;135:403–412.
  • Kovtunenko VA, Kunisch K. High precision identification of an object: optimality conditions based concept of imaging. SIAM J. Control Optim. 2014;52:773–796.
  • Conca C, Donato P. Nonhomogeneous Neumann problems in domains with small holes. RAIRO Model. Math. Anal. Numer. 1988;22:561–607.