Available online at www.sciencedirect.com [610086]

Available online at www.sciencedirect.com
ScienceDirect
Mathematics and Computers in Simulation 133 (2017) 165–174
www.elsevier.com/locate/matcom
Original articles
Homogenization results for the calcium dynamics in living cells
Claudia Timofte∗
University of Bucharest, Faculty of Physics, Bucharest-Magurele, P .O. Box MG-11, Romania
Received 29 November 2014; received in revised form 10 June 2015; accepted 28 June 2015
Available online 8 July 2015
Abstract
Via the periodic unfolding method, the effective behavior of a nonlinear system of coupled reaction–diffusion equations arising
in the modeling of the dynamics of calcium ions in living cells is analyzed. We deal, at the microscale, with two reaction–diffusion
equations governing the concentration of calcium ions in the endoplasmic reticulum and, respectively, in the cytosol, coupled
through an interfacial exchange term. Depending on the magnitude of this term, various models arise at the macroscale. In particu-
lar, we obtain, at the limit, a bidomain model. Such a model is widely used for studying the dynamics of the calcium ions, which are
recognized to be important intracellular messengers between the endoplasmic reticulum and the cytosol inside the biological cells.
c⃝2015 International Association for Mathematics and Computers in Simulation (IMACS). Published by Elsevier B.V . All rights
reserved.
Keywords: Homogenization; Calcium dynamics; The bidomain model; The periodic unfolding method
1. Introduction
Calcium is a very important second messenger in a living cell, participating in many cellular processes, such as
protein synthesis, muscle contraction, cell cycle, metabolism or apoptosis (see, for instance, [8]). Intracellular free
calcium concentrations must be very well regulated and many buffer proteins, pumps or carriers of calcium take part
in this complicated process. The finely structured endoplasmic reticulum, which is surrounded by the cytosol, is an
important multifunctional intracellular organelle involved in calcium homeostasis and many of its functions depend
on the calcium dynamics. The endoplasmic reticulum plays an important role in the metabolism of human cells. It
performs diverse functions, such as protein synthesis, translocation across the membrane and folding. This complex
and highly heterogeneous cellular structure spreads throughout the cytoplasm, generating various zones with diverse
morphology and functions. The study of the dynamics of calcium ions, acting as messengers between the endoplasmic
reticulum and the cytosol inside living cells, represents a topic of huge interest, which still requires special attention.
Many biological mechanisms involving the functions of the cytosol and of the endoplasmic reticulum are not yet
perfectly understood.
Our goal in this paper is to rigorously analyze, using the periodic unfolding method, the macroscopic behavior of a
nonlinear system of coupled reaction–diffusion equations arising in the modeling of calcium dynamics in living cells.
∗Tel.: +40 72 149 8466.
E-mail address: [anonimizat].
http://dx.doi.org/10.1016/j.matcom.2015.06.011
0378-4754/ c⃝2015 International Association for Mathematics and Computers in Simulation (IMACS). Published by Elsevier B.V . All rights
reserved.

166 C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174
Fig. 1. The reference cell Y.
We consider, at the microscale, two equations governing the concentration of calcium ions in the cytosol and, respec-
tively, in the endoplasmic reticulum, coupled through an interfacial exchange term. Depending on the magnitude of
this term, different models arise at the limit. In a particular case, we obtain, at the macroscale, a bidomain model,
which is largely used for studying the dynamics of the calcium ions in human cells. The calcium bidomain system
consists of two reaction–diffusion equations, one for the concentration of calcium ions in the cytosol and one for the
concentration of calcium ions in the endoplasmic reticulum, coupled through a reaction term. For details about the
physiological background of such a model, the reader is referred to [17]. Bidomain models arise also in other contexts,
such as the modeling of diffusion processes in partially fissured media (see [4,11] and [12]) or the modeling of the
electrical activity of the heart (see [3,2] and [19]).
Our models can serve as a tool for biophysicists to analyze the complex mechanisms involved in the calcium
dynamics in living cells, justifying in a rigorous manner some biological points of view concerning such processes.
The problem of obtaining the calcium bidomain equations using homogenization techniques was addressed by a
formal approach in [13] and by a rigorous one, based on the use of the two-scale convergence method, in [14]. Our
results constitute a generalization of some of the results contained in [13] and [14]. The proper scaling of the interfa-
cial exchange term has an important influence on the limit problem and, using some techniques from [10], we extend
the analysis from [14] to the case in which the parameter γarising in the exchange term belongs to R. Also, the tool
we use for obtaining the above mentioned macroscopic models, namely the periodic unfolding method, allows us to
treat a large class of heterogeneous media.
The layout of this paper is as follows: Section 2 is devoted to the setting of the microscopic problem. In Section 3,
we present the main convergence results, which are proven in the last section.
2. Setting of the microscopic problem
LetΩbe a bounded domain in Rn(n≥3), having a Lipschitz boundary ∂Ωformed by a finite number of con-
nected components. The domain Ωis supposed to be a periodic structure made up of two connected parts, Ωε
1and
Ωε
2, separated by an interface Γε. We assume that only the phase Ωε
1reaches the outer fixed boundary ∂Ω. Here,
εis considered to be a small positive real parameter related to the characteristic dimension of our two regions. For
modeling the dynamics of the concentration of calcium ions in a biological cell, the phase Ωε
1represents the cytosol,
while the phase Ωε
2is the endoplasmic reticulum. Let Y1be an open connected Lipschitz subset of the elementary cell
Y=(0,1)nandY2=Y\Y1(see Fig. 1). We consider that the boundary ΓofY2is locally Lipschitz and that its
intersections with the boundary of Yare reproduced identically on the opposite faces of the elementary cell. Moreover,
if we repeat Yin a periodic manner, the union of all the sets Y1is a connected set, with a locally C2boundary. Also,
we consider that the origin of the coordinate system lies in a ball contained in the above mentioned union (see [12]).
For any ε∈(0,1), let
Zε= {k∈Zn|εk+εY⊆Ω},
Kε= {k∈Zε|εk±εei+εY⊆Ω,∀i=1, . . . , n},

C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174 167
Fig. 2. A cross section of the periodic structure.
where eiare the vectors of the canonical basis of Rn. We denote (see Fig. 2)
Ωε
2=int
k∈Kε(εk+εY2)
,Ωε
1=Ω\Ωε
2
and we set
θ=Y\Y2.
Forα1, β1∈R, with 0 < α 1< β 1, we denote by M(α1, β1,Y)the collection of all the matrices A∈(L∞(Y))n×n
with the property that, for any ξ∈Rn,(A(y)ξ, ξ) ≥α1|ξ|2,|A(y)ξ| ≤β1|ξ|, almost everywhere in Y. We consider
the matrices Aε(x)=A(x/ε)defined on Ω, where A∈M(α1, β1,Y)is aY-periodic smooth symmetric matrix and
we denote the matrix AbyA1inY1and, respectively, by A2inY2.
If(0,T)is the time interval under consideration, we shall be concerned with the macroscopic behavior of the
solutions of the following microscopic system:


∂tuε−div(A1ε∇uε)=f(uε)in(0,T)×Ωε
1,
∂tvε−div(A2ε∇vε)=g(vε)in(0,T)×Ωε
2,
A1ε∇uε·ν=A2ε∇vε·νon(0,T)×Γε,
A1ε∇uε·ν=εγh(uε, vε)on(0,T)×Γε,
uε=0 on (0,T)×∂Ω,
uε(0,x)=u0(x)inΩε
1, vε(0,x)=v0(x)inΩε
2,(2.1)
where νis the unit outward normal to Ωε
1and the scaling exponent γis a given real number, related to the speed
of the interfacial exchange. As we shall see, three important cases will arise at the limit, i.e. γ=1,γ=0 and
γ= −1 (see, also, Remark 2). We assume that the initial conditions are non-negative and that the functions fandg
are Lipschitz-continuous, with f(0)=g(0)=0. We also suppose that
h(uε, vε)=hε
0(x)(vε−uε), (2.2)
where hε
0(x)=h0(x/ε)andh0=h0(y)is a real Y-periodic function in L∞(Γ), with h0(y)≥δ > 0. Besides, we
consider that
H=
Γh0(y)dσy̸=0.

168 C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174
As in [14], we can treat in a similar way the case in which the function his Lipschitz-continuous in both arguments
and is given by:
h(r,s)=h(r,s)(s−r), (2.3)
with 0 <hmin≤h(r,s)≤hmax<∞.
Throughout the paper, we denote by Ca generic fixed positive constant, whose value can change from line to line.
Since it is not easy to find an explicit solution of the well-posed microscopic problem (2.1), we need to apply a
homogenization procedure for obtaining a suitable model that describes the averaged properties of the complicated
microstructure. Using the periodic unfolding method proposed by D. Cioranescu, A. Damlamian, G. Griso, P. Donato
and R. Zaki (see [6] and [7]), we can find the asymptotic behavior of the solution of our problem. For the case γ=1,
this behavior is described by a new nonlinear system (see (3.1)), a bidomain model . So, in this case, at a macroscopic
scale, our medium can be represented by a continuous model, i.e. the superimposition of two interpenetrating contin-
uous media, the cytosol and the endoplasmic reticulum, which coexist at any point. For the other two relevant cases,
see (3.2) and (3.3).
The approach based on the periodic unfolding method allows us to work with quite general media. For our particular
geometry, we use two unfolding operators, mapping functions defined on oscillating domains into functions given on
fixed domains and, as a consequence, we do not need to use extension operators and we can deal with media possessing
less regularity than those considered usually in the literature.
It might seem that these simplified assumptions about the complex calcium dynamics inside a cell are quite strong.
However, the homogenized solution fits well with experimental data (see [17]). Also, one could argue that the period-
icity of the microstructure is not a realistic assumption and it would be interesting to work with a random microstruc-
ture. Still, such a periodic structure provides a very good description, in agreement with all the experimental findings
(see [15]).
We can deal, in a similar manner, with the more general case of a heterogeneous medium represented by a matrix
Aε
0=A0(x,x/ε)or by a matrix Dε=D(t,x/ε), under reasonable assumptions on the matrices A0andD. For
instance, we can suppose that Dis a symmetric matrix, with D, ∂tD∈L∞(0;T;L∞
per(Y))n×nand such that, for any
ξ∈Rn,(D(t,x)ξ, ξ) ≥α2|ξ|2and|D(t,x)ξ| ≤β2|ξ|, almost everywhere in (0,T)×Y, for 0 < α 2< β 2.
3. The main convergence results
In this section, we shall present the effective behavior of the solutions of the microscopic model (2.1) for the three
important cases mentioned above.
(a) Let us deal first with the case γ=1.
Theorem 1. The solution (uε, vε)of system (2.1) converges in the sense of (4.5), asε→0, to the unique solution
(u, v)of the following macroscopic problem:


θ ∂tu−div(A1∇u)−H(v−u)=θf(u)in(0,T)×Ω,
(1−θ) ∂ tv−div(A2∇v)+H(v−u)=(1−θ)g(v) in(0,T)×Ω,
u(0,x)=u0(x), v( 0,x)=v0(x)inΩ.(3.1)
Here,
H=
Γh0(y)dσy
andA1andA2are the homogenized matrices, given by:
A1
i j=
Y1
a1
i j+n
k=1a1
ik∂χ1j
∂yk
dy,
A2
i j=
Y2
a2
i j+n
k=1a2
ik∂χ2j
∂yk
dy,

C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174 169
where a1
i j=A1
i j,a2
i j=A2
i jandχ1k∈H1
per(Y1)/R, χ2k∈H1
per(Y2)/R,k=1, . . . , n, are the weak solutions of the
cell problems
−divy(A1(y)(∇yχ1k+ek))=0,y∈Y1,
A1(y)(∇yχ1k+ek)·ν=0,y∈Γ,
−divy(A2(y)(∇yχ2k+ek))=0,y∈Y2,
A2(y)(∇yχ2k+ek)·ν=0,y∈Γ.(3.2)
At a macroscopic scale, we obtain a continuous model, a so-called bidomain model , similar to those arising in the
context of the modeling of diffusion processes in partially fissured media (see [4] and [12]) or in the case of the
modeling of the electrical activity of the heart (see [3,2] and [19]). If we assume that his given by (2.3), then, at the
limit, the exchange term appearing in (3.1) is of the form |Γ|h(u, v).
(b) For γ=0, i.e. for high contact resistance, we obtain, at the macroscale, only one concentration field. So,
u=v=u0andu0is the unique solution of the following problem:

∂tu0−div(A0∇u0)=θf(u0)+(1−θ)g(u0)in(0,T)×Ω,
u0(0,x)=u0(x)+v0(x)inΩ.(3.3)
Here, the effective matrix A0is given by:
A0
i j=
Y1
a1
i j+n
k=1a1
ik∂χ1j
∂yk
dy+
Y2
a2
i j+n
k=1a2
ik∂χ2j
∂yk
dy,
in terms of the functions χ1k∈H1
per(Y1)/R, χ2k∈H1
per(Y2)/R,k=1, . . . , n, weak solutions of the local problems
−divy(A1(y)(∇yχ1k+ek))=0,y∈Y1,
A1(y)(∇yχ1k+ek)·ν=0,y∈Γ,
−divy(A2(y)(∇yχ2k+ek))=0,y∈Y2,
A2(y)(∇yχ2k+ek)·ν=0,y∈Γ.(3.4)
In this case, the exchange at the interface leads to the modification of the limiting diffusion matrix, but the insulation
is not enough strong to impose the existence of two different limit concentrations.
(c) For the case γ= −1, i.e. for very fast interfacial exchange of calcium between the cytosol and the endoplasmic
reticulum (i.e. for weak contact resistance), at the limit, we also obtain u=v=u0and, in this case, the effective
concentration field u0satisfies:
∂tu0−div(A0∇u0)=θf(u0)+(1−θ)g(u0)in(0,T)×Ω,
u0(0,x)=u0(x)+v0(x)inΩ.(3.5)
The effective coefficients are given by:
A0,i j=
Y1
a1
i j+n
k=1a1
ik∂w1j
∂yk
dy+
Y2
a2
i j+n
k=1a2
ik∂w2j
∂yk
dy,
where w1k∈H1
per(Y1)/R, w2k∈H1
per(Y2)/R,k=1, . . . , n, are the weak solutions of the cell problems


−divy(A1(y)(∇yw1k+ek))=0,y∈Y1,
−divy(A2(y)(∇yw2k+ek))=0,y∈Y2,
(A1(y)∇yw1k)·ν=(A2(y)∇yw2k)·ν, y∈Γ,
(A1(y)∇yw1k)·ν+h0(y)(w 1k−w2k)= − A1(y)ek·ν, y∈Γ.(3.6)
It is important to notice that the diffusion coefficients depend now on h0. A similar result holds true for the case in
which his given by (2.3). Let us notice that in this case the homogenized matrix is no longer constant, but it depends
on the solution u0. A similar effect was noticed in [1]. 

170 C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174
Remark 2. For simplicity, we address here only the relevant cases γ= − 1,0,1. For the case γ∈(−1,1), we get,
at the limit, the macroscopic problem (3.2), while for γ > 1, we obtain a problem similar to (3.1), but without the
exchange term H(v−u)or|Γ|h(u, v), respectively. Finally, for the case γ <−1, we obtain, at the limit, a standard
composite medium without any barrier resistance. It is worth mentioning that in this case we get w1k=w2konΓ, for
k=1, . . . , n.
Let us notice that the geometry of the domain Ωis crucial. For instance, as proven in a counterexample of H.K.
Hummel [16], if both phases are disconnected, for γ > 1, the sequence of solutions of the microscopic problem might
diverge. Also, if we deal with a different geometry, i.e. if we suppose that the phase Ωε
1is still connected, but Ωε
2is
disconnected, it follows that the homogenized matrix A2=0 and, for γ=1, the system (3.1) consists of a partial
differential equation, coupled with an ordinary differential one, which, in some particular cases, can be solved, leading
us to only one partial differential equation with memory. 
Remark 3. The conditions imposed on the nonlinear functions f,gandhcan be relaxed. For instance, we can
consider that fandgare maximal monotone graphs, verifying suitable growth conditions (see [9]). Also, as in [18,20]
and [23], we can work with more general functions h.
Remark 4. Let us notice that we can also treat the case in which the initial conditions depend on ε. A standard choice
is to consider that uε(0,x)=uε
0(x)inΩε
1,vε(0,x)=vε
0(x)inΩε
2, with uε
0(x)∈L2(Ωε
1), vε
0(x)∈L2(Ωε
2). We
assume that the extensions by zero of these functions verify uε
0⇀u0,vε
0⇀ v 0, weakly in L2(Ω). In this situation,
forγ=1 (see Theorem 1.), we get u(0,x)=θu0(x)andv(0,x)=(1−θ)v0(x), while for the other cases, namely
(b) and (c), we obtain u(0,x)=θu0(x)+(1−θ)v0(x)andu0(0,x)=θu0(x)+(1−θ)v0(x).
4. Proof of the main results
Let us introduce now the function spaces and norms we shall work with in the sequel. Let
H1
∂Ω(Ωε
1)= {v∈H1(Ωε
1)|v=0 on∂Ω∩∂Ωε
1},
V(Ωε
1)=L2(0,T;H1
∂Ω(Ωε
1)), V(Ωε
1)= {v∈V(Ωε
1)|∂tv∈L2((0,T)×Ωε
1)},
V(Ωε
2)=L2(0,T;H1(Ωε
2)), V(Ωε
2)= {v∈V(Ωε
2)|∂tv∈L2((0,T)×Ωε
2)},
with
(u(t), v(t))Ωεα=
Ωεαu(t,x)v(t,x)dx, ∥u(t)∥2
Ωεα=(u(t),u(t))Ωεα,
(u, v)Ωεα,t=t
0(u(t), v(t))Ωεαdt, ∥u∥2
Ωεα,t=(u,u)Ωεα,t,
forα=1,2. Also, let
V(Ω)=L2(0,T;H1(Ω)), V(Ω)= {v∈V(Ω)|∂tv∈L2((0,T)×Ω)},
with
(u(t), v(t))Ω=
Ωu(t,x)v(t,x)dx, ∥u(t)∥2
Ω=(u(t),u(t))Ω,
(u, v)Ω,t=t
0(u(t), v(t))Ωdt, ∥u∥2
Ω,t=(u,u)Ω,t
and
V0(Ω)= {v∈V(Ω)|v=0 on∂Ωa.e. on (0,T)},V0(Ω)=V0(Ω)∩V(Ω).
Similar spaces can be defined for Γε. We use the notation (u, v)Γε=
Γεgεuvdσx, where gεis the Riemannian
tensor on Γε.

C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174 171
The variational formulation of problem (2.1) is as follows: find (uε, vε)∈V(Ωε
1)×V(Ωε
2), with (uε(0,x), vε(0,x))
=(u0(x), v0(x))∈(L2(Ω))2and
(∂tuε(t), ϕ(t))Ωε
1+(∂tvε(t), ψ( t))Ωε
2
+(A1ε∇uε(t),∇ϕ(t))Ωε
1+(A2ε∇vε(t),∇ψ(t))Ωε
2
−εγ(h(uε(t), vε(t)), ϕ( t)−ψ(t))Γε=(f(uε(t)), ϕ( t))Ωε
1+(g(vε(t)), ψ( t))Ωε
2, (4.1)
for a.e. t∈(0,T)and any (ϕ, ψ) ∈V(Ωε
1)×V(Ωε
2).
Following the same techniques used in [14], it is not difficult to prove that (4.1) is a well-posed problem and that
uεandvεare non-negative and bounded almost everywhere.
Taking (uε, vε)as test function in (4.1), integrating with respect to time and taking into account that uεandvεare
bounded and non-negative, it follows that there exists a constant C≥0, independent of ε, such that
∥uε(t)∥2
Ωε
1+ ∥vε(t)∥2
Ωε
2+ ∥∇ uε∥2
Ωε
1,t+ ∥∇ vε∥2
Ωε
2,t+εγ(h(uε, vε),uε−vε)Γε,t≤C, (4.2)
for a.e. t∈(0,T). Also, as in [14] or [23], we can see that there exists a positive constant C≥0, independent of ε,
such that
∥∂tuε(t)∥2
Ωε
1+ ∥∂tvε(t)∥2
Ωε
2≤C, (4.3)
forγ≥1 and
∥∂tuε∥L2(0,T;H−1(Ωε
1))+ ∥∂tvε∥L2(0,T;H−1(Ωε
2))≤C, (4.4)
forγ < 1. These a priori estimates will allow us to use the periodic unfolding method and to obtain the needed
convergence results in all the above mentioned relevant cases. For retrieving the macroscopic behavior of the solution
of problem (2.1), we use two unfolding operators, Tε
1andTε
2, which transform functions defined on oscillating
domains into functions given on fixed domains (see [5,6] and [10]).
Let us prove now Theorem 1. As already mentioned, applying the same tools used in [14], it follows that (4.1)
is a well-posed problem. For γ=1, using the obtained a priori estimates and the properties of the operators Tε
1
andTε
2, it follows that there exist u∈L2(0,T;H1
0(Ω)),v∈L2(0,T;H1(Ω)),u∈L2((0,T)×Ω;H1
per(Y1)),
v∈L2((0,T)×Ω;H1
per(Y2))such that, passing to a subsequence, for ε→0, we have:


Tε
1(uε)→ustrongly in L2((0,T)×Ω,H1(Y1)),
Tε
1(∇uε) ⇀∇u+ ∇ yuweakly in L2((0,T)×Ω×Y1),
Tε
2(vε) ⇀ v weakly in L2((0,T)×Ω,H1(Y2)),
Tε
2(∇vε) ⇀∇v+ ∇ yvweakly in L2((0,T)×Ω×Y2).(4.5)
Moreover, as in [14] and [21], ∂tu∈L2(0,T;L2(Ω)), ∂ tv∈L2(0,T;L2(Ω))andu∈C0([0,T];H1
0(Ω)), v∈
C0([0,T];H1(Ω)). So, u∈V0(Ω)andv∈V(Ω).
Let us mention that, in fact, under our hypotheses, passing to a subsequence, Tε
1(uε)converges strongly to uin
Lp((0,T)×Ω×Y1), for 1 ≤p<∞. As a consequence, since the Nemytskii operator corresponding to the nonlinear
function fis continuous, it follows that f(Tε
1(uε))converges to f(u). A similar result holds true for vε.
For getting the limit problem (3.1), we take in (4.1) the admissible test functions
ϕ(t,x)=ϕ1(t,x)+εϕ2
t,x,x
ε
,
ψ(t,x)=ψ1(t,x)+εψ2
t,x,x
ε
,(4.6)
withϕ1, ψ1∈C∞
0((0,T)×Ω),ϕ2∈C∞
0((0,T)×Ω;H1
per(Y1))andψ2∈C∞
0((0,T)×Ω;H1
per(Y2)). Integrating
in time and applying in each term the corresponding unfolding operators, we obtain:

172 C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174
T
0
Ω×Y1∂tTε
1(uε)Tε
1(ϕ)dxdydt+T
0
Ω×Y2∂tTε
2(vε)Tε
2(ψ)dxdydt
+T
0
Ω×Y1Tε
1(A1ε)Tε
1(∇uε)·Tε
1(∇ϕ)dxdydt
+T
0
Ω×Y2Tε
2(A2ε)Tε
2(∇vε)·Tε
2(∇ψ)dxdydt
+T
0
Ω×Γh0(y)(Tε
1(uε)−Tε
2(vε))(Tε
1(ϕ)−Tε
2(ψ)) dxdσydt
=T
0
Ω×Y1Tε
1(f(uε))Tε
1(ϕ)dxdydt+T
0
Ω×Y2Tε
2(g(vε))Tε
2(ψ)dxdydt. (4.7)
Using the above convergence results and Lebesgue’s convergence theorem, we can pass to the limit in (4.7) (see, for
details, [6,9,21] and [22]).
Thus, we get:
T
0
Ω×Y1∂tuϕ1dxdydt+T
0
Ω×Y2∂tvψ1dxdydt
+T
0
Ω×Y1A1(y)(∇u+ ∇ yu)·(∇ϕ1+ ∇ yϕ2)dxdydt
+T
0
Ω×Y2A2(y)(∇v+ ∇ yv)·(∇ψ1+ ∇ yψ2)dxdydt
+T
0
Ω×Γh0(y)(u−v)(ϕ 1−ψ1)dxdσydt
=T
0
Ω×Y1f(u)ϕ1dxdydt+T
0
Ω×Y2g(v)ψ 1dxdydt. (4.8)
Using standard density arguments, it follows that (4.8) holds true for any ϕ1∈L2(0,T;H1
0(Ω)),ψ1∈
L2(0,T;H1(Ω)),ϕ2∈L2((0,T)×Ω;H1
per(Y1))andψ2∈L2((0,T)×Ω;H1
per(Y2)). This is just the weak
formulation of the limit problem (3.1). Indeed, if we take ϕ1=0 and ψ1=0, we easily get the cell problems (3.2) and
u=n
k=1∂u
∂xkχ1k,v=n
k=1∂v
∂xkχ2k. (4.9)
Then, if we recall that θ= |Y1|, 1−θ= |Y2|, taking ϕ2=0 and, respectively, ψ2=0, and using (4.9), we obtain ex-
actly (3.1). We notice that passing to the limit, with ε→0, in the initial conditions, we get u(0,x)=u0(x), v( 0,x)=
v0(x), for all x∈Ω.
Since uandvare uniquely determined (see [14]), the above convergences for the microscopic solutions hold for
the whole sequence and this ends the proof of Theorem 1. 
Remark 5. It is worth mentioning here that in the case in which γ > 1, we have
εγ(h(uε, vε), vε−uε)Γε=εγ−1ε (h(uε, vε), vε−uε)Γε→0,
and, thus, the exchange term H(v−u)(or, respectively, |Γ|h(u, v)) disappears in the limit. 
Let us treat now the other two relevant cases. By applying Tε
1andTε
2in (4.2), we observe that
∥Tε
1(uε)−Tε
2(vε)∥L2((0,T)×Ω×Γ)≤Cε1−γ
2.
Thus, for the case γ=0, we have, at the macroscale,
u=v=u0∈V0(Ω).

C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174 173
Taking similar tests function to those used in the case γ=1, i.e.
ϕ(t,x)=φ(t,x)+εϕ1
t,x,x
ε
,
ψ(t,x)=φ(t,x)+εϕ2
t,x,x
ε
,
withφ∈C∞
0((0,T)×Ω),ϕ1∈C∞
0((0,T)×Ω;H1
per(Y1))andϕ2∈C∞
0((0,T)×Ω;H1
per(Y2)), by unfolding,
passing to the limit and using density arguments, we get
T
0
Ω×Y1∂tu0φdxdydt+T
0
Ω×Y2∂tu0φdxdydt
+T
0
Ω×Y1A1(y)(∇u0+ ∇ yu)·(∇φ+ ∇ yϕ1)dxdydt
+T
0
Ω×Y2A2(y)(∇u0+ ∇ yv)·(∇φ+ ∇ yϕ2)dxdydt
=T
0
Ω×Y1f(u0)φdxdydt+T
0
Ω×Y2g(u0)φdxdydt,
forφ∈L2(0,T;H1
0(Ω)),ϕ1∈L2((0,T)×Ω;H1
per(Y1))andϕ2∈L2((0,T)×Ω;H1
per(Y2)), which, since
u=n
k=1∂u0
∂xkχ1k,v=n
k=1∂u0
∂xkχ2k,
withχ1kandχ2ksolutions of the cell problems (3.4), leads immediately to the macroscopic problem (3.3). 
For the case γ= −1, we still obtain at the macroscale
u=v=u0∈V0(Ω).
Moreover, in this case, following the techniques from [10], one can prove that
Tε
1(uε)−Tε
2(vε)
ε⇀u−vweakly in L2((0,T)×Ω×Γ). (4.10)
Hence, using the same tests function as those used in the case γ=0, by unfolding and passing to the limit, we get
T
0
Ω×Y1∂tu0φdxdydt+T
0
Ω×Y2∂tu0φdxdydt
+T
0
Ω×Y1A1(y)(∇u0+ ∇ yu)·(∇φ+ ∇ yϕ1)dxdydt
+T
0
Ω×Y2A2(y)(∇u0+ ∇ yv)·(∇φ+ ∇ yϕ2)dxdydt
+T
0
Ω×Γh0(u−v)(ϕ1−ϕ2)dxdσydt
=T
0
Ω×Y1f(u0)φdxdydt+T
0
Ω×Y2g(u0)φdxdydt,
forφ∈L2(0,T;H1
0(Ω)),ϕ1∈L2((0,T)×Ω;H1
per(Y1))andϕ2∈L2((0,T)×Ω;H1
per(Y2)), which, using the fact
that
u=n
k=1∂u0
∂xkw1k,v=n
k=1∂u0
∂xkw2k,
withw1kandw2ksolutions of the cell problems (3.6), leads to the limit problem (3.5). 

174 C. Timofte / Mathematics and Computers in Simulation 133 (2017) 165–174
Remark 6. It is worth mentioning here that for γ <−1,
Tε
1(uε)−Tε
2(vε)
ε
L2((0,T)×Ω×Γ)≤Cε−1−γ
2→0.
Then, as in [10], it follows that u=von(0,T)×Ω×Γand we obtain, at the limit, a standard composite medium
without any barrier resistance. 
5. Conclusions
Via the periodic unfolding method, the asymptotic behavior of the solution of a system of coupled partial differ-
ential equations arising in the modeling of calcium dynamics in living cells was analyzed. At the microscale, two
reaction–diffusion equations for the concentration of calcium ions in the cytosol and, respectively, in the endoplasmic
reticulum, coupled through an interfacial exchange term, were considered. Depending on the scaling of the exchange
term, three important cases were considered. The advantage provided by our approach is that, avoiding the use of
extension operators, it allows us to deal with quite general media.
Acknowledgments
The author is grateful to the anonymous referees for their valuable comments and suggestions, which improved the
content and the presentation of this paper.
References
[1] G. Allaire, Z. Habibi, Second order corrector in the homogenization of a conductive-radiative heat transfer problem, Discrete Contin. Dyn.
Syst. Ser. B 18 (1) (2013) 1–36.
[2] M. Amar, D. Andreucci, P. Bisegna, R. Gianni, On a hierarchy of models for electrical conduction in biological tissues, Math. Meth. Appl.
Sci. 29 (2006) 767–787.
[3] M. Amar, D. Andreucci, R. Gianni, Evolution and memory effects in the homogenization limit for electrical conduction in biological tissues,
Math. Models Methods Appl. Sci. 14 (2004) 1261–1295.
[4] G.I. Barenblatt, Y .P. Zheltov, I.N. Kochina, Basic concepts in the theory of seepage of homogeneous liquids in fissured rocks (strata), Prikl.
Mat. Mekh. 24 (1960) 852–864.
[5] D. Cioranescu, A. Damlamian, P. Donato, G. Griso, R. Zaki, The periodic unfolding method in domains with holes, SIAM J. Math. Anal. 44
(2012) 718–760.
[6] D. Cioranescu, A. Damlamian, G. Griso, The periodic unfolding method in homogenization, SIAM J. Math. Anal. 40 (2008) 1585–1620.
[7] D. Cioranescu, P. Donato, R. Zaki, Asymptotic behavior of elliptic problems in perforated domains with nonlinear boundary conditions,
Asymptot. Anal. 53 (2007) 209–235.
[8] H. Coe, M. Michalak, Calcium binding chaperones of the endoplasmic reticulum, Gen. Physiol. Biophys. 28 (2009) F96–F103.
[9] C. Conca, J.I. D ´ıaz, C. Timofte, Effective chemical processes in porous media, Math. Models Methods Appl. Sci. 13 (2003) 1437–1462.
[10] P. Donato, K.H. Le Nguyen, R. Tardieu, The periodic unfolding method for a class of imperfect transmission problems, J. Math. Sci. 176
(2011) 891–927.
[11] H.I. Ene, J.L. Auriault, Macroscopic modelling of heat transfer in composites with interfacial thermal barrier, Int. J. Heat Mass Transfer 37
(1994) 2885–2892.
[12] H.I. Ene, D. Polisevski, Model of diffusion in partially fissured media, Z. Angew. Math. Phys. 53 (2002) 1052–1059.
[13] P. Goel, J. Sneyd, A. Friedman, Homogenization of the cell cytoplasm: the calcium bidomain equations, Multiscale Model. Simul. 5 (4) (2006)
1045–1062.
[14] I. Graf, M. Peter, J. Sneyd, Homogenization of a nonlinear multiscale model of calcium dynamics in biological cells, J. Math. Anal. Appl.
419 (2014) 28–47.
[15] C. Henriquez, Simulating the electrical behavior of cardiac tissue using the bidomain model, Crit. Rev. Biomed. Eng. 21 (1993) 1–77.
[16] H.K. Hummel, Homogenization for heat transfer in polycrystals with interfacial resistances, Appl. Anal. 75 (2000) 403–424.
[17] J. Keener, J. Sneyd, Mathematical Physiology, Springer, 2008.
[18] A. Muntean, M. Neuss-Radu, A multiscale Galerkin approach for a class of nonlinear coupled reaction–diffusion systems in complex media,
J. Math. Anal. Appl. 371 (2010) 705–718.
[19] M. Pennacchio, G. Savar ´e, P.C. Franzone, Multiscale modeling for the bioelectric activity of the heart, SIAM J. Math. Anal. 37 (2005)
1333–1370.
[20] M. Radu-Neuss, W. J ¨ager, Effective transmission conditions for reactions–diffusion processes in domains separated by an interface, SIAM J.
Math. Anal. 39 (2007) 687–720.
[21] C. Timofte, Multiscale analysis in nonlinear thermal diffusion problems in composite structures, Cent. Eur. J. Phys. 8 (2010) 555–561.
[22] C. Timofte, Multiscale analysis of diffusion processes in composite media, Comput. Math. Appl. 66 (2013) 1573–1580.
[23] C. Timofte, Homogenization results for ionic transport in periodic porous media, Comput. Math. Appl. 68 (2014) 1024–1031.