\documentclass[reqno]{amsart}
\usepackage{hyperref}

\AtBeginDocument{{\noindent\small {\em 
Electronic Journal of Differential Equations}, 
Vol. 2007(2007), No. 97, pp. 1--29.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu
or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu  (login: ftp)}
\thanks{\copyright 2007 Texas State University - San Marcos.}
\vspace{9mm}}

\begin{document}
\title[\hfilneg EJDE-2007/97\hfil Derivation of  models]
{Derivation of models of compressible miscible displacement in
partially fractured reservoirs}

\author[C. Choquet\hfil EJDE-2007/97\hfilneg]
{Catherine Choquet}

\address{Catherine Choquet \newline
 Laboratoire Analyse Topologie Probabilit\'es, CNRS UMR 6632 \\
 Universit\'e P. C\'ezanne \\
 Facult\'e des Sciences et Techniques de Saint-J\'er\^ome, Case Cour A \\
 13397 Marseille Cedex 20,  France}
\email{c.choquet@univ-cezanne.fr}

\thanks{Submitted May 29, 2007. Published July 2, 2007.}
\subjclass[2000]{76S05, 35K55, 35B27, 76M50}
\keywords{Miscible compressible displacement; porous medium;
\hfill\break\indent partially fractured reservoir;
 double porosity; homogenization}

\begin{abstract}
 We derive rigorously homogenized models for the displacement of one
 compressible miscible fluid by another in  fractured porous media.
 We denote by $\epsilon$ the characteristic size of the heterogeneity
 in the medium.  A parameter $\alpha \in [0,1]$ characterizes the
 cracking degree of the rock.
 We carefully define an adapted microscopic model which is scaled by
 appropriate powers of $\epsilon$.
 We then study its limit as $\epsilon \to 0$.
 Assuming a totally fractured or a partially fractured medium, we obtain
 two effective macroscopic limit models.
 The first one is a double porosity model. The second one is of single
 porosity type but it still contains some effects due to the partial
 storage in the matrix part. The convergence is shown using
 two-scale convergence techniques.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{remark}[theorem]{Remark}
\allowdisplaybreaks

\section{Introduction}

Most of the natural  reservoirs are characterized by the existence
of a system of highly conductive fissures together with a large
number of matrix blocks. Because the fractures can rapidly
distribute pollution over vast areas, they are perceived as
controlling the water quality of the entire aquifer, although
their storativity is usually significantly smaller than that of
the surrounding matrix. Therefore, fractured aquifers are
considered vulnerable to pollution process. One can assume that
such reservoirs possess two distinct porous structures. The
problem of flow through a fractured environment is thus primarily
a problem of flow through a dual-porosity system. Its two
components are hydraulically interconnected and can not be treated
separately. The degree of interconnection between these two flow
systems defines the character of the entire flow domain and is a
function of the hydraulic properties of each of them.


Two main approaches are usually used for the modelling of such
reservoirs. The first one treats the system as a global porous
medium with averaged porosity and permeability. It leads to a
single porosity model. The second one uses the concept of double
porosity introduced by Barenblatt et al \cite{Bar}. Water in the
matrix is considered practically immobile. The less permeable part
of the rock contributes as global sink or source terms for the
transported solutes in the fracture. These two types of model have
been rigorously derived using homogenization tools for some
displacement problems. Assuming  that the reservoir is a porous
medium with highly oscillating porosity and permeability with
respect to a parameter $\epsilon$ and passing to the limit
$\epsilon \to 0$, one obtains a single porosity model (see for
instance \cite{Cho1}). Modelling the reservoir with a periodic
structure controlled by a parameter $\epsilon >0$ wich represents
the size of each block of the matrix,  scaling the equations of
the flow in the matrix by  appropriate powers of $\epsilon$ to
represent the discontinuities of the flow,  and letting
$\epsilon \to 0$, one gets a double porosity model (see for instance
\cite{Dou,Cho2}).


By the way, the former derivation is based on the assumption of a
totally fissured medium, ignoring the heterogeneity of a natural
reservoir. The matrix of cells may also be connected so that some
flow occurs directly within the cell matrix. This is the case of a
{\it partially fissured medium} considered in \cite{DS} and
\cite{Sho}. The authors introduce two flows in the microscopic
model for the matrix. The first one is the slow scale flow usually
introduced for the double porosity model. The second one is a
global flow within the matrix. The homogenization process then
leads to a macroscopic model containing both single porosity and
double porosity characteristics.

The authors in \cite{DS} and \cite{Sho} consider the flow of
single phase  fluid. The aim of the present paper is to perform
the same work for a miscible and compressible displacement in a
partially fractured medium. In Section 2, we thus derive a new
microscopic model with two flows in the matrix. A special
attention is devoted to the respect of the miscibility assumption.
This model is consistent with the totally fractured case and with
the non fractured case. The proportion between the rapidly varying
and the slowly varying part of the flow is specified by a
parameter $\alpha \in [0,1]$. We derive in Section 3 uniform
estimates on the microscopic solutions. Choosing $\alpha=0$ or
$\alpha >0$, we let the scaling parameter $\epsilon$  tend to zero
in Sections 4 and 5 and we get the associated  macroscopic models.
We mainly use  two-scale convergence arguments.


Contrary to the case studied in \cite{DS,Sho}, we show that the
double  porosity part of the model almost disappears as soon as a
direct flow occurs in the matrix. It emphasizes in particular the
role of the dispersion tensor which models all the velocities
heterogeneity at the microscopic level. It is characteristic of a
miscible flow.

From a physical viewpoint this contribution aims to give a better
understanding of the transport processes in a more or less
fissured medium. From a more mathematical viewpoint, it gives a
unified approach for the homogenization in a fractured and in a
non fractured medium.

\section{The microscopic system}

\subsection{Decomposition of the flows}

We consider the displacement of two miscible compressible species
in a fractured porous reservoir $\Omega$. The domain $\Omega$
thus consists of two subdomains, the fissures $\Omega_f$  and the
matrix $\Omega_m$.

We begin by describing the displacement in the fracture $\Omega_f$
where we do not decompose the flow. We follow the lines of
\cite{B,Pe,DR}. We assume identical compressibilities  for the
species. We thus consider that the density of the fluid is of the
form
$$
 \rho = \rho_o \, e^{p_f}, \quad \rho_o >0,
$$
where $p_f$ is the pressure in $\Omega_f$. Let $f_i$ be the
concentration of one of the two species of the mixture  and let
$v_i$ be its velocity in the flow. The conservation of mass of
this component is expressed by the following equation in
$\Omega_f$,
$$
\partial_t (\phi_f \rho f_i)+ \mathop{\rm div} (\rho f_i v_i)
=\rho q \hat{f_i},
$$
where $\phi_f$ is the porosity of the fracture and $q\hat{f_i}$ is a
source term such that $\hat{f_1}+\hat{f_2}=1$.
We then define the average velocity of the flow
$\underline{v}_f = (f_1v_1+f_2v_2)/(f_1+f_2)=f_1v_1+f_2v_2$.
Due to the Darcy law $\underline{v}_f$ is given by
$$ \underline{v}_f = - (k_f/\mu(f_1)) \nabla p_f.$$
We neglect the gravitational terms for sake of clarity in the estimates below.
The permeability of the fracture is $k_f$.
The viscosity $\mu$ is a nonlinear function depending on one of the two
concentrations in the mixture.
We cite for instance the Koval model \cite{Ko} where $\mu$ is defined for
$c \in (0,1)$ by
$\mu(c)=\mu(0) (1+(M^{1/4}-1)c )^{-4}$, the constant $M=\mu(0)/\mu(1)$
being the mobility ratio.
Setting $j_i=\rho f_i (v_i-\underline{v}_f)$, the latter equation becomes

\begin{equation}
 \partial_t (\phi_f \rho f_i)
+ \mathop{\rm div} (\rho f_i \underline{v}_f) +\mathop{\rm
div}(j_i) =\rho q \hat{f_i}. \label{2.1}
\end{equation}
The tensor $j_i$ models the effects of the heterogeneity of the velocities
in the mixture.
Analogous to Fick's law this dispersive flux is considered proportional
 to the concentration gradient
$j_i = \rho \mathcal{D}(\underline{v}_f) \nabla f_i$, where the
dispersion tensor is

\begin{equation}
\mathcal{D}(\underline{u})
=\phi_f \bigl(D_m  Id + D_p(\underline{u})  \bigr)
= \phi_f \bigl(  D_m  Id +
|\underline{u}| \left( \alpha_l  \mathcal{E}(\underline{u}) +
  \alpha_t  (Id -  \mathcal{E}(\underline{u}) ) \right)
\bigr) ,
\label{2.2}
\end{equation}
where $\mathcal{E}(\underline{u})_{ij} = \underline{u}_i
\underline{u}_j / |\underline{u}|^2  $, $\alpha_l$ and $\alpha_t$
are the longitudinal and transverse  dispersion constants and
$D_m$ is the molecular diffusion. For the usual rates of flow,
these reals are such that $\alpha_l \ge \alpha_t \ge D_m > 0$.
Such a dispersive tensor is characteristic of a miscible model.
Using the definition of the density $\rho$ in  (\ref{2.1}),
dividing by $\rho>0$, and using the classical assumption of weak
compressibility to neglect the terms containing $\underline{v}_f
\cdot \nabla p_f$, we get
\begin{equation}
\phi_f \partial_t f_i +\phi_f f_i \partial_t p_f + \mathop{\rm
div} ( f_i \underline{v}_f) -\mathop{\rm
div}(\mathcal{D}(\underline{v}_f) \nabla f_i) =q \hat{f_i}.
\label{2.3}
\end{equation}
Bearing in mind $f_1+f_2=1$ and summing up the later relation for
$i=1,2$, we model the conservation of the total mass by
\begin{equation}
\phi_f \partial_t p_f + \mathop{\rm div}( \underline{v}_f) =q,
\quad \underline{v}_f = - \frac{k_f}{\mu(f_1)} \nabla p_f.
\label{2.4}
\end{equation}
The flow in $\Omega_f$ is then completely modelled by (\ref{2.4})
coupled with
\begin{equation}
\phi_f \partial_t f_1 + \underline{v}_f  \cdot \nabla f_1
-\mathop{\rm div}(\mathcal{D}(\underline{v}_f) \nabla f_1) =
q(\hat{f_1}-f_1). \label{2.5}
\end{equation}


We now perform a similar calculation in the matrix part of the
domain. But, following \cite{DS},  we assume that the flow of each
component  is made of two parts. The first component accounts for
the global diffusion in the pore system. The second one
corresponds to high frequency spatial variations which lead to
local storage in the matrix. The $i$th concentration in the matrix
$m_i$ is then given by
$$
m_i=\alpha c_i + \beta  C_i, \quad 0 \le \alpha < 1, \; \alpha +
\beta =1.
$$
Note that the value of the parameter $\alpha$ may be for instance
given by experimental data on samples of porous media and by
stochastic reconstruction (see \cite{Adl} and the references
therein). The partial concentrations $c_i$ and $C_i$ satisfy
\begin{equation}
 \partial_t (\phi \rho c_i)
+ \mathop{\rm div} (\rho c_i v_i) =\rho q \hat{c}_i, \quad
 \partial_t (\phi \rho C_i)
+ \mathop{\rm div} (\rho C_i V_i) =\rho q \hat{C}_i, \label{2.6}
\end{equation}
where $\phi$ is the porosity of the matrix, $\rho=\rho_o e^p$ is
the  density of the mixture and $p$ is the pressure in the matrix.
We also define two Darcy velocities
$\underline{v}=(c_1v_1+c_2v_2)/(c_1+c_2)$ and
$\underline{V}=(C_1V_1+C_2V_2)/(C_1+C_2)$. We assume that they
both obey a Darcy law depending on the pressure $p$. The rapidly
varying and the slowly varying components are distinguished by two
different permeabilities $k$ and $K$ satisfying
\begin{equation}
\underline{v} = - (k/\mu(m_1)) \nabla p, \quad
\underline{V} = - (K/\mu(m_1)) \nabla p.
\label{2.7}
\end{equation}
Introducing these Darcy velocities  in (\ref{2.6}), we obtain
\begin{equation}
\begin{gathered}
 \partial_t (\phi \rho c_i) + \mathop{\rm div} (\rho c_i
\underline{v}) + \mathop{\rm div}( \rho c_i (v_i-\underline{v}))
=\rho q \hat{c}_i,\\
 \partial_t (\phi \rho C_i)
 + \mathop{\rm div} (\rho C_i \underline{V})
+ \mathop{\rm div} (\rho C_i (V_i-\underline{V})) =\rho q
\hat{C}_i.
\end{gathered} \label{2.8}
\end{equation}
We now use the classical dispersive characteristic of a miscible displacement.
To this aim, we define a global average velocity
$\underline{\mathcal{V}} = (\alpha(c_1+c_2)\underline{v}
+ \beta (C_1+C_2)\underline{V}) / ( \alpha(c_1+c_2) + \beta (C_1+C_2))$,
that is
\begin{equation}
\underline{\mathcal{V}}
= \alpha(c_1+c_2)\underline{v} + \beta (C_1+C_2)\underline{V}.
\label{2.9}
\end{equation}
The relative velocities are assumed to satisfy a Fick's law given
by $ c_i (v_i - \underline{v})+c_i( \underline{v} -
\underline{\mathcal{V}}) = c_i (v_i - \underline{\mathcal{V}}) =
\mathcal{D}(\underline{\mathcal{V}}) \nabla c_i$ and $ C_i (V_i -
\underline{V})+c_i( \underline{V} - \underline{\mathcal{V}}) = C_i
(V_i - \underline{\mathcal{V}}) =
\mathcal{D}(\underline{\mathcal{V}}) \nabla C_i$, the tensor
$\mathcal{D}$ being  of the form (\ref{2.2}). Using the definition
of $\rho$ and neglecting the quadratic velocities terms, we then
write (\ref{2.8}) in the following form for $i=1,2$.
\begin{equation}
 \phi \partial_t d_i
+ \phi d_i \partial_t p + \mathop{\rm div}(\underline{\mathcal{V}}
d_i) - \mathop{\rm div}(\mathcal{D}(\underline{\mathcal{V}})
\nabla d_i ) =q \hat{d}_i, \quad d_i=c_i   \text{ or }   C_i.
\label{2.10}
\end{equation}
Finally we use $m_1+m_2=1$ to get a relation for the conservation
of the total mass $ \phi \partial_t p + \mathop{\rm
div}(\underline{\mathcal{V}}) =q.$ The flow in the matrix part is
thus governed by the following set of equations:
\begin{gather}
 \phi \partial_t p + \mathop{\rm div}(\underline{\mathcal{V}})
=q, \quad \underline{\mathcal{V}} = \alpha(c_1+c_2) \underline{v}
+ (1- \alpha (c_1+c_2)) \underline{V},
\label{2.11} \\
 \underline{v} = - \frac{k}{\mu(m_1)} \nabla p, \quad
\underline{V} = - \frac{K}{\mu(m_1)} \nabla p,
\label{2.12} \\
 \phi \partial_t C_1
+\underline{\mathcal{V}} \cdot \nabla C_1 - \mathop{\rm
div}(\mathcal{D}(\underline{\mathcal{V}}) \nabla C_1 ) =q
(\hat{C_1}-C_1), \label{2.13} \\
\phi \partial_t c_i +\underline{\mathcal{V}} \cdot \nabla c_i -
\mathop{\rm div}(\mathcal{D}(\underline{\mathcal{V}}) \nabla c_i )
=q (\hat{c}_i-c_i), \quad i=1,2. \label{2.14}
\end{gather}
Note that choosing $\alpha=1$ and thus $\beta=0$ and $k=K$,
this model is exactly the same as the one derived in the fracture.

\subsection{The scaled microscopic system}

Since we neglect the gravitational terms, we describe the
far-field repository by  a domain $\Omega \subset \mathbb{R}^2$ with a
periodic structure, controlled by a parameter $\epsilon >0$ which
represents the size of each block of the matrix. Note that all the
results of the paper remain true in a domain $\Omega$ of $\mathbb{R}^3$
(see Remark \ref{rmk1} below for the minor modifications of the
proof). The $\mathcal{C}^1$ boundary of $\Omega$ is $\Gamma$ and
$\nu$ is the corresponding exterior normal. As in \cite{Sho}, the
standard period ($\epsilon =1$) is a cell $Y$ consisting of a
matrix block $Y_m$ of external $\mathcal{C}^1$ boundary $\partial
Y_m$ and of a fracture domain $Y_f$. We assume that $\vert Y \vert
=1$. The $\epsilon$-reservoir consists of copies $\epsilon Y$
covering $\Omega$. The two subdomains of $\Omega$ are defined by
$$
\Omega_f^{\epsilon}= \Omega \cap \bigl\lbrace \cup_{\xi \in \mathcal{A}}
 \epsilon (Y_f + \xi )\bigr\rbrace ,
\quad \Omega_m^{\epsilon}= \Omega \cap \bigl\lbrace \cup_{\xi \in
\mathcal{A}} \epsilon ( Y_m + \xi )\bigr\rbrace ,
$$
where
$\mathcal{A}$ is an appropriate infinite lattice. The
fracture-matrix interface is denoted by $\Gamma_{fm}^\epsilon =
\partial \Omega_f^{\epsilon} \cap \partial \Omega_m^{\epsilon}
\cap \Omega$ and $\nu_{fm}$ is the corresponding unit normal
pointing out   $\Omega_f^{\epsilon}$. An example of unit cell in
the case of a partially fractured medium is given in Fig. 1.
See  \cite{Acer} for a more complicated admissible 2D-structure.

\begin{figure}[ht]
 \begin{center}
\setlength{\unitlength}{.8mm}
\begin{picture}(60,64)(0,0)
\put(0,0){\line(1,0){60}}
\put(0,0){\line(0,1){60}}
\put(60,0){\line(0,1){60}}
\put(0,60){\line(1,0){60}}
\qbezier(0,40)(20,40)(20,60)
\qbezier(40,60)(40,40)(60,40)
\qbezier(0,20)(20,20)(20,0)
\qbezier(40,0)(40,20)(60,20)
\put(8,50){$Y_f$} \put(50,50){$Y_f$}
\put(8,8){$Y_f$} \put(50,8){$Y_f$}
\put(29,29){$Y_m$}
\put(28,62){$\Gamma_{mm}$}
\put(16,16){$\Gamma_{fm}$}
\end{picture} 
 \end{center}
\caption{Unit cell of a partially fractured medium}
\end{figure}



The major difference between a partially fractured structure and a
totally  fractured one is that the matrix block $Y_m$ is not
completely surrounded by the fracture domain $Y_f$. Let us denote
by $J=(0,T)$ the time interval of interest. To homogenize the
reservoir, we shall let tend to zero the size $\epsilon$ of the
cells.



Our starting point consists of the equations derived in the latter
subsection. As we assume a periodic structure in the reservoir,
the porosities
$(\phi_f^\epsilon(x),\phi^{\epsilon}(x))=(\phi_f(\frac{x}{\epsilon}),
\phi (\frac{x}{\epsilon }))$ and the permeabilities
$(k_f^\epsilon(x),k^{\epsilon}(x))=(k_f(\frac{x}{\epsilon}),
 k(\frac{x}{\epsilon}))$ of the fracture and of the matrix are periodic
of period $(\epsilon Y_f,\epsilon Y_m)$. These  quantities are
assumed to be smooth and bounded, but globally they are
discontinuous across $\Gamma_{fm}^{\epsilon}$. A slowly varying
flow and a rapidly varying flow occur in the matrix
$\Omega_m^{\epsilon}$. The equations for the rapidly varying flow
will be scaled by appropriate powers of $\epsilon$ to conserve the
flow between the matrix and the fractures as $\epsilon \to 0$
({\it cf} \cite{Dou,DS}). Scaling  (\ref{2.4})-(\ref{2.5}) and
(\ref{2.11})-(\ref{2.14}), we get
\begin{gather}
 \phi_f^\epsilon  \partial_t f_1^\epsilon +
\underline{v}_f^\epsilon \cdot \nabla f_1^\epsilon - \mathop{\rm
div}(\mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon )
= q (\hat{f}_1-f_1^\epsilon) \quad \text{in}\  \Omega_f^\epsilon
\times J,
\label{2.15} \\
  \phi_f^\epsilon \partial_t p_f^\epsilon
+\mathop{\rm div}(\underline{v}_f^\epsilon) = q, \
\underline{v}_f^\epsilon = -
\frac{k_f^\epsilon}{\mu(f_1^\epsilon)} \nabla p_f^\epsilon \quad
\text{in}\  \Omega_f^\epsilon  \times J,
\label{2.16}  \\
 \phi^\epsilon \partial_t C_1^\epsilon +
\underline{\mathcal{V}}^\epsilon \cdot \nabla C_1^\epsilon -
\mathop{\rm
div}(\mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon )
\nabla C_1^\epsilon ) = q (\hat{C}_1-C_1^\epsilon) \quad
\text{in}\  \Omega_m^\epsilon \times J,
\label{2.17} \\
 \phi^\epsilon \partial_t c_1^\epsilon +
\underline{\mathcal{V}}^\epsilon \cdot \nabla c_1^\epsilon -
\mathop{\rm
div}(\mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon )
\nabla c_1^\epsilon ) = q (\hat{c}_1-c_1^\epsilon) \quad \text{in
}  \Omega_m^\epsilon \times J,
\label{2.18} \\
 \phi^\epsilon \partial_t c_2^\epsilon +
\underline{\mathcal{V}}^\epsilon \cdot \nabla c_2^\epsilon -
\mathop{\rm
div}(\mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon )
\nabla c_2^\epsilon ) = q (\hat{c}_2-c_2^\epsilon) \quad
\text{in}\  \Omega_m^\epsilon \times J,
\label{2.19} \\
 \phi^\epsilon \partial_t p^\epsilon +
\mathop{\rm div}(\underline{\mathcal{V}}^\epsilon  ) =q, \quad
    \underline{\mathcal{V}}^\epsilon=  \underline{\mathcal{V}}_s^\epsilon+  \epsilon \underline{\mathcal{V}}_h^\epsilon,
\quad \text{in}\  \Omega_m^\epsilon \times J,
\label{2.20} \\
 \underline{\mathcal{V}}_s^\epsilon = -\alpha
(c_1^\epsilon+c_2^\epsilon) \frac{k^\epsilon}{\mu(m_1^\epsilon)}
\nabla p^\epsilon, \;
  \underline{\mathcal{V}}_h^\epsilon = -(1-\alpha (c_1^\epsilon+c_2^\epsilon) ) \frac{\epsilon k^\epsilon}{\mu(m_1^\epsilon)} \nabla p^\epsilon,
\label{2.21}
\end{gather}
where $m_1^\epsilon= \alpha c_1^\epsilon+\beta C_1^\epsilon$, $0
\le \alpha < 1$, $\alpha + \beta=1$. The flow in the fissures is
described by \eqref{2.15}-\eqref{2.16}. The matrix behavior is
described by \eqref{2.17}-\eqref{2.21}. In particular,
\eqref{2.17} governs the slowly varying component while
(\ref{2.18})-(\ref{2.19}) governs the high frequency varying ones.
The  source term $q$ is  a nonnegative function of $L^2(\Omega
\times J)$ and
$$
\alpha \hat c_1 + \beta \hat C_1=\hat f_1, \quad
0 \le \hat f_1 \le 1, \quad
\hat c_1 + \hat c_2 =1.
$$
 We assume that the porosities and the symmetric permeability tensors satisfy
$$
0 < \phi_- \le \phi_f(x),\ \phi(x) \le \phi_-^{-1}, \quad
k_- |\xi|^2 \le k_f(x)\xi \cdot \xi, \  k(x)\xi \cdot \xi
\le k_-^{-1} |\xi|^2,
$$
$k_- >0$, a.e. in $\Omega$, for all $\xi \in \mathbb{R}^2$.
The viscosity $\mu \in W^{1,\infty}(0,1)$ is such that
$$
0 < \mu_- \le \mu(x) \le \mu_+ \quad \forall x \in (0,1).
$$
The tensor $\mathcal{D}$ is already defined in (\ref{2.2}). The
tensor $\mathcal{D}^\epsilon$ has a similar structure but its
diffusive part $(\alpha + \beta \epsilon^2)D_m Id$ contains the
same proportions of slowly and rapidly varying flows than the
matrix. The main property of these tensors is
\begin{equation}
\begin{gathered}
\mathcal{D}(\underline{v}_f^\epsilon) \xi \cdot \xi \ge \phi_-
(D_m + \alpha_t \vert \underline{v}_f^\epsilon \vert  ) \vert \xi
\vert^2  , \quad  \forall \xi \in \mathbb{R}^2,
\\
  \mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon ) \xi
\cdot \xi \ge \phi_- (D_m (\alpha +\beta \epsilon^2) + \alpha_t
\vert \underline{\mathcal{V}}_s^\epsilon + \epsilon^2
\underline{\mathcal{V}}_h^\epsilon \vert  ) \vert \xi \vert^2  ,
\quad  \forall \xi \in \mathbb{R}^2.
\end{gathered} . \label{2.22}
\end{equation}
The model is completed by the following boundary and initial conditions.
We begin by the continuity relations across  the interface
$\Gamma_{fm}^\epsilon \times J$.
\begin{gather}
\beta\mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon
\cdot \nu_{fm} =  \mathcal{D}^\epsilon(
\underline{\mathcal{V}}^\epsilon)  \nabla C_1^\epsilon \cdot
\nu_{fm},
\label{2.23} \\
\alpha\mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon
\cdot \nu_{fm} =  \mathcal{D}^\epsilon(
\underline{\mathcal{V}}^\epsilon)  \nabla c_1^\epsilon \cdot
\nu_{fm},
\label{2.24}\\
 \alpha\mathcal{D}(\underline{v}_f^\epsilon ) \nabla
(1-f_1^\epsilon )\cdot \nu_{fm} = -
\alpha\mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon
\cdot \nu_{fm} =  \mathcal{D}^\epsilon(
\underline{\mathcal{V}}^\epsilon)  \nabla c_2^\epsilon \cdot
\nu_{fm},
\label{2.25}\\
 f_1^\epsilon = \alpha c_1^\epsilon + \beta C_1^\epsilon ,
\label{2.26}\\
 \underline{v}_f^\epsilon \cdot \nu_{fm} =
\underline{\mathcal{V}}^\epsilon \cdot \nu_{fm}, \quad
p_f^\epsilon = p^\epsilon. \label{2.27}
\end{gather}
We add a zero flux condition out of the full domain $\Omega$.
\begin{gather}
 \mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon \cdot
\nu =0\  \text{on} \  \partial \Omega_f^\epsilon \cap \Gamma,
\label{2.28}\\
 \mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon ) \nabla
C_1^\epsilon \cdot \nu =
\mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon ) \nabla
c_1^\epsilon \cdot \nu =
\mathcal{D}^\epsilon(\underline{\mathcal{V}}^\epsilon ) \nabla
c_2^\epsilon \cdot \nu =0 \  \text{on} \  \partial
\Omega_m^\epsilon \cap \Gamma,
\label{2.29}\\
 \underline{v}_f^\epsilon \cdot \nu =0 \  \text{on} \  \partial
\Omega_f^\epsilon \cap \Gamma, \ \underline{\mathcal{V}}^\epsilon
\cdot \nu =0 \  \text{on} \  \partial \Omega_m^\epsilon \cap
\Gamma. \label{2.30}
\end{gather}
The initial conditions in $\Omega$ are the following.
\begin{gather}
(f_1^\epsilon(x,0),C_1^\epsilon(x,0),c_1^\epsilon(x,0),c_2^\epsilon(x,0)
)=(\chi_f^\epsilon f_1^o(x),C_1^o(x),c_1^o(x),c_2^o(x) ) ,
\label{2.31}\\
 p_f^\epsilon(x,0)= \chi_f^\epsilon(x) p^o(x), \
p^\epsilon(x,0)=\chi_m^\epsilon(x) p^o(x) . \label{2.32}
\end{gather}
We assume that $p^o$ belongs to $H^1(\Omega)$, and that $(f_1^o,C_1^o,c_1^o,c_2^o) \in L^\infty(\Omega) \times (L^2(\Omega_m^\epsilon))^3$ satisfies

\begin{gather}
 0 \le f_1^o(x) \le 1 \  \text{a.e.\  in}\  \Omega,
\label{2.33}\\
\alpha c_1^o(x)+ \beta C_1^o(x) = \chi_m^\epsilon f_1^o (x) , \ 0
\le C_\alpha \le c_1^o(x)+c_2^o(x) \le 1 \  \text{a.e.\  in}\
\Omega_m^\epsilon. \label{2.34}
\end{gather}
The constant $0 \le C_\alpha <1$ is introduced to prevent a
degeneration in the limit study of the pressure equation in
$\Omega_m^\epsilon$ (see Section 5).


\subsection{Variational formulation and existence for the
 microscopic model}

We now state an existence result for the problem
\eqref{2.15}-\eqref{2.21}, \eqref{2.23}-\eqref{2.32}. The usual
equations modelling a miscible and compressible displacement in
porous media are of the form (\ref{2.4})-(\ref{2.5}). The
existence of weak solutions for this problem is proved in
\cite{Cho3}. In the present paper, the decomposition of the flow
in the matrix part of the domain induces two additional
difficulties. The first one is a new coupling between
concentrations and pressure due to the term $\alpha
(c_1^\epsilon+c_2^\epsilon)$ in \eqref{2.21}. But from a
mathematical viewpoint this difficulty is similar to the one due
to the concentration dependent viscosity in the Darcy law. The
second novelty occurs at the interface $\Gamma_{fm}^\epsilon$. One
has to link one concentration  $f_1^\epsilon$ in the fracture with
three concentrations $(C_1^\epsilon,c_1^\epsilon,c_2^\epsilon)$ in
the matrix. Thus, following \cite{Sho}, we introduce appropriate
concentrations spaces for the problem. Let $H^\epsilon$ be the
Hilbert space $H^\epsilon=L^2(\Omega_f^\epsilon) \times
L^2(\Omega_m^\epsilon) \times L^2(\Omega_m^\epsilon)$ with the
inner product
\begin{align*}
& \bigl( [u_f,u_m,U_m],[\psi_f,\psi_m,\Psi_m]
\bigr)_{H^\epsilon}\\
&= \int_{\Omega_f^\epsilon} u_f(x) \, \psi_f(x) \, dx
 +  \int_{\Omega_m^\epsilon} u_m(x) \, \psi_m(x) \,
dx +  \int_{\Omega_m^\epsilon} U_m(x) \, \Psi_m(x) \, dx.
\end{align*}
Let $\gamma_j^\epsilon \, : H^1(\Omega_j^\epsilon) \to L^2(\partial \Omega_j^\epsilon)$ be the usual trace map and $\chi_j^\epsilon$ be the characteristic function associated with $\Omega_j^\epsilon$, $j=f,m$.
Let $V^\epsilon$ be the following Banach space
\begin{align*}
V^\epsilon &= H^\epsilon \cap  \big\{ (u_f,u_m,U_m) \in
H^1(\Omega_f^\epsilon) \times  H^1(\Omega_m^\epsilon)  \times
H^1(\Omega_m^\epsilon)  ; \\
&\quad \gamma_f^\epsilon u_f = \alpha \gamma_m^\epsilon u_m +
\beta \gamma_m^\epsilon U_m \  \text{on}\  \Gamma_{fm}^\epsilon
\big\}
\end{align*}
endowed with the norm
\begin{align*}
\Vert  (u_f,u_m,U_m) \Vert_{V^\epsilon}
& = \Vert \chi_f^\epsilon
u_f \Vert_{L^2(\Omega)} + \Vert \chi_m^\epsilon u_m
\Vert_{L^2(\Omega)} + \Vert \chi_m^\epsilon U_m
\Vert_{L^2(\Omega)}\\
& \quad + \Vert \chi_f^\epsilon \nabla u_f
\Vert_{(L^2(\Omega))^2} + \Vert \chi_m^\epsilon \nabla u_m
\Vert_{(L^2(\Omega))^2}  + \Vert \chi_m^\epsilon \nabla U_m
\Vert_{(L^2(\Omega))^2} .
\end{align*}
The introduction of similar spaces for the pressure is useless because
we only use one pressure variable in the matrix part of the domain.
Note that the pair of spaces $(H^\epsilon,V^\epsilon)$  possesses the
same ``good'' properties as $(L^2(\Omega),H^1(\Omega))$.
In particular, we have the compact embedding $V^\epsilon \subset H^\epsilon$.
Thus, adapting the proof of \cite{Cho3} to the present piecewise structure,
one can state the following existence result.


\begin{theorem} \label{thm1}
Let $0<\epsilon <1$. There exists a solution
$(p_f^{\epsilon},p^\epsilon,f_1^\epsilon,c_1^{\epsilon},
C_1^{\epsilon},c_2^{\epsilon})$ of Problem
\eqref{2.15}-\eqref{2.21}, \eqref{2.23}-\eqref{2.32} in the
following sense.

\noindent (i) The pressure part $(p_f^\epsilon,p^{\epsilon})$
belongs
 to  $L^{2}(J;H^1(\Omega_f^\epsilon)) \times L^{2}(J;H^1(\Omega_m^\epsilon))$
 and is a weak solution of \eqref{2.16}, \eqref{2.20}--\eqref{2.21},
\eqref{2.27}, \eqref{2.32}.
 Indeed, for any function $\psi \in \mathcal{C}^1(J;H^1(\Omega))$,
\begin{equation}
\begin{aligned}
&- \int_{\Omega \times J} (\chi_f^\epsilon \phi_f^\epsilon p^\epsilon_f
 + \chi_m^\epsilon \phi^\epsilon p^\epsilon) \partial_t \psi \\
&+ \int_{\Omega \times J} \bigl( \chi_f^\epsilon
  \frac{k_f^\epsilon}{\mu(f_1^\epsilon)}  \nabla p_f^\epsilon
  + \chi_m^\epsilon (\alpha (c_1^\epsilon +c_2^\epsilon)(1-\epsilon^2)
  + \epsilon^2  )\frac{k^\epsilon}{\mu(m_1^\epsilon)}
  \nabla p_\epsilon \bigr) \cdot \nabla \psi\\
&=  - \int_\Omega (\chi_f^\epsilon \phi_f^\epsilon
  + \chi_m^\epsilon \phi^\epsilon ) p^o \psi(x,0)
+ \int_{\Omega \times J} q  \psi  .
\end{aligned} \label{2.35}
\end{equation}

\noindent(ii) The concentration part
$(f_1^\epsilon,c_1^{\epsilon},C_1^{\epsilon},c_2^{\epsilon})$ is such that
$( f_1^\epsilon,c_1^{\epsilon},C_1^{\epsilon}) \in L^2(J;V^\epsilon)
\cap H^1(J;(V^\epsilon)')$ and
 $  c_2^{\epsilon} \in L^2(J;H^1(\Omega_m^\epsilon))
\cap H^1(J;(H^1(\Omega_m^\epsilon))')$.
It satisfies for any $(d_f,d_1,D_1)  \in L^2(J;V^\epsilon)$ and any
 $d_2  \in L^2(J;H^1(\Omega_m^\epsilon))$ the following relations.
\begin{equation}
\begin{aligned}
& \int_{\Omega_f^\epsilon \times J} \phi_f^\epsilon \partial_t f_1^\epsilon\,d_f
  + \int_{\Omega_m^\epsilon \times J} \phi^\epsilon \partial_t c_1^\epsilon\,d_1
  + \int_{\Omega_m^\epsilon \times J} \phi^\epsilon \partial_t C_1^\epsilon
  D_1
  + \int_{\Omega_f^\epsilon \times J}( \underline{v}_f^\epsilon \cdot
 \nabla f_1^\epsilon)\, d_f
 \\
& + \int_{\Omega_m^\epsilon \times J} \underline{\mathcal{V}}^\epsilon \cdot ( d_1 \nabla c_1^\epsilon
  + D_1 \nabla C_1^\epsilon)
 + \int_{\Omega_f^\epsilon \times J} \mathcal{D}( \underline{v}_f^\epsilon ) \nabla f_1^\epsilon \cdot \nabla d_f
 \\
& + \int_{\Omega_m^\epsilon \times J} \mathcal{D}^\epsilon( \underline{\mathcal{V}}^\epsilon ) \nabla c_1^\epsilon \cdot \nabla d_1
 + \int_{\Omega_m^\epsilon \times J} \mathcal{D}^\epsilon( \underline{\mathcal{V}}^\epsilon ) \nabla C_1^\epsilon \cdot \nabla D_1
 \\
& = \int_{\Omega_f^\epsilon \times J} q \, (\hat f_1-f_1^\epsilon)\, d_f
 + \int_{\Omega_m^\epsilon \times J} q \, (\hat c_1-c_1^\epsilon)\, d_1
 + \int_{\Omega_m^\epsilon \times J} q \, (\hat C_1-C_1^\epsilon)\, D_1 ,
\end{aligned} \label{2.36}
\end{equation}
and
\begin{equation}
\begin{aligned}
&\int_{\Omega_m^\epsilon \times J} \phi^\epsilon \partial_t c_2^\epsilon
\, d_2
+ \int_{\Omega_m^\epsilon \times J}( \underline{\mathcal{V}}^\epsilon
\cdot \nabla c_2^\epsilon)\, d_2  \\
& + \int_{\Omega_m^\epsilon \times J} \mathcal{D}^
\epsilon( \underline{\mathcal{V}}^\epsilon ) \nabla c_2^\epsilon
\cdot \nabla d_2
 - \int_{\partial \Omega_m^\epsilon \times J} (\mathcal{D}^\epsilon
( \underline{\mathcal{V}}^\epsilon ) \nabla c_2^\epsilon \cdot \nu_m)\,
 \gamma_m^\epsilon d_2 \\
&=  \int_{\Omega_m^\epsilon \times J} q \, (1-c_2^\epsilon)\, d_2  .
\end{aligned} \label{2.37}
\end{equation}
Furthermore, the following maximum principles hold:
\begin{gather}
 0 \le f_1^\epsilon(x,t) \le \hat f_1 \quad \text{a.e.  in }
  \Omega_f^\epsilon \times J, \nonumber \\
 0 \le m_1^\epsilon(x,t) \le \hat f_1 \quad  \text{a.e.  in }
   \Omega_m^\epsilon \times J,  \label{2.38} \\
 0 \le C_\alpha \le c_1^\epsilon(x,t) + c_2^\epsilon(x,t) \le 1\quad
  \text{a.e.  in }  \Omega_m^\epsilon \times J.  \label{2.39}
\end{gather}
\end{theorem}

\begin{proof}
The proof of this existence result follows the lines of
\cite{Cho3} and is  based on a fixed point approach. The necessary
a priori estimates are the same as the uniform ones derived in
Section 3 below. We thus do not detail the proof in the present
paper. We only give the details for the crucial maximum principles
(\ref{2.38})-(\ref{2.39}). In view of the method developed in
\cite{Cho3}, we can assume in what follows that the Darcy
velocities are regularized so that $\underline{v}_f^\epsilon$ and
$\underline{\mathcal{V}}^\epsilon$ belong to $L^\infty(\Omega
\times J)$. We first consider the problem of the left-hand sides
of estimates (\ref{2.38}). We note that the function $\alpha
c_1^\epsilon + \beta C_1^\epsilon=m_1^\epsilon$ satisfies the
following system in $\Omega_m^\epsilon \times J$.
 \begin{gather}
  \phi^\epsilon \partial_t m_1^\epsilon
 + \underline{\mathcal{V}}^\epsilon \cdot \nabla m_1^\epsilon
 - \mathop{\rm div}(\mathcal{D}^\epsilon (\underline{\mathcal{V}}^\epsilon)
  \nabla m_1^\epsilon )
 = q (\hat f_1 - m_1^\epsilon ) ,  \label{2.40}
\\
 m_1^\epsilon = f_1^\epsilon, \quad
 \mathcal{D}^\epsilon (\underline{\mathcal{V}}^\epsilon)
\nabla m_1^\epsilon \cdot \nu_{fm}
 = (\alpha^2+\beta^2) \mathcal{D}(\underline{v}_f^\epsilon)
  \nabla f_1^\epsilon \cdot \nu_{fm}\quad \text{on } \Gamma_{fm}^\epsilon \times J,
  \label{2.41}
\\
\mathcal{D}^\epsilon (\underline{\mathcal{V}}^\epsilon) \nabla
m_1^\epsilon \cdot \nu = 0 \quad  \text{on } (\partial
\Omega_{m}^\epsilon  \cap \Gamma ) \times J,
  \label{2.42}
\\
m_1^\epsilon(x,0)=\alpha  c_1^o(x) + \beta C_1^o(x)
= \chi_m^\epsilon(x) f_1^o(x) \quad \text{in }Ê\;  \Omega_m^\epsilon.
   \label{2.43}
\end{gather}
For any function $f$, we denote by $(f)^-$ the function
$(f)^-=\sup(0,-f)$. We now multiply  \eqref{2.15} by
$(\alpha^2+\beta^2)(f_1^\epsilon)^-$ (respectively
(\ref{2.40}) by $(m_1^\epsilon)^-$) and we integrate over
$\Omega_f^\epsilon$ (respectively $\Omega_m^\epsilon$). Noting
that $(f_1^\epsilon)^-=(m_1^\epsilon)^-$ on
$\Gamma_{fm}^\epsilon$, we sum up the resulting relations to kill
the non zero boundary terms. We get
\begin{equation}
\begin{aligned}
& \frac{d}{2dt} \int_{\Omega} \bigl(
\chi_f^\epsilon (\alpha^2+\beta^2) \phi_f^\epsilon \vert (f_1^\epsilon)^-
 \vert^2 + \chi_m^\epsilon \phi^\epsilon \vert (m_1^\epsilon )^-\vert^2
\bigr)
 - \int_{\Omega_m^\epsilon}  (\underline{\mathcal{V}}^\epsilon \cdot
\nabla m_1^\epsilon ) \, (m_1^\epsilon)^-
 \\
&- (\alpha^2+\beta^2)  \int_{\Omega_f^\epsilon} ( \underline{v}_f^\epsilon
\cdot \nabla f_1^\epsilon)\, (f_1^\epsilon)^-
-  (\alpha^2+\beta^2) \int_{\Omega_f^\epsilon} \mathcal{D}
(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon \cdot \nabla (f_1^\epsilon)^-
 \\
&- \int_{\Omega_m^\epsilon}  \mathcal{D}^\epsilon
 (\underline{\mathcal{V}}^\epsilon) \nabla m_1^\epsilon \cdot
\nabla (m_1^\epsilon )^-
+  \int_{\Omega} q \hat f_1  \bigl(
\chi_f^\epsilon (\alpha^2+\beta^2) (f_1^\epsilon)^-
+ \chi_m^\epsilon (m_1^\epsilon )^- \bigr)
 \\
&+ \int_{\Omega} q \, \bigl(  \chi_f^\epsilon
(\alpha^2+\beta^2) \vert (f_1^\epsilon)^- \vert^2+ \chi_m^\epsilon
\vert (m_1^\epsilon )^- \vert^2 \bigr) =0.
\end{aligned} \label{2.44}
\end{equation}
Using $-\nabla f \cdot \nabla (f)^-= \nabla (f)^- \cdot \nabla (f)^-$
and the basic properties (\ref{2.22}) of the tensors
$(\mathcal{D},\mathcal{D}^\epsilon)$, we note that
\begin{gather*}
-\int_{\Omega_f^\epsilon} \mathcal{D}(\underline{v}_f^\epsilon )
\nabla f_1^\epsilon \cdot \nabla (f_1^\epsilon)^- \, dx
 \ge \int_{\Omega_f^\epsilon} \phi_- (D_m + \alpha_t \vert
\underline{v}_f^\epsilon\vert  ) \, \vert \nabla (f_1^\epsilon)^-\vert^2
\, dx ,
\\
\begin{aligned}
&- \int_{\Omega_m^\epsilon}  \mathcal{D}^\epsilon
(\underline{\mathcal{V}}^\epsilon) \nabla m_1^\epsilon \cdot
\nabla (m_1^\epsilon )^- \, dx \\
& \ge \int_{\Omega_m^\epsilon} \phi_-
(D_m(\alpha + \beta \epsilon^2) + \alpha_t \vert
\underline{\mathcal{V}}_s^\epsilon + \epsilon^2
\underline{\mathcal{V}}_h^\epsilon \vert ) \, \vert \nabla
(m_1^\epsilon )^-\vert^2 \, dx.
\end{aligned}
\end{gather*}
The convective terms in (\ref{2.44}) are then estimated as follows using
the Cauchy-Schwarz and Young inequalities.
\begin{gather*}
 \Bigl\vert  \int_{\Omega_f^\epsilon} ( \underline{v}_f^\epsilon  \cdot
\nabla f_1^\epsilon)\, (f_1^\epsilon)^-  \Bigr\vert
\le
\int_{\Omega_f^\epsilon}  \frac{ \alpha_t}{2} \vert
\underline{v}_f^\epsilon\vert  \, \vert \nabla (f_1^\epsilon)^-\vert^2
+ C \Vert \underline{v}_f^\epsilon\Vert_\infty \int_{\Omega_f^\epsilon}
\vert (f_1^\epsilon)^-\vert^2   ,
\\
\Bigl\vert  \int_{\Omega_m^\epsilon} (\underline{\mathcal{V}}^\epsilon \cdot
 \nabla m_1^\epsilon ) \, (m_1^\epsilon )^- \Bigr\vert
\le
\int_{\Omega_m^\epsilon}  \frac{ \alpha D_m}{2} \, \vert
\nabla (m_1^\epsilon )^- \vert^2
+ \frac{C}{\alpha} \Vert  \underline{\mathcal{V}}^\epsilon \Vert_\infty^2
\int_{\Omega_m^\epsilon}  \vert (m_1^\epsilon )^-\vert^2  .
\end{gather*}
All the terms  of (\ref{2.44}) containing the source $q$ are nonnegative.
The previous estimates then lead to
\begin{align*}
& \frac{\phi_-}{2} \frac{d}{dt} \int_{\Omega}
\bigl( (\alpha^2+\beta^2) \chi_m^\epsilon \vert (f_1^\epsilon)^- \vert^2
+ \chi_m^\epsilon \vert (m_1^\epsilon)^-\vert^2 \bigr) \, dx
+ (\alpha^2+\beta^2) \phi_-  \int_{\Omega_f^\epsilon}  (D_m
\\
&+\frac{ \alpha_t}{2}  \vert \underline{v}_f^\epsilon\vert  ) \,
\vert \nabla (f_1^\epsilon)^-\vert^2   dx
 + \phi_-  \int_{\Omega_m^\epsilon} (D_m(\frac{\alpha}{2}
+ \beta \epsilon^2) + \alpha_t \vert \underline{\mathcal{V}}_s^\epsilon
+\epsilon^2 \underline{\mathcal{V}}_h^\epsilon \vert )
\, \vert \nabla  (m_1^\epsilon )^-\vert^2  dx
\\
&\le  C \int_{\Omega} \bigl( (\alpha^2+\beta^2) \chi_f^\epsilon
\vert (f_1^\epsilon)^-\vert^2 + \chi_m^\epsilon
\vert (m_1^\epsilon)^-\vert^2 \bigr) \, dx .
\end{align*}
Using the Gronwall lemma with Hypotheses (\ref{2.33})-(\ref{2.34})
on the initial data, we conclude that $(f_1^\epsilon)^-(x,t)=0$
a.e. in $\Omega_f^\epsilon \times J$ and $(m_1^\epsilon
)^-(x,t)=0$ a.e. in $\Omega_m^\epsilon \times J$, that is the
first part of (\ref{2.38}). The second part is proved similarly,
multiplying   \eqref{2.15} by $(\alpha^2+\beta^2)(\hat
f_1-f_1^\epsilon)^-$ (respectively   (\ref{2.40}) by $(\hat
f_1-m_1^\epsilon )^-$).

We now justify  (\ref{2.39}).
To this aim, we consider the problem satisfied by
$(c_1^\epsilon+c_2^\epsilon)$ in $\Omega_m^\epsilon \times J$:
\begin{gather}
 \phi^\epsilon \partial_t (c_1^\epsilon+c_2^\epsilon)
 +  \underline{\mathcal{V}}^\epsilon \cdot \nabla
 (c_1^\epsilon+c_2^\epsilon)
 -\mathop{\rm div}(\mathcal{D}^\epsilon \nabla  (c_1^\epsilon+c_2^\epsilon) )
 =
 q (1-c_1^\epsilon-c_2^\epsilon), \label{2.45} \\
\mathcal{D}^\epsilon( \underline{\mathcal{V}}^\epsilon)
\nabla  (c_1^\epsilon+c_2^\epsilon) \cdot \nu =0 \quad
\text{on }  \partial \Omega_m^\epsilon \times J ,
  \label{2.46} \\
(c_1^\epsilon+c_2^\epsilon)(x,0)=c_1^o(x)+c_2^o(x) \quad \text{in }  \Omega.
 \label{2.47}
 \end{gather}
 Multiplying (\ref{2.45}) by $(c_1^\epsilon+c_2^\epsilon-C_\alpha)^-$,
integrating over $\Omega_m^\epsilon$ and using the same tools as in the
 previous part of this proof, we state that
$c_1^\epsilon(x,t)+c_2^\epsilon(x,t) \ge C_\alpha \ge 0$ almost everywhere
in $\Omega_m^\epsilon \times J$.
 We now detail the proof of the second part of (\ref{2.39}).
 We multiply (\ref{2.45}) by $(1-c_1^\epsilon-c_2^\epsilon)^-$ and we
integrate over $\Omega_m^\epsilon$.
 Using the basic properties (\ref{2.22}) of the tensor $\mathcal{D}^\epsilon$
we get
 \begin{align*}
& \frac{\phi_- }{2} \frac{d}{dt} \int_{\Omega_m^\epsilon}  \vert (1-  c_1^\epsilon -c_2^\epsilon )^-\vert^2
+ \phi_- \int_{\Omega_m^\epsilon}   D_m(\alpha + \beta \epsilon^2) \vert \nabla (  1-c_1^\epsilon - c_2^\epsilon )^- \vert^2
 \\
&+ \int_{\Omega_m^\epsilon}  (\underline{\mathcal{V}}^\epsilon \cdot  \nabla (  c_1^\epsilon +c_2^\epsilon ) )\, (1 - c_1^\epsilon - c_2^\epsilon )^-
- \int_{\Omega_m^\epsilon} q  (1-  c_1^\epsilon -c_2^\epsilon ) \,
(1 -  c_1^\epsilon - c_2^\epsilon )^-
\le 0 .
\end{align*}
 The fourth term of the left-hand side is nonnegative.
 The convective term is estimated with the Cauchy-Schwarz and Young
inequalities as follows.
 \begin{align*}
 &  \Bigl\vert \int_{\Omega_m^\epsilon}  (\underline{\mathcal{V}}^\epsilon \cdot  \nabla (  c_1^\epsilon +c_2^\epsilon ) )\, (1 - c_1^\epsilon - c_2^\epsilon )^-  \Bigr\vert
 \\
 &  = \Bigl\vert \int_{\Omega_m^\epsilon}  (\underline{\mathcal{V}}^\epsilon \cdot  \nabla ( 1- c_1^\epsilon -c_2^\epsilon )^- )\, (1 - c_1^\epsilon - c_2^\epsilon )^-  \Bigr\vert
 \\
 & \le  \int_{\Omega_m^\epsilon}  \frac{\phi_- \alpha D_m}{2}\,
 \vert \nabla ( 1- c_1^\epsilon -c_2^\epsilon )^- \vert^2
+ C \Vert \underline{\mathcal{V}}^\epsilon \Vert_\infty^2
\int_{\Omega_m^\epsilon}  \vert (1 - c_1^\epsilon - c_2^\epsilon )^-
 \vert^2.
\end{align*}
 Using these estimates in the first relation, we get
\begin{align*}
 & \frac{\phi_- }{2}  \frac{d}{dt} \int_{\Omega_m^\epsilon}
 \vert (1-  c_1^\epsilon -c_2^\epsilon )^-\vert^2 \, dx
+ \phi_- \int_{\Omega_m^\epsilon}   D_m(\frac{\alpha}{2}
+ \beta \epsilon^2)\vert \nabla (  1-c_1^\epsilon - c_2^\epsilon )^-
\vert^2 \, dx
 \\
& \le C \int_{\Omega_m^\epsilon} \vert (1 -
c_1^\epsilon - c_2^\epsilon )^- \, dx \vert^2\, dx.
\end{align*}
The latter relation and the Gronwall lemma combined with Assumption
(\ref{2.34}) let us conclude that
$ (1 - c_1^\epsilon - c_2^\epsilon )^-(x,t)=0$ and thus
$c_1^\epsilon(x,t)+c_2^\epsilon(x,t) \le 1 $ a.e. in
$\Omega_m^\epsilon \times J$.
 This completes the proof of the theorem.
 \end{proof}

\section{Uniform estimates}

We begin by stating the following properties of the pressure solutions of
the problem \eqref{2.16}, \eqref{2.20}-\eqref{2.21}, \eqref{2.27},
(\ref{2.30}), \eqref{2.32}.


\begin{lemma} \label{lem1}
The pressure satisfies the following uniform estimates
\begin{gather*}
 \Vert p_f^{\epsilon}\Vert_{L^\infty(J;L^2(\Omega_f^\epsilon ))}
 + \Vert p_f^{\epsilon}\Vert_{L^2(J;H^1(\Omega_f^\epsilon ))} \le C,
 \\
 \Vert \underline{v}_f^\epsilon \Vert_{(L^2(J;L^2(\Omega_f^\epsilon)))^2}
\le C,
 \\
 \Vert p^{\epsilon}\Vert_{L^\infty(J;L^2(\Omega_m^\epsilon ))}
 \le C,
\\
\Vert \alpha^{1/2} (c_1^\epsilon+c_2^\epsilon)^{1/2} \nabla p^\epsilon \Vert_{(L^2(J;L^2(\Omega_m^\epsilon)))^2}
 + \Vert \epsilon  \nabla p^\epsilon \Vert_{(L^2(J;L^2(\Omega_m^\epsilon)))^2}
 \le C,
 \\
 \Vert \underline{\mathcal{V}}_s^\epsilon
\Vert_{(L^2(J;L^2(\Omega_m^\epsilon)))^2} \le C,
\\\
\Vert \underline{\mathcal{V}}_h^\epsilon
\Vert_{(L^2(J;L^2(\Omega_m^\epsilon)))^2} \le C.
\end{gather*}
Furthermore the time derivative
$(\chi_f^\epsilon \phi_f^\epsilon \partial_t p_f^\epsilon
+ \chi_m^\epsilon \phi^\epsilon \partial_t p^\epsilon)$ is
uniformly bounded in $L^2(J;(H^1(\Omega))')$.
\end{lemma}



\begin{proof}
The estimates are derived from integration by parts.
We multiply   \eqref{2.16} by $p_f^{\epsilon}$ and integrate over
$\Omega_f^\epsilon \times J$.
We multiply   \eqref{2.20} by $p^{\epsilon}$ and integrate over
$\Omega_m^\epsilon \times J$.
 Summing up the  resulting relations, we obtain
\begin{align*}
& \frac{1}{2} \int_{\Omega_f^\epsilon}\phi_f^\epsilon\,
  \vert p_f^{\epsilon}\vert^2 \,  dx
+ \frac{1}{2} \int_{\Omega_m^\epsilon}\phi^{\epsilon}\,
  \vert p^{\epsilon}\vert^2 \,  dx
+\int_{\Omega_f^\epsilon \times J} \frac{k_f^\epsilon}{\mu(f_1^\epsilon)}
\nabla p_f^\epsilon \cdot \nabla p_f^\epsilon \, dx\,dt
\\
&+\int_{\Omega_m^\epsilon \times J}  \bigl( \alpha (c_1^\epsilon+c_2^\epsilon)
 (1-\epsilon^2) + \epsilon^2 \bigr)\, \frac{k^\epsilon}{\mu(m_1^\epsilon)}
  \nabla p^\epsilon \cdot \nabla p^\epsilon \, dx\,dt
\\
&= \frac{1}{2} \int_{\Omega} ( \chi_f^\epsilon \phi_f^\epsilon(x)
 + \chi_m^\epsilon \phi^{\epsilon}(x))\,  \vert p^{o}(x)\vert^2  \,  dx
 +\int_{\Omega \times J} q \, (\chi_f^\epsilon p_f^{\epsilon}
 + \chi_m^\epsilon p^{\epsilon})\,  dx\,dt.
\end{align*}
Applying the Cauchy-Schwarz and Young inequalities with the properties
of $\phi_f^\epsilon$, $\phi^\epsilon$, $k_f^\epsilon$, $k^\epsilon$ and
 $\mu$ in the latter relation, we get
\begin{align*}
& \frac{\phi_-}{2} \int_{\Omega_f^\epsilon} \vert p_f^{\epsilon}\vert^2 \,  dx
+ \frac{\phi_-}{2} \int_{\Omega_m^\epsilon}  \vert p^{\epsilon}\vert^2 \,  dx
+ \frac{k_-}{\mu_+} \int_{\Omega_f^\epsilon \times J}  \vert \nabla p_f^\epsilon \vert^2 \, dx\,dt
\\
&+ \frac{k_-}{\mu_+}
 \int_{\Omega_m^\epsilon \times J} \bigl(
  \alpha (c_1^\epsilon+c_2^\epsilon)\, \vert \nabla p^\epsilon \vert^2
+   \epsilon^2 \, (1- \alpha (c_1^\epsilon+c_2^\epsilon))\,  \vert \nabla p^\epsilon \vert^2
\bigr)\, dx\,dt
\\
&\le C\bigl( \Vert p^o \Vert_{L^2(\Omega)}, \Vert q \Vert_{L^2
 (\Omega \times J)} \bigr)
+ \int_{\Omega_f^\epsilon \times J} \vert p_f^{\epsilon}\vert^2 \,  dx\,dt
+ \int_{\Omega_m^\epsilon \times J} \vert p^{\epsilon}\vert^2 \,  dx\,dt.
\end{align*}
Using the Gronwall lemma, we prove the desired estimates. The
result on the time derivatives then follows straightforward  from
 \eqref{2.16}, \eqref{2.20}-\eqref{2.21}.
\end{proof}

We can now establish the following results concerning the concentrations
functions $(f_1^{\epsilon}, C_1^\epsilon,c_1^\epsilon,c_2^\epsilon)$.

\begin{lemma} \label{lem2}
\begin{itemize}
\item[(i)] The functions $(f_1^{\epsilon},C_1^\epsilon,
c_1^\epsilon,c_2^\epsilon)$
are uniformly bounded in the space
$L^{\infty}(J;L^2(\Omega_f^\epsilon)) \times (L^{\infty}
(J;L^2(\Omega_f^\epsilon)))^3$ and are such that
\begin{gather*}
 0 \le f_1^{\epsilon}(x,t) \le \hat f_1 \le 1 \quad
\text{ almost   everywhere   in }   \Omega_f^\epsilon \times J ,\\
 0 \le \alpha c_1^{\epsilon}(x,t) + \beta C_1^\epsilon(x,t) \le \hat f_1
\le 1 \quad \text{ almost   everywhere  in   }     \Omega_m^\epsilon \times J \\
 0 \le  c_1^{\epsilon}(x,t) + c_2^\epsilon(x,t) \le 1 \quad
\text{ almost   everywhere  in   }     \Omega_m^\epsilon \times
J ;
\end{gather*}
\item[(ii)] the sequence
 $( (D_m^{1/2}+\alpha_t^{1/2} \vert \underline{v}_f^\epsilon \vert^{1/2})
 \nabla  f_1^{\epsilon})$ is uniformly bounded in
$(L^2(\Omega_f^\epsilon \times J ))^2$;

\item[(iii)] for $i=1,2$, the diffusive terms
$\alpha^{1/2} (1+ (c_1^\epsilon+c_2^\epsilon)^{1/2} \vert \nabla p^\epsilon
\vert^{1/2} )\nabla c_i^\epsilon$ and
$\epsilon ( 1 +  \vert \epsilon \nabla p^\epsilon \vert^{1/2} )
\nabla c_i^\epsilon $
are uniformly  bounded in $(L^2(\Omega_m^\epsilon \times J ))^2$.
The same estimates hold for $C_1^\epsilon$.
\end{itemize}
\end{lemma}


\begin{proof}
The maximum principles of {\it (i)} are a direct consequence of
the construction of the solution $(f_1^{\epsilon},
C_1^\epsilon,c_1^\epsilon,c_2^\epsilon)$ in Theorem \ref{thm1}. We
write the variational formulation (\ref{2.36}) with the test
function $(d_f,d_1,D_1)=(f_1^\epsilon,c_1^\epsilon,C_1^\epsilon)$.
We get
\begin{equation}
\begin{aligned}
&\frac{1}{2} \int_{\Omega_f^\epsilon} \phi_f^\epsilon \vert f_1^\epsilon \vert^2 \, dx
+  \frac{1}{2} \int_{\Omega_m^\epsilon} \phi^\epsilon ( \vert   c_1^\epsilon \vert^2 + \vert  C_1^\epsilon \vert^2 ) \, dx
+ \int_{\Omega_f^\epsilon \times J} \mathcal{D}(\underline{v}_f^\epsilon ) \nabla f_1^\epsilon \cdot \nabla f_1^\epsilon \, dx\,dt
 \\
&+ \int_{\Omega_m^\epsilon \times J}  (\mathcal{D}^\epsilon (\underline{\mathcal{V}}^\epsilon) \nabla  c_1^\epsilon  \cdot  \nabla   c_1^\epsilon +  \mathcal{D}^\epsilon (\underline{\mathcal{V}}^\epsilon) \nabla  C_1^\epsilon  \cdot  \nabla   C_1^\epsilon )\, dx\,dt
 \\
&+  \int_{\Omega_f^\epsilon \times J} ( \underline{v}_f^\epsilon  \cdot \nabla f_1^\epsilon)\, f_1^\epsilon \, dx\,dt
+ \int_{\Omega_m^\epsilon \times J} \underline{\mathcal{V}}^\epsilon \cdot ( c_1^\epsilon\,  \nabla   c_1^\epsilon  +   C_1^\epsilon\, \nabla   C_1^\epsilon )  \, dx\,dt
 \\
&+  \int_{\Omega \times J} q \,(\chi_f^\epsilon  \vert f^\epsilon_1\vert^2
 + \chi_m^\epsilon\,  (\vert c^\epsilon_1\vert^2
 + \vert C^\epsilon_1\vert^2 ) ) \, dx\,dt\\
&=
\int_{\Omega  \times J} q  \hat f_1 f^\epsilon_1 \, dx\,dt
+ \int_{\Omega_m^\epsilon  \times J} q  (\hat c_1 c^\epsilon_1
 +  \hat C_1 C^\epsilon_1  ) \, dx\, dt\\
&\quad + \frac{1}{2} \int_{\Omega} \bigl(  \phi_f^\epsilon \vert f_1^o \vert^2
+\phi^\epsilon ( \vert   c_1^o \vert^2 + \vert  C_1^o \vert^2 )\bigr) dx.
\end{aligned} \label{3.1}
\end{equation}
The convective terms in (\ref{3.1}) are estimated as follows using the
Cauchy-Schwarz and Young inequalities.
In the fractured part, we write
\[
\Bigl\vert  \int_{\Omega_f^\epsilon \times J} ( \underline{v}_f^\epsilon
 \cdot \nabla f_1^\epsilon)\, f_1^\epsilon \, dx \Bigr\vert
\le
\int_{\Omega_f^\epsilon \times J}  \frac{ \alpha_t}{2} \vert \underline{v}_f^\epsilon\vert  \, \vert \nabla f_1^\epsilon \vert^2  \, dx
+ C \Vert f_1^\epsilon \Vert_\infty^2 \int_{\Omega_f^\epsilon}  \vert \underline{v}_f^\epsilon \vert  \, dx,
\]
where $0 \le f_1^\epsilon(x,t) \le 1$ a.e. in $\Omega_f^\epsilon
\times J$  and $\underline{v}_f^\epsilon$ is uniformly bounded in
$(L^1(\Omega_f^\epsilon \times J))^2$ thanks to Lemma \ref{lem1}.
In the matrix part, the work is more difficult because we do not
have an estimate for $c_1^\epsilon$ and $C_1^\epsilon$ in
$L^\infty(\Omega_m^\epsilon \times J)$. We thus get firstly
\begin{align*}
 \Bigl\vert \int_{\Omega_m^\epsilon \times J} \underline{\mathcal{V}}^\epsilon
 \cdot ( c_1^\epsilon\,  \nabla   c_1^\epsilon
 +   C_1^\epsilon\, \nabla   C_1^\epsilon )  \Bigr\vert
&\le
\int_{\Omega_m^\epsilon \times J}  \frac{ \alpha_t}{2} \vert
 \underline{\mathcal{V}}_s^\epsilon +\epsilon^2
 \underline{\mathcal{V}}_h^\epsilon  \vert \, (\vert  \nabla
  c_1^\epsilon \vert^2 + \vert \nabla C_1^\epsilon \vert^2 )
\\
&\quad + C \int_{\Omega_f^\epsilon}( \vert
  \underline{\mathcal{V}}_s^\epsilon
  \vert +\vert \underline{\mathcal{V}}_h^\epsilon \vert) \, (\vert
c_1^\epsilon \vert^2+ \vert C_1^\epsilon \vert^2 )\, dx .
\end{align*}
In such a matrix domain, the Gagliardo-Nirenberg inequality reads
$$
\Vert u \Vert_{L^4(\Omega_m^\epsilon)} \le C \Vert u
\Vert_{L^2(\Omega_m^\epsilon)}^{1/2} \Vert (\alpha
+ \beta \epsilon )u \Vert_{H^1(\Omega_m^\epsilon)}^{1/2},
\quad \forall u \in H^1(\Omega_m^\epsilon).
$$
The second term of the right-hand side of the latter relation is
then treated as follows using the Gagliardo-Nirenberg inequality
and Lemma \ref{lem1}.
\begin{align*}
&\int_{\Omega_f^\epsilon} \vert  \underline{\mathcal{V}}^\epsilon \vert \,
(\vert  c_1^\epsilon \vert^2+ \vert C_1^\epsilon \vert^2 ) \\
&\le \frac{k_+}{\mu_-}
\int_{\Omega_f^\epsilon} \bigl( \alpha(c_1^\epsilon+c_2^\epsilon)
(1-\epsilon)+ \epsilon \bigr) \vert \nabla p^\epsilon\vert \,
\bigl( \vert c_1^\epsilon\vert^2 + \vert c_2^\epsilon\vert^2 \bigr)
\\
&\le C \Bigl(  \int_{\Omega_f^\epsilon}
 \bigl( \alpha^2(c_1^\epsilon+c_2^\epsilon)^2
 + \epsilon^2 (1-\alpha(c_1^\epsilon+c_2^\epsilon) \bigr)^2
  \vert \nabla p^\epsilon\vert^2  \Bigr)^{1/2}
\\
&\quad \times  \bigl( \Vert c_1^\epsilon \Vert_{L^4(\Omega_m^\epsilon )}^2
  + \Vert C_1^\epsilon \Vert_{L^4(\Omega_m^\epsilon )}^2  \bigr) \\
&\le C  \bigl( \Vert c_1^\epsilon \Vert_{L^4(\Omega_m^\epsilon )}^2
 + \Vert C_1^\epsilon \Vert_{L^4(\Omega_m^\epsilon )}^2  \bigr)
\\
&\le C  \bigl( \Vert c_1^\epsilon \Vert_{L^2(\Omega_m^\epsilon )}
\Vert (\alpha + \beta \epsilon)c_1^\epsilon \Vert_{H^1(\Omega_m^\epsilon)}
+ \Vert C_1^\epsilon \Vert_{L^2(\Omega_m^\epsilon )}
 \Vert (\alpha + \beta \epsilon)C_1^\epsilon \Vert_{H^1(\Omega_m^\epsilon)}
 \bigr)
\\
& \le \frac{C}{\delta} \bigl( \Vert c_1^\epsilon
 \Vert_{L^2(\Omega_m^\epsilon )}^2 + \Vert C_1^\epsilon
 \Vert_{L^2(\Omega_m^\epsilon )}^2 \bigr)
+ \phi_-\delta  \int_{\Omega_f^\epsilon} (\alpha + \epsilon^2)
 D_m \bigl( \vert \nabla c_1^\epsilon \vert^2
 + \vert \nabla C_1^\epsilon \vert^2 \bigr)  ,
\end{align*}
for any $\delta >0$.
The  last term in the left-hand side of (\ref{3.1}) is nonnegative.
Using the latter estimates, the Cauchy-Schwarz and Young inequalities
for the right-hand side source terms  and the basic properties (\ref{2.22})
 of the tensors $\mathcal{D}$ and $\mathcal{D}^\epsilon$, it follows
from (\ref{3.1}) that
\begin{align*}
&\frac{\phi_-}{2} \int_{\Omega}(  \chi_f^\epsilon \vert f_1^\epsilon \vert^2
 + \chi_m^\epsilon  (\vert   c_1^\epsilon \vert^2+ \vert C_1^\epsilon \vert^2))
 \, dx
+ \phi_- \int_{\Omega_f^\epsilon  \times J}  (D_m +\frac{ \alpha_t}{2}
 \vert \underline{v}_f^\epsilon\vert  ) \, \vert \nabla f_1^\epsilon \vert^2
  \, dx\,dt
\\
&+ \phi_- \int_{\Omega_m^\epsilon  \times J}  ((\alpha + \epsilon^2 )(1-\delta)D_m +\frac{ \alpha_t}{2}  \vert  \underline{\mathcal{V}}_s^\epsilon + \epsilon^2 \underline{\mathcal{V}}_h^\epsilon \vert  ) \, \bigl( \vert \nabla c_1^\epsilon \vert^2  + \vert C_1^\epsilon \vert^2 \bigr) \, dx\,dt
\\
&\le C + C  \int_{\Omega_f^\epsilon \times J}  \vert f_1^\epsilon \vert^2
 \, dx\,dt
 +  \frac{C}{\delta} \int_{\Omega_m^\epsilon  \times J}
( \vert c_1^\epsilon \vert^2  +  \vert  C_1^\epsilon \vert^2 )\, dx\,dt.
\end{align*}
We choose $0 < \delta < 1$.
Using the Gronwall lemma yields to the result for $f_1^\epsilon$,
$c_1^\epsilon$ and $C_1^\epsilon$.
Once we know the estimate for $c_1^\epsilon$, we obtain similar ones
for $c_2^\epsilon$ by multiplying   (\ref{2.18}) by $c_1^\epsilon$,
(\ref{2.19}) by $c_2^\epsilon$, integrating over $\Omega_m^\epsilon$
and summing up the results to kill the terms on $\Gamma_{fm}$.
 Our claim is proved.
\end{proof}


\begin{remark} \label{rmk1} \rm
This proof fails if we consider a domain $\Omega \subset \mathbb{R}^3$ instead of
$\mathbb{R}^2$ because of the use of the Gagliardo-Nirenberg inequality.
Nevertheless we can get the same result using minor modifications.
The simplest way is to add a $L^\infty$-estimate for $c_1^\epsilon$ and
$C_1^\epsilon$.
To this aim, we note that the function
$d_1^\epsilon= \frac{\alpha}{\beta}C_1^\epsilon-c_1^\epsilon$ is solution
of the following problem in $\Omega_m^\epsilon \times J$.
\begin{gather*}
\phi^\epsilon \partial_t d_1^\epsilon +
\underline{\mathcal{V}}^\epsilon \cdot \nabla d_1^\epsilon
-\mathop{\rm div}(\mathcal{D}^\epsilon(
\underline{\mathcal{V}}^\epsilon ) \nabla d_1^\epsilon )
= q (\frac{\alpha}{\beta}\hat C_1 -\hat c_1 -d_1^\epsilon),\\
 \mathcal{D}^\epsilon( \underline{\mathcal{V}}^\epsilon ) \nabla
d_1^\epsilon \cdot \nu =0 \quad  \text{on }   \partial
\Omega_m^\epsilon \times J, \\
 d_1^\epsilon(x,0)= \frac{\alpha}{\beta}C_1^o(x)-c_1^o(x)
 \quad \text{in } \Omega_m^\epsilon.
\end{gather*}
Similar arguments to the ones used in Theorem \ref{thm1} together
with an additional assumption on $\frac{\alpha}{\beta}C_1^o-c_1^o$
lead to a maximum principle for
$\frac{\alpha}{\beta}C_1^\epsilon-c_1^\epsilon$. Combining it with
Lemma \ref{lem2} (ii), we get a $L^\infty$-bound for
$c_1^\epsilon$ and $C_1^\epsilon$.
\end{remark}

\section{The macroscopic model: the case of a totally fractured media
($\alpha=0$)}

We assume in this section that $\alpha=0$ and then $\beta=1$.
We thus do not  consider the concentrations variables
$(c_1^\epsilon,c_2^\epsilon)$.
And we define   global pressure and concentration functions
$\theta^\epsilon$ and $\xi^\epsilon$ by
\begin{equation}
\theta^\epsilon =\begin{cases}
 p_f^\epsilon & \text{in }   \Omega_f^\epsilon \times J, \\
 p^\epsilon   & \text{in }   \overline{\Omega_m^\epsilon} \times J , \\
\end{cases}
 \qquad \xi^\epsilon = \begin{cases}
 f_1^\epsilon & \text{in }   \Omega_f^\epsilon \times J, \\
 C_1^\epsilon & \text{in }   \overline{\Omega_m^\epsilon} \times J, \\
\end{cases}
\label{4.1}
\end{equation}
and the global porosity $\Phi^\epsilon$, tensor of permeability
$K^\epsilon$ and new diffusion tensor $\mathcal{D}^\epsilon $ by
$ \Phi^\epsilon = \chi_f^\epsilon  \phi_f^\epsilon
+ \chi_m^\epsilon \phi^\epsilon$,
$K^\epsilon = \chi_f^\epsilon   k_f^\epsilon  +\chi_m^\epsilon  \epsilon^2
 k^\epsilon$ and
 $ \mathcal{D}^\epsilon = \chi_f^\epsilon  \mathcal{D}
+\chi_m^\epsilon \mathcal{D}^\epsilon$.

\subsection{The limit double porosity model}

The aim of this section is to derive rigorously  the double porosity model
described below.
We obtain a macroscopic fracture system driven by equations in
$\Omega \times J$, similar to the microscopic ones:
\begin{gather}
 \overline{\phi_f}^{Y_f}  \partial_t P_f - \mathop{\rm div}
\underline{\mathcal{V}}_f = q - \int_{Y_m}\phi(y) \,   \partial_tp
\,  dy , \quad \underline{\mathcal{V}}_f
=-\frac{\overline{K}_f}{\mu(F_1)}  \nabla P_f ,
\label{4.2} \\
\begin{aligned}
&  \overline{\phi_f}^{Y_f}   \partial_t
F_1 + \underline{\mathcal{V}}_f    \cdot \nabla F_1 - \mathop{\rm
div}( \overline{\mathcal{D}}_f(\nabla P_f,\mu(F_1))  \nabla F_1 )
+q\,  |Y_f| \,  F_1 \\
&= q\hat f_1  -q \int_{Y_m} C_1 \,  dy
- \int_{Y_m} \phi(y)\,  \partial_t C_1\,   dy
+ \int_{Y_m}  \frac{k(y)}{\mu(C_1)}  \nabla_y p \cdot \nabla_y C_1 \,  dy.
\end{aligned}\label{4.3}
\end{gather}
The matrix plays the role of a source and produces the additional
right-hand side source-like terms.
The homogenization process  gives for this macroscopic level an
effective porosity $\overline{\phi_f}^{Y_f}$, an effective rock
permeability $\overline{K}_f$ and an homogenized tensor of diffusion
$\overline{\mathcal{D}}_f$  defined by
\begin{gather}
\overline{\phi_f}^{Y_f} = \int_{Y_f} \phi_f(y)  dy ,
\label{4.4} \\
\overline{K}_{f_{ij}}=\int_{Y_f} k_f(y) (\nabla_y  v^i(y)+e^i )
\cdot (\nabla_y v^j(y)+e^j )dy
\quad 1\le i,j \le2,
\label{4.5} \\
\overline{\mathcal{D}}_{f_{ij}}(\nabla P_f,\mu(F_1))=\int_{Y_f} D_f (\nabla_y   w^i+e^i )
\cdot (\nabla_y w^j+e^j )dy,
\label{4.6}
\end{gather}
$ D_f=\mathcal{D}\Bigl( \frac{K_f(y)}{\mu(F_1)} \nabla P_f \Bigr)$
and
$K_f(y)_{ij}= k_f(y)\left(\nabla_y  v^i(y)+e^i\right) \cdot \left(\nabla_y v^j(y)+e^j\right)$.
 The functions $(v^i)_{1 \le i \le 2}$ and $(w^i)_{1 \le i \le 2}$ are
respectively solutions of the cell problems (\ref{4.7}) and (\ref{4.8}) below.
\begin{gather}
\begin{gathered}
 -\mathop{\rm div}_y ( k_f(y)
(\nabla_y v^i(y)+e^i ))=0 \quad \text{in }  Y_f,\\
 k_f(y) (\nabla_y  v^i(y)+e^i )\cdot \nu_y=0 \quad
\text{on }  \Gamma_{fm} , \; y \mapsto v^i(x,y)\;
Y\text{-periodic},
\end{gathered}
\label{4.7}
\\
\begin{gathered}
 -\mathop{\rm div}_y ( D_f(x,y) (\nabla_y
w^i(x,y)+e^i ))=0 \quad \text{in } Y_f, \\
 D_f(x,y) (\nabla_y  w^i(x,y)+e^i )\cdot \nu_y=0 \quad
\text{on }  \Gamma_{fm} , \; y \mapsto w^i(x,y)\;
Y\text{-periodic},
\end{gathered}\label{4.8}
\end{gather}
where $e^j$ is the unit vector in the $j$-th direction.
On the other hand, to each $x \in \Omega$ corresponds a matrix block,
driven by equations in $\{x\} \times Y_m \times J$ which give the new
source terms:
\begin{gather}
 \phi(y) \,  \partial_tp +\mathop{\rm div}_y
\underline{\mathcal{V}}_h = q   , \quad \underline{\mathcal{V}}_h
= - \frac{k(y)}{\mu(C_1)} \nabla_y p,
\label{4.9} \\
 \phi(y)\,  \partial_t
C_1 + \underline{\mathcal{V}}_h \cdot \nabla_y C_1 - \mathop{\rm
div}_y( \mathcal{D}(y,\underline{\mathcal{V}}_h)\,  \nabla_y C_1)
+q\,  C_1 = q \hat f_1 . \label{4.10}
\end{gather}
The equations \eqref{4.2}-(\ref{4.3}), (\ref{4.9})-(\ref{4.10}) are
provided with the following initial and boundary conditions
\begin{gather}
 \overline{K}_f   \nabla P_f \cdot \nu = 0 ,
\quad
\overline{\mathcal{D}}_f   \nabla F_1 \cdot \nu = 0
 \quad  \text{on }  \Gamma \times J ,
\quad p=P_f , \; C_1=F_1 \quad  \text{on }  \Gamma_{fm},
\label{4.11} \\
 P_f(x,0) = p(x,y,0) = p^o(x) , \quad
F_1(x,0) =  C_1(x,y,0) = f_{1}^o(x) \quad  \text{in }
  \Omega \times Y_m   . \label{4.12}
\end{gather}
We claim and prove the following convergence result.


\begin{theorem} \label{thm2}

As the scaling parameter $\epsilon$ tends to zero, the microscopic
model \eqref{2.15}--\eqref{2.17}, \eqref{2.20}--\eqref{2.21},
\eqref{2.23}, \eqref{2.26}--\eqref{2.32} with $\alpha =0$ converges
to the double porosity macroscopic model \eqref{4.2}--\eqref{4.12}.
\end{theorem}

We have captured the interactions between the local and the global
scales.
This   model is consistent with the double porosity formulations of
the engineering literature, {\it cf} \cite{Bar}.
But in \cite{Bar}, the exchanges between fractures and blocks are
assumed quasi-stationary.
Without this hypothesis, we obtain additional memory terms.

The proof of the homogenization process will be carried out by using the
 two-scale convergence  introduced by G.Nguetseng in \cite{Ng} and
developed by Allaire in \cite{A}.
We recall the basic definition and properties of this concept.

\begin{proposition} \label{prop1}
A sequence of functions $(v^{\epsilon})$ bounded in $L^2(\Omega \times J )$
two-scale converges to a limit  $v^o(x,y,t)$ belonging to
$L^2(\Omega  \times  Y \times J )$,
$v^\epsilon {\stackrel{2}{\rightharpoonup}}  v^o$,  if
$$
\lim_{\epsilon \to 0}\int_{\Omega \times J}v^{\epsilon}(x,t)\,
 \Psi (x,x/\epsilon,t) \,  dx\,dt
=\int_{\Omega \times J}\int_{Y}v^o(x,y,t)\,  \Psi(x,y,t)\,  dx\,dy\,dt,
$$
for any test function  $\Psi(x,y,t)$, Y-periodic in the second variable,
satisfying
$$
\lim_{\epsilon \to 0} \int_{\Omega \times J} \vert \Psi (x,x/\epsilon,t)
\vert^2\,  dx\,dt
= \int_{\Omega \times J} \int_Y |\Psi(x,y,t)|^2\,  dx\,dy\,dt .
$$
(i) From each bounded sequence $ (v^{\epsilon} )$ in  $L^2 (\Omega \times J )$
one can extract a subsequence which  two-scale converges.

\noindent (ii)
Let $ (v^{\epsilon} )$ be a bounded sequence in  $L^2(J;H^1 (\Omega  ))$
which converges weakly to $v$ in  $L^2(J;H^1 (\Omega  ))$.
Then  $ v^{\epsilon}   {\stackrel{2}{\rightharpoonup}} v$  and there exists
a function $v^1 \in L^2 (\Omega \times J;H_{per}^1 (Y ) )$ such that,
up to a subsequence,
$ \nabla v^{\epsilon}  {\stackrel{2}{\rightharpoonup}}  \nabla v (x,t )
+ \nabla_y v^1 (x,y,t )$.

\noindent(iii) Let $ (v^{\epsilon} )$ be a bounded sequence in
$L^2 (\Omega \times J )$ with $(\epsilon \nabla  v^{\epsilon} )$ bounded
in $(L^2 (\Omega \times J ))^2$.
 Then, there exists a function $v^o \in L^2 (\Omega \times  J;H_{per}^1 (Y ) )$
such that, up to a subsequence,
$ v^{\epsilon}   {\stackrel{2}{\rightharpoonup}} v^o$ and
$\epsilon \nabla  v^{\epsilon}  {\stackrel{2}{\rightharpoonup}}
\nabla_y v^o (x,y,t )$.
\end{proposition}

To exploit the {\it a priori} estimates obtained in the fractured part
$\Omega_f^\epsilon$,  we need to extend the functions $p_f^\epsilon$
and $f_1^\epsilon$ to the whole domain $\Omega$.
To this aim, following \cite{Acer}, we claim that there exists three
constants $k_i=k_i(Q_f) >0 $, $i=1,2,3$, and a linear and continuous
extension operator
$ \Pi^\epsilon  : H^1(\Omega_f^\epsilon) \to H_{loc}^1(\Omega) $
 such that $\Pi^\epsilon v =v$ a.e. in $\Omega_f^\epsilon$ and
$$
 \int_{\Omega (\epsilon k_1)} |\Pi^\epsilon v|^2\,  dx  \le k_2
\int_{\Omega_f^\epsilon}|v|^2\,  dx ,
 \quad
 \int_{\Omega (\epsilon k_1)} |\nabla (\Pi^\epsilon  v)|^2\,  dx
\le  k_3 \int_{\Omega_f^\epsilon}|\nabla v|^2\,  dx
$$
for all $v \in H^1(\Omega_f^\epsilon)$, with $\Omega (\epsilon
k_1) =\lbrace x \in \Omega   \mid \mathop{\rm dist}(x, \Gamma) >
\epsilon k_1 \rbrace $. To avoid dealing with boundary layers, we
make the following additional assumption on the structure of the
domain $\Omega$:
$$
\Omega_m^\epsilon =\Omega (\epsilon k_1) \cap \bigl\lbrace
  \cup_{k \in \mathbb{Z}^2} \epsilon \left(Y_m +k \right) \bigr\rbrace
\quad \text{and} \quad \Omega_f^\epsilon = \Omega \setminus
\overline{\Omega_m^\epsilon}.
$$
 For any subset $\Omega' \subset \subset \Omega$,  we get with
Lemmas \ref{lem1} and \ref{lem2} $ \int_{\Omega' \times J} |\nabla
( \Pi^\epsilon  p_f^\epsilon ) |^2  dx\,dt \le C$ and $
\int_{\Omega' \times J} |\nabla ( \Pi^\epsilon f_1^\epsilon ) |^2
dx\,dt \le C $, if $\epsilon <\mathop{\rm
dist}(\Omega',\Gamma)/k_1$. Thus we can state convergence results
in any subset $\Omega' \subset \subset \Omega$ and conclude with a
density argument. But, for sake of simplicity, we prefer assume
that the blocks are removed in an $\epsilon k_1$-neighborhood of
$\Gamma$. We then get
$$
\int_{\Omega \times J} |\nabla ( \Pi^\epsilon  p_f^\epsilon ) |^2\,  dx\,dt
\le C , \quad \int_{\Omega \times J} |\nabla ( \Pi^\epsilon
f_1^\epsilon ) |^2\,  dx\,dt \le C .
$$
We  ensure  the existence of functions $P_f \in L^2(J;H^1 (\Omega ) )$,
 $P_f^1 \in L^2 (\Omega \times J;H_{per}^1 (Y ) )$,
 $p \in L^2 (\Omega \times J;H_{per}^1 (Y ) )$,
 $F_1 \in L^2(J;H^1 (\Omega ))$,
 $F^1_1 \in L^2 (\Omega \times J;H_{per}^1 (Y ) )$
 and $C_1 \in L^2 (\Omega \times J;H_{per}^1 (Y ) )$
 such that, up to subsequences not relabeled for convenience, we have
 the following convergence.
 \begin{gather*}
 \Pi^\epsilon p_f^\epsilon
\rightharpoonup P_f  \text{ weakly  in }  L^2(J;H^1(\Omega)), \\
\nabla (\Pi^\epsilon p_f^\epsilon)\
{\stackrel{2}{\rightharpoonup}}   \nabla P_f(x,t) + \nabla_y
P_f^1(x,y,t) ,
\\
 \theta^\epsilon  \stackrel{2}{\rightharpoonup}   p(x,y,t),
\quad \epsilon \nabla \theta^\epsilon \ {\stackrel{2}{\rightharpoonup}}\
\nabla_y p(x,y,t) ,
\\
 \Pi^\epsilon f_1^\epsilon
\rightharpoonup F_1 \quad \text{weakly  in }  L^2(J;H^1(\Omega)),
 \nabla (\Pi^\epsilon f_1^\epsilon ) \
{\stackrel{2}{\rightharpoonup}}\  \nabla F_1(x,t) + \nabla_y
F^1_1(x,y,t) ,
\\
 \xi^\epsilon\
{\stackrel{2}{\rightharpoonup}}\   C_1(x,y,t), \quad \epsilon
\nabla \xi^\epsilon \ {\stackrel{2}{\rightharpoonup}}\   \nabla_y
C_1(x,y,t).
\end{gather*}
We also assert that
$ \Phi^\epsilon \  {\stackrel{2}{\rightharpoonup}}\  \Phi(y)
= \chi_f(y)  \phi_f(y) + \chi_m(y) \phi(y)$,
\[
K^\epsilon \  {\stackrel{2}{\rightharpoonup}}\  K(y)
= \chi_f(y)  k_f(y) + \chi_m(y)  k(y),
\]
and that we can consider that $\Phi^\epsilon$ and $K^\epsilon$ are admissible
test functions for the two-scale convergence.
We give in the following lemma some additional compactness result.


\begin{lemma} \label{lem3}
The sequence  $(\Pi^\epsilon f_1^\epsilon)$ is sequentially compact
 in $L^2(\Omega \times J)$.
\end{lemma}

\begin{proof} We begin by proving that the sequence
$(\Phi^\epsilon \partial_t \xi^\epsilon )$ is uniformly bounded in
$L^2(J;(H^2(\Omega))')$. Let $g \in L^2(J;H^2(\Omega))$. Equations
\eqref{2.15} and \eqref{2.17} give
\begin{align*}
&\langle \Phi^\epsilon  \partial_t\xi^\epsilon , g
\rangle_{L^2(J;(H^2(\Omega))') \times L^2(J;H^2(\Omega)) }
\\
&=  \int_{\Omega \times J} (\frac{K^\epsilon}{\mu(\xi^\epsilon)} \nabla \theta^\epsilon \cdot \nabla \xi^\epsilon) \,  g
-\int_{\Omega \times J}\mathcal{D}^\epsilon  \nabla \xi^\epsilon \cdot
\nabla g
+\int_{\Omega \times J}q\, (\hat f_1-  \xi^\epsilon)\,  g .
\end{align*}
 In view of the previous lemmas, we have
\begin{gather*}
\begin{aligned}
&\Bigl\vert \int_{\Omega \times J} (\frac{K^\epsilon}{\mu(\xi^\epsilon)}
  \nabla\theta^\epsilon  \cdot    \nabla \xi^\epsilon)\, g\,  dx\,dt
\Bigr\vert \\
&\le \frac{C}{\mu_-}\Vert  \vert K^\epsilon \nabla\theta^\epsilon\vert^{1/2}
\nabla \xi^\epsilon \Vert_{(L^2(\Omega \times J))^2}
 \Vert  \vert K^\epsilon \nabla\theta^\epsilon \vert^{1/2}
\Vert_{L^4(\Omega \times J)} \Vert  g \Vert_{L^{4}(\Omega \times J)} \\
&\le C \Vert g \Vert_{L^2(J;H^2(\Omega))} ,
\end{aligned}\\
  \Bigl\vert \int_{\Omega \times J}\mathcal{D}^\epsilon \nabla \xi^\epsilon
\cdot   \nabla g \,  dx\,dt \Bigr\vert
\le C \Vert \nabla g \Vert_{(L^2(J;L^4(\Omega)))^2}
\le C \Vert g \Vert_{L^2(J;H^2(\Omega))}  ,
\\
 \Bigl\vert \int_{\Omega \times J}q\,  (\hat f_1-  \xi^\epsilon)\,  g\,
dx\,dt\Bigr\vert
\le C \Vert g \Vert_{L^2(\Omega \times J)}.
\end{gather*}
Then
$ \vert \langle \Phi^\epsilon  \partial_t\xi^\epsilon , g
\rangle \vert \le C \Vert g \Vert_{L^2(J;H^2(\Omega))}$
and  $(\phi_f^\epsilon \partial_t ( \Pi^\epsilon  f_1^\epsilon))$ is
uniformly bounded in $L^2(J;(H^2(\Omega))')$.
A   compactness argument of Aubin's type ensures that
$(\phi_f^\epsilon ( \Pi^\epsilon  f_1^\epsilon))$ is compact in
$L^2(J;(H^1(\Omega))')$.
We thus can pass to the limit in the product
$ \langle \phi_f^\epsilon( \Pi^\epsilon  f_1^\epsilon) ,
 \Pi^\epsilon  f_1^\epsilon  \rangle_{L^2(J;(H^1(\Omega))')
\times L^2(J;H^1(\Omega))}$.
Since $\phi_f^\epsilon(x) \ge \phi_- >0$ almost everywhere in $\Omega$,
it follows that
$(\Pi^\epsilon  f_1^\epsilon)$ is compact in $L^2(\Omega \times J)$.
\end{proof}

Now, we have the first tools to study the behavior  of the microscopic
 system as $\epsilon$ tends to zero.
We begin by the pressure equations.

\subsection{The Pressure Problem}

We want now to pass to the limit in System \eqref{2.16},
\eqref{2.20}--\eqref{2.21}, (\ref{2.22}), (\ref{2.30}),
\eqref{2.32}  recalled below.
\begin{gather*}
 \Phi^\epsilon  \partial_t \theta^\epsilon - \mathop{\rm div}
((1/\mu(\xi^\epsilon)) K^\epsilon  \nabla \theta^\epsilon ) = q
\quad \text{in }  \Omega \times J,
\\
 K^\epsilon   \nabla \theta^\epsilon \cdot \nu =0 \quad  \text{on }
\Gamma \times J, \quad \theta^\epsilon (x,0)= p^o(x)  \quad \text{in }
  \Omega .
\end{gather*}
 As in \cite{A}, we multiply the first equation by a test function in
 the form
$\Psi(x,t)+\epsilon  \Psi_1 (x,x  /   \epsilon,t ) +\psi(x,x / \epsilon,t)$,
with $\Psi \in \mathcal{D} (\Omega \times J  )$,  $\Psi_1 \in \mathcal{D}
(\Omega \times J;C_{per}^{\infty} (Y ) )$  and
$\psi \in \mathcal{D} (\Omega \times   J;C_{per}^{\infty} (Y ) )$ such that
 $\psi (x,y,t)=0$ for $y \in Y_f$.
Integrating over $\Omega \times J$, we obtain
\begin{align*}
&- \int_{\Omega \times J}\Phi^{\epsilon} (x)\, \theta^{\epsilon}
\partial_t \bigl(\Psi  (x,t )
+\epsilon \,  \Psi_1(x,x/\epsilon ,t ) +\psi(x,x/\epsilon ,t)\bigr) \, dx\,dt
\\
&+ \int_{\Omega }\Phi^{\epsilon} (x)\, p^{o}(x) \bigl(\Psi  (x,0 )
 +\epsilon   \Psi_1(x,x/\epsilon ,0)
+\psi(x,x / \epsilon ,0)\bigr) \, dx
\\
& \quad
+\int_{\Omega \times J}\frac{K^{\epsilon}(x)}{\mu(\xi^\epsilon)}  \nabla  \theta^{\epsilon}
\cdot
\bigl(\nabla \Psi +\epsilon   \nabla_x \Psi_1^\epsilon
+\nabla_y \Psi_1^\epsilon + \nabla_x \psi^\epsilon
+ \frac{1} {\epsilon}  \nabla_y \psi^\epsilon \bigr)\, dx\,dt
\\
& \quad
=\int_{\Omega \times J} q
\bigl(\Psi(x,t) +\epsilon  \Psi_1(x,x / \epsilon ,t) + \psi(x, x / \epsilon,t)
\bigr)\, dx\,dt.
\end{align*}
Letting $\epsilon \to 0$, we get
\begin{align*}
&- \int_{\Omega \times J}\int_{Y_f} \phi_f(y)\,    P_f \, \partial_t \Psi(x,t )
-\int_{\Omega  \times J} \int_{Y_m} \phi(y)\,  p  \,  \partial_t(\Psi(x,t) + \psi(x,y,t))
\\
& + \int_{\Omega }\int_{Y} ( \chi_f(y)\phi_f(y)+\chi_m(y)\phi(y))\,   p^o \,  (\Psi(x,0) + \chi_m(y) \psi(x,y,0))
\\
& +\int_{\Omega \times J}\int_{Y_f} \frac{k_f(y)}{\mu(F_1)}  ( \nabla P_f + \nabla_y P_f^1 )\cdot  ( \nabla \Psi + \nabla_y \Psi_1)
\\
& +\lim_{\epsilon \to 0} \int_{\Omega_m^\epsilon \times J}
 \frac{ \epsilon k^\epsilon}{\mu(C_1^\epsilon)}  \nabla p^\epsilon
  \cdot \nabla_y\psi(x,\frac{x}{\epsilon},t) \\
&=\int_{\Omega \times J \times Y}q\,  (\Psi +\chi_m(y)\psi) .
\end{align*}
By density, this equality holds true for any functions
$ (\Psi  ,\Psi_1 ,\psi  ) \in H_o^1 (\Omega \times J ) \times L^2
(\Omega \times J;H_{per}^1 (Y ) ) \times L^2(\Omega \times J;H_{per}^1(Y_m))$.
Another integration by parts shows (taking successively the test functions
$\psi$, $\Psi_1$ and $\Psi$ equal to zero) that it is a variational
formulation of the following two-scale homogenized system in
$\Omega \times J$:
\begin{gather*}
\overline{\phi_f}^{Y_f}    \partial_tP_f -\mathop{\rm
div} \Bigl( \frac{1}{\mu(F_1)} \int_{Y_f} k_f(y)(\nabla P_f
+\nabla_y P_f^1 )dy\Bigr) =q -\int_{Y_m}\phi(y)\, \partial_t p
\,  dy ,
\\
 -\mathop{\rm div}_y \Bigl( \frac{k_f(y)}{\mu(F_1)}
(\nabla P_f +\nabla_y P_f^1) \Bigr)
=0 \quad \text{in }   Y_f  , \\
 \phi(y) \,  \partial_t p + \mathop{\rm div}_y
\underline{\mathcal{V}}_h = q  \quad \text{in }   Y_m , \quad
\text{where }  - \epsilon\, \frac{k^\epsilon}{\mu(C_1^\epsilon)}
\nabla p^\epsilon  \ {\stackrel{2}{\rightharpoonup}}\
\underline{\mathcal{V}}_h
\\
 k_f(y) (\nabla  P_f +\nabla_y P_f^1)\cdot  \nu_y =0 \quad
\text{on }  \Gamma_{fm}  , \\
 k_f(y) (\nabla P_f + \nabla_y P_f^1)
\cdot \nu = 0 \quad \text{on }  \Gamma ,
\\
 P_f(x,0)= p(x,y,0) = p^o(x) \quad  \text{in }  \Omega
\times Y_m.
\end{gather*}
Now we  eliminate the function $P_f^1$ in the former system.
We use  the solution $(v^i)_{1 \le i \le 2}$ of the cell problem (\ref{4.7}) and  the
homogenized permeability tensor  $\overline{K}_f$ defined by (\ref{4.5}).
Through the relation
$ P_f^1(x,y,t)=\sum_{i=1}^2  \partial_{x_i}P_f(x,t)\,  v^i(y)$,
we recover the following homogenized system.
\begin{gather}
 \overline{\phi_f}^{Y_f} \,    \partial_tP_f - \mathop{\rm div}
\Bigl(\frac{\overline{K}_f}{\mu(C_1)}  \nabla P_f \Bigr) = q -
\int_{Y_m}\phi \,  \partial_t p \, dy \quad \text{in}\  \Omega
\times J,
\label{4.13} \\
  \phi \,  \partial_t p + \mathop{\rm
div}_y \underline{\mathcal{V}}_h = q  \quad \text{in} \  \Omega
\times Y_m \times J,
\label{4.14}\\
  \overline{K}_f   \nabla P_f  \cdot \nu
= 0 \ \text{on}\ \partial \Omega \times J, \  P_f(x,0)= p(x,y,0) =
p^o(x) \  \text{in}\  \Omega \times Y_m . \label{4.15}
\end{gather}
We now determine the nonexplicit limit $\underline{\mathcal{V}}_h$
using a dilation operator.


\subsection{Introduction of an Appropriate Dilation Operator}

Due to the nonlinearities and to the strong coupling of the
problem, the two-scale convergence does not provide an explicit
form for the source-like terms which model the influence of the
matrix at the macroscopic level. To overcome this difficulty, we
follow an idea of \cite{Dou}, that can also be compared with the
unfolding method of Cioranescu et al \cite{Cio}. Firstly, we
recall that in this section we assume no direct flow between the
matrix cells. We thus can assume that the component s of
$\Omega_m^\epsilon$ are strictly separated. For each $\epsilon
>0$, we then define a dilation operator $\widetilde{\  \cdot \  }$
mapping measurable functions on $\Omega_m^\epsilon \times J$ to
measurable functions on $\Omega \times Y_m \times J$ by
$$
\widetilde{u}(x,y,t)=u(c^\epsilon(x)+\epsilon\, y,t) \quad
\text{for}\  y \in Y_m,\  (x,t)\in \Omega \times J,
$$
where
$c^\epsilon(x)$ denotes the lattice translation point of the
$\epsilon$-cell domain containing $x$. This dilation allows to
reach directly the scale of the standard cell. Since
$\Omega_m^\epsilon = \Omega(\epsilon k_1) \cap \left( \cup_{k  \in
\mathbb{Z}^2} \epsilon (Y_m +k) \right)$, we recall that $c^\epsilon(x) =
\epsilon  k$ for $x \in \epsilon (Y_m +k)$. Thus the function
$\widetilde{u}$ is constant in $x$ on each block $ \epsilon (Y_m
+k)$, $ k \in \mathbb{Z}^2$, of $\Omega$. We extend this operator from
$Y_m$ to $\cup_{k}(Y_m +k)$ periodically. The dilation has the
following properties ({\it cf} \cite{Dou}). Any  function $u \in
L^2(J;H^1(\Omega_m^\epsilon))$ satisfies
$$
\Vert \widetilde{u} \Vert_{L^2(\Omega_m^\epsilon  \times J \times
Y_m)} = \Vert u \Vert_{L^2(\Omega  \times J)}, \quad \nabla_y
\widetilde{u} =\epsilon\,  \widetilde{\nabla_x u} \  \text{a.e.\
in}\  \Omega  \times J \times Y_m .
$$
Moreover, if $(v,w) \in (L^2(0,T;H^1(\Omega_m^\epsilon)))^2$, then
\begin{gather*}
  \left(  \widetilde{v},\widetilde{w} \right)_{L^2(\Omega \times J \times Y_m )}
= \left( v,w \right)_{L^2(\Omega_m^\epsilon \times J)} ,
\\
\left(  \widetilde{v},w \right)_{L^2(\Omega \times J \times Y )}
= \left( v,\widetilde{w} \right)_{L^2(\Omega \times J \times Y )},
\\
\left\Vert \nabla_y \widetilde{v} \right\Vert_{(L^2(\Omega \times
J  \times Y_m ))^2} =\epsilon \bigl\Vert \widetilde{ \nabla_x v }
\bigr\Vert_{(L^2(\Omega_m^\epsilon \times J ))^2}  .
\end{gather*}
The following result makes the link between two-scale convergence and
weak convergence of dilated sequences.
We refer to  \cite{BLM} for its detailed proof.

\begin{proposition} \label{prop2}
If $(v^\epsilon)$ is a bounded sequence of $L^2(\Omega_m^\epsilon
\times J)$ such that $\widetilde{v^\epsilon}$ converges weakly to
$v^o$ in  $L^2(\Omega \times J ; L_{per}^2( Y_m) )$ and
$\chi_m^\epsilon   v^\epsilon$ two-scale converges to $v$, then we
have $v^o=v$ a.e.  in  $\Omega  \times J \times Y_m$.
\end{proposition}

We now determine the limit behavior of the dilated solutions $\widetilde C_1^\epsilon $
and $\widetilde p^\epsilon $. The outline of this study is the following. We find the
equations satisfied by the dilated  solutions in $\Omega \times
Y_m \times J$. Some estimates lead to the convergence of $(\widetilde p^\epsilon ,\widetilde C_1^\epsilon )$
to $(p,C_1)$    as $\epsilon \to 0$. On the other hand, for each
fixed $k \in \mathbb{Z}^2$, we note that the restriction in
$\epsilon(Y_m+k)$ of $(\widetilde p^\epsilon ,\widetilde C_1^\epsilon )$ is independent of $x$. Thus we
define corresponding functions  $(\widetilde p^\epsilon_k (y,t),\widetilde C^\epsilon_{1k}(y,t))$. We get
enough compactness results to find the equations satisfied by
their limit $(p_k,C_{1k})$ for each $k \in \mathbb{Z}^2$. Then, by a
density argument, the limit equations for $(p,C_1)$ are deduced.


We begin with finding the equations satisfied by the dilated solutions
$\widetilde C_1^\epsilon $ and $\widetilde p^\epsilon $.
We define a test function $\hat{\psi}$ by
\begin{gather*}
\hat{\psi}(x,z,t) = \begin{cases}
 \psi \bigl( (z - c^\epsilon(x))/\epsilon , t \bigr)
 &\text{for }   z \in \epsilon Y_m + c^\epsilon(x) \\
 0 & \text{for }   z \not \in \epsilon  Y_m + c^\epsilon(x)  ,
\end{cases}
\end{gather*}
for any function $\psi \in L^2(J;H_o^1(Y_m))$.
We multiply   \eqref{2.20} by $\hat{\psi}$ and integrate over
$\Omega_m^\epsilon$.
Since $ \Omega_m^\epsilon = \cup_{x \in \Omega} \left(\epsilon
 Y_m + c^\epsilon(x) \right)$
and $ (\epsilon   Y_m + c^\epsilon(x_1) ) \cap (\epsilon
 Y_m +c^\epsilon(x_2) ) = \emptyset$ for $x_1 \ne x_2$,
we get for almost every $x \in \Omega_m^\epsilon$
\begin{align*}
&\int_{ \epsilon  Y_m + c^\epsilon(x)  } \!\! \bigl(
\phi^\epsilon(z)\, \partial_t   p^\epsilon)(z,t)\,
\hat{\psi}(x,z,t) + \epsilon^2
\frac{k^\epsilon(z)}{\mu(C_1^\epsilon)}  \nabla p^\epsilon(z,t)
\cdot \nabla_z \hat{\psi}(x,z,t) \bigr) \, dz
\\
& = \int_{ \epsilon  Y_m + c^\epsilon(x)  } q\,
\hat{\psi}(x,z,t)\,  dz .
\end{align*}
With the change of variable $z \mapsto \epsilon  y + c^\epsilon(x)
= \epsilon   ( y + k  )$, $k \in \mathbb{Z}^2$, we recover the variational
formulation of the  equation
\begin{equation}
\phi(y)\,   \partial_t \widetilde p^\epsilon   +\mathop{\rm div}_y  \widetilde q^\epsilon  =\widetilde{q},
\quad \widetilde q^\epsilon =-(1/\mu(\widetilde C_1^\epsilon ))k(y)  \nabla_y \widetilde p^\epsilon  . \label{4.16}
\end{equation}
We proceed in the same way for the concentration   \eqref{2.17}.
We get
\begin{equation}
\phi(y)\,  \partial_t \widetilde C_1^\epsilon  + \widetilde q^\epsilon  \cdot \nabla_y \widetilde C_1^\epsilon  -\mathop{\rm
div}_y  (\mathcal{D}(y,\widetilde q^\epsilon )\nabla_y \widetilde C_1^\epsilon  ) + \widetilde{q}\,  \widetilde C_1^\epsilon  =
\widetilde{q} \hat f_1 . \label{4.17}
\end{equation}
This system is provided with the following boundary and initial conditions.
\begin{gather}
 \widetilde p^\epsilon  = \widetilde{P_f^\epsilon}\quad \text{and}\quad
  \widetilde C_1^\epsilon  = \widetilde{f_1^\epsilon} \quad \text{in }  H^{1/2}(\partial Y_m)
\text{ for }  (x,t) \in \Omega \times J,
\label{4.18}\\
 \widetilde p^\epsilon (x,y,0) =\widetilde{p^o}(x,y), \quad \widetilde C_1^\epsilon (x,y,0)=\widetilde{C_1^o}(x,y) \quad
\text{in } \Omega \times Y_m .
\label{4.19}
\end{gather}
Using the basic properties of the dilation operator and the estimates
of Section 3, we claim  that, for subsequences not relabeled for convenience,
the following convergence take place.
\begin{gather*}
  \widetilde p^\epsilon   \rightharpoonup  p , \  \widetilde C_1^\epsilon  \rightharpoonup   C_1 \quad
\text{weakly  in }  L^2(\Omega \times J \times Y_m) ,
\\
 \nabla_y \widetilde p^\epsilon  \rightharpoonup \nabla_y p, \  \nabla_y \widetilde C_1^\epsilon  \rightharpoonup \nabla_y C_1
 \quad \text{weakly  in } (L^2(\Omega \times J \times Y_m ))^2.
\end{gather*}
We then choose a fixed $k \in \mathbb{Z}^2$.
For each fixed $\epsilon >0$, we  define the functions $\widetilde p^\epsilon_k $ and
$\widetilde C^\epsilon_{1k}$ in $Y_m \times J$ by
\begin{gather*}
\widetilde p^\epsilon_k (y,t) = \begin{cases}
 \widetilde p^\epsilon (x,y,t)_{/x \in \epsilon (Y_m +k)}
& \text{if $k$ is   such  that } \epsilon (Y_m +k) \cap \Omega \ne \emptyset , \\
 0 &\text{otherwise,}
\end{cases}
 \\
\widetilde C^\epsilon_{1k}(y,t) =\begin{cases}
 \widetilde C_1^\epsilon (x,y,t)_{/x \in \epsilon (Y_m +k)}
&\text{if  k  is  such  that }  \epsilon (Y_m +k) \cap \Omega \ne \emptyset , \\
 0 &\text{otherwise.}
\end{cases}
\end{gather*}
Roughly speaking, $k=(k_1,k_2) \in \mathbb{Z}^2$ is such that
$\epsilon(Y_m +k) \cap \Omega \ne \emptyset$ if $k_i < |\Omega|_i / \epsilon$,
 $i=1,2$ (where $|\Omega|_i$ denotes the value of the measure of $\Omega$
in the $i$-th direction).
For each $\epsilon >0$ such that $\epsilon (Y_m +k) \cap \Omega \ne \emptyset$,
 $(\widetilde p^\epsilon_k ,\widetilde C^\epsilon_{1k})$ is a solution of   Pb.  (\ref{4.16})-(\ref{4.19}) in
$Y_m \times J$.
Furthermore, since any $f \in L^2(\Omega \times J)$ satisfies
$$
\Vert \widetilde{f_k} \Vert_{L^2(Y_m \times J)}
= \frac{1}{\epsilon |Y_m|} \Vert \widetilde{f} \Vert_{L^2( \epsilon(Y_m+k)
\times Y_m \times J)}
\le  \frac{1}{\epsilon |Y_m|} \Vert \widetilde{f} \Vert_{L^2( \Omega
\times J  \times Y_m )} ,
$$
we have enough regularity properties to get with  (\ref{4.16})-(\ref{4.19})
some  estimates for $\widetilde p^\epsilon_k $ and $\widetilde C^\epsilon_{1k}$.
They are  analogous to those obtained for $p_f^\epsilon$ and $f_1^\epsilon$
in the fracture.
We claim that
\begin{gather*}
\Vert \widetilde p^\epsilon_k  \Vert_{L^\infty( J;L_{per}^2(Y_m)) \cap L^2(  J;H_{per}^1(Y_m))}
\le C ,
\quad  0 \le \widetilde C^\epsilon_{1k}(y,t) \le 1 \quad  \text{a.e.  in }  Y_m \times J ,
\\
\bigl\Vert \bigl( D_m + \alpha_t  |\nabla_y \widetilde p^\epsilon_k  | \bigr)^{1/2}
\nabla_y \widetilde C^\epsilon_{1k} \bigr\Vert_{(L^2( J \times Y_m ))^2}
+  \Vert \partial_t \widetilde C^\epsilon_{1k} \Vert_{L^2( J;(H_{per}^2(Y_m))')} \le C ,
\end{gather*}
where $C$ is a generic constant, independent of $\epsilon$ and $k$.
Using in particular compactness arguments of Aubin's type for the
concentration $\widetilde C^\epsilon_{1k}$, we deduce of these latter estimates that,
up to subsequences not relabeled for convenience, the following convergence
 take place.
\begin{gather*}
\widetilde p^\epsilon_k  \rightharpoonup p_k \quad \text{weakly  in }  L^2(J;H_{per}^1(Y_m)) , \\
\widetilde C^\epsilon_{1k} \rightharpoonup C_{1k} \quad \text{weakly in } L^2( J;H_{per}^1(Y_m))
\text{ and   a.e.   in }  Y_m \times J ,
\end{gather*}
for some limit functions $p_k \in L^2( J;H_{per}^1(Y_m))$ and
$C_k \in L^2( J;H_{per}^1(Y_m))$.
We can then easily pass to the limit  $\epsilon \to 0$ in the pressure
equation (\ref{4.16}) satisfied by $\widetilde p^\epsilon_k $.
We get
\begin{gather}
 \phi(y)\,  \partial_t p_k - \mathop{\rm div}_y  (
(1/\mu(C_{1k})) k(y)  \nabla_y p_k ) =q \quad \text{in } Y_m
\times J . \label{4.20}
\end{gather}
This equation is satisfied for each $k \in \mathbb{Z}^2$. It reminds to
make the link with the limit $(p,C_1)$ of $(\widetilde p^\epsilon ,\widetilde C_1^\epsilon )$ by density
arguments. On the first hand, for each $k \in \mathbb{Z}^2$, there exists
$\epsilon(k)>0$ such that   $\epsilon(Y_m +k) \cap \Omega \ne
\emptyset$ and then $\widetilde p^\epsilon_k (y,t) = \widetilde p^\epsilon (x,y,t)_{/  x \in \epsilon
(Y_m+k)}$ for any $ \epsilon < \epsilon(k)$. Since $\lim_{\epsilon
\to 0} \vert \epsilon (Y_m +k) \vert = 0$, the subset $ \epsilon
(Y_m +k) \cap \Omega $ tends when $\epsilon \to 0$ to a point
$\{x_k\} \subset \Omega$. Then $p_k(y,t)=p(x_k,y,t)$. In the same
way, we get $C_{1k}(y,t)=C_1(x_k,y,t)$. On the other hand,
$\lim_{\epsilon \to 0} \vert \Omega  \setminus ( ( \cup_{k \in
\mathbb{Z}^2} \epsilon (Y_m +k) ) \cap    \Omega ) \vert = 0$. Thus the
set $ \left( \cup_{k \in \mathbb{Z}^2} \{x_k \} \right)$ is dense in
$\Omega$. Since the functions $(p(x_k,y,t), C_1(x_k,y,t))$ satisfy
(\ref{4.20}) for any $k \in \mathbb{Z}^2$, we conclude by density that $p$
is a solution of (\ref{4.9}).

To end this subsection, we add the following result of strong
convergence.

\begin{lemma} \label{lem4}
The sequence $\mu(C_1^\epsilon)^{-1/2} \epsilon \nabla p^\epsilon$ satisfies
$$
 \lim_{\epsilon \to 0} \Vert \mu(C_1^\epsilon)^{-1/2} \epsilon
\nabla p^\epsilon \Vert_{(L^2(\Omega_m^\epsilon \times J ))^2} =
\Vert \mu(C_1)^{-1/2}  \nabla_y p \Vert_{(L^2(\Omega  \times J
\times Y_m ))^2}.
$$
\end{lemma}

 Using the terminology of \cite{A},
$\mu(C_1^\epsilon)^{-1/2} \epsilon \nabla p^\epsilon$ is said strongly
two-scale converging to  $\mu(C_1)^{-1/2}  \nabla_y p$.
 For any bounded sequence $v^\epsilon$ of $(L^2(\Omega \times J))^2$
such that $v^\epsilon \  {\stackrel{2}{\rightharpoonup}}\ v$ and any
test admissible function $\psi$, we can assert that
\begin{align*}
& \lim_{\epsilon \to 0} \int_{\Omega \times J} ( \mu(C_1^\epsilon)^{-1/2}
\epsilon \nabla p^\epsilon \cdot v^\epsilon ) \, \psi(x,x / \epsilon,t)
 \, dx\,dt
\\
& =\int_{\Omega \times J} \int_{Y_m} (\mu(C_1)^{-1/2}  \nabla_y p \cdot v)
\, \psi(x,y,t) \, dx\,dy\,dt .
\end{align*}


\begin{proof}
We multiply   (\ref{4.16}) by $\widetilde p^\epsilon$ and we integrate over $\Omega \times (0,t) \times Y_m=\Omega_t \times Y_m$, for $t \in J$.
Letting $\epsilon \to 0$ and using   (\ref{4.9}), we  write
\begin{align*}
& \lim_{\epsilon \to 0} \Bigl(
\frac{d}{2dt} \int_{\Omega \times Y_m} \phi(y) \, \vert \widetilde p^\epsilon  \vert^2
- \frac{1}{2} \int_{\Omega \times Y_m} \phi(y) \, \vert {\widetilde p^o}
\vert^2
+ \int_{\Omega_t \times Y_m} \frac{k(y)}{\mu(\widetilde C_1^\epsilon )} \nabla_y \widetilde p^\epsilon  \cdot \nabla_y \widetilde p^\epsilon 
\Bigr)
 \\
&= \lim_{\epsilon \to 0} \int_{\Omega_t \times Y_m} \widetilde{q}\, \widetilde p^\epsilon 
= \int_{\Omega_t \times Y_m} q\, p
\\
&= \frac{d}{2dt} \int_{\Omega \times Y_m} \phi(y) \, \vert p \vert^2
- \frac{1}{2} \int_{\Omega \times Y_m} \phi(y) \, \vert {\widetilde p^o}
\vert^2
+ \int_{\Omega_t \times Y_m} \frac{k(y)}{\mu(C_1)} \nabla_y p \cdot
 \nabla_y p .
\end{align*}
Bearing in mind that $k$ is a symmetric  definite positive tensor, we
conclude in particular from the letter relation that
\begin{align*}
&\lim_{\epsilon \to 0} \Bigl\Vert \frac{1}{\mu(\widetilde C_1^\epsilon )^{1/2}}  \nabla_y \widetilde p^\epsilon \Bigr\Vert_{(L^2(\Omega  \times J \times Y_m))^2}
=
 \Bigl\Vert \frac{1}{\mu(C_1)^{1/2}}  \nabla_y p \Bigr\Vert_{(L^2(\Omega  \times J \times Y_m ))^2}
 \\
&= \lim_{\epsilon \to 0} \Bigl\Vert \frac{1}{\mu(\widetilde C_1^\epsilon )^{1/2}} \widetilde{\epsilon \nabla p^\epsilon }\Bigr\Vert_{(L^2(\Omega  \times J \times Y_m))^2}
= \lim_{\epsilon \to 0} \Bigl\Vert \ensuremath{ \widetilde{ \frac{1}{\mu(C_1^\epsilon)^{1/2}}\epsilon \nabla p^\epsilon }}\Bigr\Vert_{(L^2(\Omega  \times J \times Y_m))^2}
\\
&= \lim_{\epsilon \to 0} \Bigl\Vert  \frac{1}{\mu(C_1^\epsilon)^{1/2}}\epsilon \nabla p^\epsilon \Bigr\Vert_{(L^2(\Omega_m^\epsilon  \times J ))^2}.
\end{align*}
Lemma \ref{lem4} is proved.
\end{proof}

\subsection{Corrector for the flux function and concentration problem}

Now we study  the behavior  of the concentration problem
\eqref{2.15}-\eqref{2.17} as $\epsilon$ tends to zero. We recall
it is
\begin{gather*}
 \Phi^\epsilon\,  \partial_t  \xi^\epsilon -
(1/\mu(\xi^\epsilon)) K^\epsilon  \nabla \theta^\epsilon \cdot
\nabla \xi^\epsilon - \mathop{\rm div}  ( \mathcal{D}^\epsilon
\nabla \xi^\epsilon) +q\,  \xi^\epsilon
=q \hat{f_1} , \\
 \mathcal{D}^\epsilon   \nabla \xi^\epsilon \cdot \nu = 0   ,
\quad \xi^\epsilon (x,0)=f_1^o(x)   .
\end{gather*}
The starting point of the asymptotic study is similar to the one used
for the pressure problem in Subsection 4.1.
We multiply the latter equation by a test function
$\Psi (x,t )+\epsilon   \Psi_1 (x, x/\epsilon,t ) + \psi
 (x,x/ \epsilon ,t  )$ and we integrate over $\Omega \times J$.
 Letting $\epsilon \to 0$, we obtain the limit variational formulation
corresponding to \eqref{2.15}-\eqref{2.17}.
The difficulty is due to the nonlinearities involving the Darcy velocity
in the convective terms and in the dispersive ones.

We begin by deriving the matrix equation (\ref{4.10}).
Indeed we have obtained sufficient compactness results in Subsection 4.2
to pass to the limit via the dilated equation (\ref{4.17}).
On the one hand, using the a.e. convergence of $\widetilde C^\epsilon_{1k}$ to $C_{1k}$,
we get by density
\begin{align*}
& \phi(y) \, \partial_t C_1
- \mu(C_1)^{-1/2} k(y) \bigl( \lim_{\epsilon \to 0} \mu(\widetilde C_1^\epsilon )^{-1/2} \nabla_y \widetilde p^\epsilon  \bigr) \cdot \nabla_y C_1
\\
& -\mathop{\rm div}_y \bigl( \mu(C_1)^{-1/2} \mathcal{D}
\bigl( \lim_{\epsilon \to 0} \mu(\widetilde C_1^\epsilon )^{-1/2} \nabla_y \widetilde p^\epsilon  \bigr)
\nabla_y C_1 \bigr) \\
&= q\, (\hat{f_1} - C_1).
\end{align*}
On the other hand, we note that Lemma \ref{lem4} gives the strong
convergence of $\widetilde{ \mu(C_1^\epsilon)^{-1/2} \epsilon
\nabla p^\epsilon } = \mu(\widetilde C_1^\epsilon )^{-1/2} \nabla_y \widetilde p^\epsilon $ to
$\mu(C_1)^{-1/2} \nabla_y p$ in $(L^2(\Omega \times J \times
Y_m))^2$. It is sufficient to explicit all the limits in the
former equation. We obtain the limit equation (\ref{4.10}).

Passing to the limit in the fractured part is less obvious.
We are going to derive a corrector for the Darcy velocity in the fracture.
Our aim is to find a function
$\underline{\mathcal{V}}_o(x,y,t) \in (L^2(\Omega \times J \times Y))^2$
such that $\underline{\mathcal{V}}_o(x,x/\epsilon,t)$ is an admissible
test function for the two-scale convergence and such that
\[
\lim_{\epsilon \to 0} \Vert \chi_f^\epsilon \bigl( -(K^\epsilon/
\mu(\xi^\epsilon)) \nabla \theta^\epsilon
- \underline{\mathcal{V}}_o(x,x/\epsilon,t) \bigr) \Vert_{(L^2(\Omega
 \times J))^2}=0.
\]
Then we will pass to the limit in
$\chi_f^\epsilon \underline{\mathcal{V}}_o(x,x/\epsilon,t) \cdot
\nabla \xi^\epsilon$ instead of $\chi_f^\epsilon (K^\epsilon/
\mu(\xi^\epsilon))\nabla \theta^\epsilon \cdot \nabla \xi^\epsilon$;
respectively we will pass to the limit in $\chi_f^\epsilon
\mathcal{D}(\underline{\mathcal{V}}_o(x,x/\epsilon,t))
\nabla \xi^\epsilon$ instead of
$\chi_f^\epsilon \mathcal{D}(K^\epsilon/\mu(\xi^\epsilon))\nabla
\theta^\epsilon) \nabla \xi^\epsilon$.
We begin with the following result.

\begin{lemma} \label{lem5}
The following convergence hold.
\begin{align*}
&\lim_{\epsilon \to 0}\int_{\Omega \times J} \frac{K^{\epsilon}}
{\mu(\xi^\epsilon)} \nabla  \theta^\epsilon \cdot \nabla \theta^\epsilon \,
  dx\,dt \\
&=\int_{\Omega \times J}\int_{Y_f} \frac{k_f(y)}{\mu(F_1)}
(\nabla P_f +\nabla_y P_f^1) \cdot (\nabla P_f\\
& \quad +\nabla_y P_f^1) \,  dx\,dy\,dt + \int_{\Omega \times
J}\int_{Y_m} \frac{k(y)}{\mu(C_1)} \nabla_y p \cdot \nabla_y p\,
dx\,dy\,dt.
\end{align*}
\end{lemma}

\begin{proof}
Let us consider  the following energy equation for $t \in J$.
We denote $\Omega_t$ the set $\Omega \times (0,t)$.
\[
\frac{1}{2}\int_\Omega \Phi^\epsilon ( \theta^\epsilon(x,t)^2
-   p^o(x)^2 ) dx
+ \int_{\Omega_t} \frac{K^\epsilon}{\mu(\xi^\epsilon)}
 \nabla \theta^\epsilon \cdot \nabla \theta^\epsilon   \,dx\,ds
=\int_{\Omega_t} q\,  \theta^\epsilon   \,dx\,ds .
\]
In view of the two-scale convergence of $\theta^\epsilon$ we have
$$
\lim_{\epsilon \to 0} \int_{\Omega_t} q\,  \theta^\epsilon \,  \,dx\,ds
= \int_{\Omega_t}\int_Y q \left( P_f(x,t) +\chi_m(y)  p(x,y,t)
\right)dx\,dy\ ds,
$$
and so, in view of the variational formulation obtained in the former
subsection,
\begin{align*}
&\lim_{\epsilon \to 0} \Bigl(
\frac{1}{2}\int_\Omega \Phi^\epsilon   {\theta^\epsilon(x,t)}^2\,  dx
 + \int_{\Omega_t} \frac{K^\epsilon}{\mu(\xi^\epsilon)}
\nabla \theta^\epsilon \cdot \nabla \theta^\epsilon  \,  \,dx\,ds
 \Bigr) \\
&=\frac{1}{2}\int_\Omega \int_Y \Phi ( P_f + \chi_m p )^2(x,y,t)\, dx\,dy
+\int_{\Omega_t}\int_{Y_f} \frac{k_f(y)}{\mu(F_1)}  (\nabla P_f +\nabla_y P_f^1)  \cdot  (\nabla P_f
\\
&\quad +\nabla_y P_f^1) \,  dx\,dy  ds
+ \int_{\Omega_t}\int_{Y_m} \frac{k(y)}{\mu(C_1)}\  \nabla_y p
\cdot \nabla_y p \,  dx\,dy ds,
\end{align*}
where we recall that
$\Phi(y)=\chi_f(y)\phi_f(y)+\chi_m(y)\phi(y)$. The limit of each
term in the left-hand side of the last relation is larger than the
corresponding two-scale limit in the right-hand side. Thus
equality holds for each contribution and Lemma \ref{lem5} is
proved.
\end{proof}

Now, comparing the results of Lemmas \ref{lem4} and \ref{lem5}, we
conclude that
\begin{align*}
&\lim_{\epsilon \to 0} \int_{\Omega \times J} \chi_f^\epsilon
\frac{K^\epsilon}{\mu(\xi^\epsilon)} \nabla \theta^\epsilon \cdot
\nabla \theta^\epsilon \, dx\,dt
\\
&=  \int_{\Omega \times J} \int_{Y_f} \frac{k_f(y)}{\mu(F_1)}
(\nabla P_f +\nabla_y P_f^1)  \cdot  (\nabla P_f
 +\nabla_y P_f^1) \,  dx\,dy\,dt .
 \end{align*}
We recall that $K^\epsilon$ is a symmetric definite positive
tensor and that it is considered as an admissible test function
for the two-scale convergence. Furthermore we have showed in Lemma
\ref{lem3} that $\chi_f^\epsilon( \Pi^\epsilon f_1^\epsilon)$
converges almost everywhere in $\Omega \times J$ to $F_1$. The
latter relation then leads to
$$
\lim_{\epsilon \to 0} \Bigl\Vert \chi_f^\epsilon \frac{K^\epsilon}
{\mu(\xi^\epsilon)} \nabla \theta^\epsilon
\Bigr\Vert_{(L^2(\Omega \times J))^2}
=
\Bigl\Vert \chi_f(y) \frac{k_f(y)}{\mu(F_1)} (\nabla P_f
+ \nabla_y P_f^1) \Bigr\Vert_{(L^2(\Omega \times J \times Y))^2}.
$$
This relation is sufficient to assert the following result.


\begin{proposition} \label{prop3}
Let us define
$$ \underline{\mathcal{V}}_o(x,y,t)=
 -\chi_f(y) \frac{k_f(y)}{\mu(F_1(x,t))}  (\nabla P_f(x,t)
 +\nabla_y P_f^1(x,y,t) ) .
$$
 Assume that $ \underline{\mathcal{V}}_o$ is an admissible test function
for the two-scale convergence, that is
$ \lim_{\epsilon \to 0} \int_{\Omega \times J} \vert
\underline{\mathcal{V}}_0(x,x/\epsilon,t) \vert^2  dx\,dt
\le   \int_{\Omega \times J}\int_Y \vert \underline{\mathcal{V}}_0(x,y,t)
\vert^2  dxdy\,dt $.
The function $ \underline{\mathcal{V}}_0(x,x/\epsilon,t)$ is a corrector
for the flux function in the following sense.
$$
\lim_{\epsilon \to 0} \Bigl\Vert - \chi_f^\epsilon
\frac{K^\epsilon}{\mu(\xi^\epsilon)} \nabla \theta^\epsilon
-  \underline{\mathcal{V}}_0(x,\frac{x}{\epsilon},t)\Bigr\Vert_{(L^2
(\Omega \times J))^2} = 0.
$$
\end{proposition}

We then can substitute $\underline{\mathcal{V}}_0(x,x/\epsilon,t)$
to $\chi_f^\epsilon({K^\epsilon}/\mu(\xi^\epsilon)) \nabla
\theta^\epsilon $ to pass to the limit in the different equations.
We get  (\ref{4.3}). This completes the proof of Theorem
\ref{thm2}.


\section{The macroscopic model: the case of a partially fractured
medium ($\alpha >0$)}

In the present section, we derive rigorously the homogenized model
corresponding to the partially fractured model
\eqref{2.15}-\eqref{2.21}, \eqref{2.23}-\eqref{2.32} when $\alpha>
0$. The limit model is described in Theorem \ref{thm3} at the end
of the section.

We begin by recalling the estimates derived in Section 3.
\begin{gather*}
 \Vert  \chi_f^\epsilon p_f^\epsilon \Vert_{L^\infty(J;L^2(\Omega))} + \Vert \chi_f^\epsilon  \nabla p_f^\epsilon \Vert_{(L^2(\Omega \times J))^2}
\le C ,
\\
 \Vert  \chi_m^\epsilon p^\epsilon \Vert_{L^\infty(J;L^2(\Omega))} + \Vert \chi_m^\epsilon \alpha(c_1^\epsilon+c_2^\epsilon) \nabla p^\epsilon \Vert_{(L^2(\Omega \times J))^2}
\le C ,
\\
 \Vert \chi_f^\epsilon f_1^\epsilon \Vert_{L^\infty(\Omega \times J)}
+ \Vert \chi_f^\epsilon (1+ \vert \nabla p_f^\epsilon \vert) \nabla f_1^\epsilon \Vert_{(L^2(\Omega \times J))^2}
\le C,
 \\
 \Vert \chi_m^\epsilon C_1^\epsilon \Vert_{L^\infty(J;L^2(\Omega))}
+ \Vert \chi_m^\epsilon (1+ \vert \underline{\mathcal{V}}_s^\epsilon \vert) \nabla C_1^\epsilon \Vert_{(L^2(\Omega \times J))^2}
\le C,
\\
 \Vert \chi_m^\epsilon c_i^\epsilon \Vert_{L^\infty(J;L^2(\Omega))}
+ \Vert \chi_m^\epsilon (1+ \vert \underline{\mathcal{V}}_s^\epsilon \vert) \nabla c_i^\epsilon \Vert_{(L^2(\Omega \times J))^2}
\le C, \  i=1,2.
\end{gather*}
 For the rest of this article, we assume that the pressure equation
\eqref{2.20} is not degenerate, that is the constant $C_\alpha$
defined in (\ref{2.34}) satisfies
$$
C_\alpha >0.
$$
The more technical degenerate case when $C_\alpha =0$ is postponed
to a forthcoming  paper. In view of (\ref{2.39}) in Theorem
\ref{thm1}, we now have $ 0 < C_\alpha \le
c_1^\epsilon(x,t)+c_2^\epsilon(x,t) \le 1$ a.e.  in $\Omega \times
J$ and $\Vert \chi_m^\epsilon \nabla p^\epsilon \Vert_{(L^2(\Omega
\times J))^2}
 \le C $.
Thus there exist functions $p \in L^\infty(J;L^2(\Omega))
\cap L^2(J; H^1(\Omega))$,   $p^1 \in L^2(\Omega \times J;H^1_{per}(Y))$,
 $(f_1,C_1,c_1,c_2) \in (L^\infty(J;L^2(\Omega)) \cap L^2(J; H^1(\Omega)))^4$
and $(f_1^1,C_1^1,c_1^1,c_2^1) \in (L^2(\Omega \times J;H^1_{per}(Y)))^4$
such that, up to extracted subsequences, as $\varepsilon \to 0$,
\begin{gather*}
 \theta^\epsilon=\chi_f^\epsilon p_f^\epsilon + \chi_m^\epsilon p^\epsilon \  {\stackrel{2}{\rightharpoonup}}\   p,\quad
 \nabla \theta^\epsilon\  {\stackrel{2}{\rightharpoonup}}\  \nabla p(x,t) + \nabla_y p^1(x,y,t)),
\\
 \chi_f^\epsilon f_1^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_f(y) f_1,\quad
\chi_f^\epsilon \nabla f_1^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_f(y) (\nabla f_1(x,t) + \nabla_y f_1^1(x,y,t)),
\\
 \chi_m^\epsilon C_1^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_m(y) C_1,\quad
\chi_m^\epsilon \nabla C_1^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_m(y) (\nabla C_1(x,t) + \nabla_y C_1^1(x,y,t)),
\\
 \chi_m^\epsilon c_i^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_m(y) c_i,\quad
\chi_m^\epsilon \nabla c_i^\epsilon\  {\stackrel{2}{\rightharpoonup}}\   \chi_m(y) (\nabla c_i(x,t) + \nabla_y c_i^1(x,y,t)),\  i=1,2.
 \end{gather*}
We begin with the following result linking the limit concentrations
$f_1$ and $m_1=\alpha c_1+\beta C_1$.


 \begin{lemma} \label{lem6}
 The concentrations $f_1(x,t)$ and $m_1(x,t)=\alpha c_1(x,t)+\beta C_1(x,t)$
are equal almost everywhere in $\Omega \times J$.
 \end{lemma}

\begin{proof}
 Let $c^\epsilon = \chi_f^\epsilon f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon
\in L^2(J;H^1(\Omega))$.
 It satisfies $\gamma_f^\epsilon c^\epsilon = \gamma_f^\epsilon f_1^\epsilon
= \gamma_m^\epsilon m_1^\epsilon = \gamma_m^\epsilon c^\epsilon$ and
$\epsilon \nabla c^\epsilon = \epsilon \chi_f^\epsilon \gamma f_1^\epsilon
+ \epsilon \chi_m^\epsilon m_1^\epsilon \in (L^2(\Omega \times J))^2$.
 We know that $c^\epsilon\  {\stackrel{2}{\rightharpoonup}}\  \chi_f(y) f_1
+ \chi_m(y) m_1$ and $\epsilon \nabla c_\epsilon
 {\stackrel{2}{\rightharpoonup}}  0$.
 For any $\underline{\Psi} \in (\mathcal{C}_o^\infty(\Omega;\mathcal{C}_{per}^\infty(Y)))^2$ we write
 $$
\int_\Omega  \epsilon \nabla c^\epsilon \cdot \underline{\Psi}
(x,\frac{x}{\epsilon}) \, dx
 = - \int_\Omega c^\epsilon \Bigl( \epsilon \mathop{\rm div}_x
\underline{\Psi}(x,\frac{x}{\epsilon}) + \mathop{\rm div}_y
\underline{\Psi}(x,\frac{x}{\epsilon}) \Bigr) \, dx.
 $$
 We take the two-scale limits on both sides.
 We get
 \begin{align*}
  0& =- \int_\Omega \int_Y (\chi_f(y) f_1(x,t)
   + \chi_m(y) m_1(x,t) ) \, \mathop{\rm div}_y \underline{\Psi}(x,y) \, dx\,dy
 \\
 & =- \int_\Omega \int_{\partial Y_f} f_1(x,t) \underline{\Psi}(x,s)
  \cdot \nu_f \, dx\,ds
 - \int_\Omega \int_{\partial Y_m} m_1(x,t) \underline{\Psi}(x,s)
 \cdot \nu_m \, dx\,ds .
 \end{align*}
 This proves that $f_1(x,t)=m_1(x,t)$ for $s \in \partial Y_f \cap \partial
Y_m = \Gamma_{fm}$ and thus
 $f_1(x,t)=m_1(x,t)$ a.e. in $\Omega \times J$.
 \end{proof}

We then claim and prove the following compactness result.

\begin{lemma} \label{lem7}
The sequences $(\chi_f^\epsilon f_1^\epsilon)$,
$(\chi_m^\epsilon m_1^\epsilon)$ and
$(\chi_m^\epsilon (c_1^\epsilon+c_2^\epsilon))$ are sequentially compact
in $L^2(\Omega \times J)$.
\end{lemma}

\begin{proof} On the one hand, let $\psi \in
L^2(J;H^2(\Omega))$. We multiply \eqref{2.15} by
$(\alpha^2+\beta^2)\chi_f^\epsilon \psi$ and (\ref{2.40}) by
$\chi_m^\epsilon \psi$. We integrate over $\Omega \times J$ and
sum up the results. Using the same type of arguments than in the
proof of Lemma \ref{lem3}, we conclude that the sequence
$\partial_t (\phi_f^\epsilon(\alpha^2+\beta^2)\chi_f^\epsilon
f_1^\epsilon+\phi^\epsilon\chi_m^\epsilon m_1^\epsilon)$ is
uniformly bounded in $L^2(J;(H^2(\Omega))^\prime)$. Since
$(\phi_f^\epsilon(\alpha^2+\beta^2)\chi_f^\epsilon
f_1^\epsilon+\phi^\epsilon\chi_m^\epsilon m_1^\epsilon)$ is
uniformly bounded in $L^\infty(\Omega \times J)$, a standard
argument of Aubin's type  proves that
$(\phi_f^\epsilon(\alpha^2+\beta^2)\chi_f^\epsilon
f_1^\epsilon+\phi^\epsilon\chi_m^\epsilon m_1^\epsilon)$ lies in a
compact subset of  $L^2(J;(H^1(\Omega))^\prime)$. Therefore, there
is $\xi \in L^2(J;(H^1(\Omega))^\prime)$, such that, up to an
extracted subsequence,
$$
\phi_f^\epsilon(\alpha^2+\beta^2)\chi_f^\epsilon f_1^\epsilon
+\phi^\epsilon\chi_m^\epsilon m_1^\epsilon \to \xi \quad \text{in }
   L^2(J;(H^1(\Omega))^\prime)   \text{ as }  \epsilon \to 0.
$$
Two-scale convergence arguments show that
\[
\xi=(\alpha^2+\beta^2) (\int_{Y_f}\phi_f (y) dy) f_1
+ (\int_{Y_m}\phi (y) dy) m_1,
\]
 where $m_1=\alpha c_1+\beta C_1=f_1$
 by Lemma \ref{lem6}.

On the other hand, the sequence
$(\chi_f^\epsilon f_1^\epsilon+\chi_m^\epsilon m_1^\epsilon)$ is uniformly
 bounded in space $L^2(J;H^1(\Omega))$.
We thus can pass to the limit in the product
$ \langle \phi_f(\alpha^2+\beta^2)\chi_f^\epsilon f_1^\epsilon
+\phi^\epsilon\chi_m^\epsilon m_1^\epsilon, \chi_f^\epsilon f_1^\epsilon
+\chi_m^\epsilon m_1^\epsilon \rangle_{(H^1(\Omega))' \times H^1(\Omega)}$
 as follows.
\begin{align*}
& \lim_{\epsilon \to 0} \Bigl(
\bigl\langle \chi_f^\epsilon (\alpha^2+\beta^2) \phi_f^\epsilon f_1^\epsilon + \chi_m^\epsilon \phi^\epsilon m_1^\epsilon,  \chi_f^\epsilon f_1^\epsilon \bigr\rangle
+ \bigl\langle \chi_f^\epsilon (\alpha^2+\beta^2) \phi_f^\epsilon f_1^\epsilon
\\
&+ \chi_m^\epsilon \phi^\epsilon m_1^\epsilon, \chi_m^\epsilon m_1^\epsilon \bigr\rangle
\Bigr)
=\bigl\langle \bigl( (\alpha^2+\beta^2) \int_{Y_f}\phi_f (y) dy  + \int_{Y_m}\phi (y) dy \bigr) f_1, \vert Y_f\vert f_1 \bigr\rangle
\\
&+ \bigl\langle \bigl( (\alpha^2+\beta^2) \int_{Y_f}\phi_f (y) dy  + \int_{Y_m}\phi (y) dy \bigr) f_1, \vert Y_m\vert m_1 \bigr\rangle
\\
&=\bigl\langle \bigl( (\alpha^2+\beta^2) \int_{Y_f}\phi_f (y) dy  + \int_{Y_m}\phi (y) dy \bigr) f_1, f_1 \bigr\rangle .
\end{align*}
As a consequence we have
\begin{align*}
& \lim_{\epsilon \to 0} \bigl\langle \bigl( (\alpha^2+\beta^2)\chi_f^\epsilon \phi_f^\epsilon + \chi_m^\epsilon \phi^\epsilon \bigr) (\chi_f^\epsilon  f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon -f_1 ),
\chi_f^\epsilon  f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon -f_1 \bigr\rangle
\\
&= \lim_{\epsilon \to 0} \Bigl(
  \bigl\langle \bigl( (\alpha^2+\beta^2)\chi_f^\epsilon \phi_f^\epsilon + \chi_m^\epsilon \phi^\epsilon \bigr) (\chi_f^\epsilon  f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon  ),
\chi_f^\epsilon  f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon  \bigr\rangle
\\
& \quad -2 \bigl\langle \bigl( (\alpha^2+\beta^2)\chi_f^\epsilon \phi_f^\epsilon + \chi_m^\epsilon \phi^\epsilon \bigr) (\chi_f^\epsilon  f_1^\epsilon + \chi_m^\epsilon m_1^\epsilon  ),
f_1  \bigr\rangle
\\
& \quad + \bigl\langle \bigl( (\alpha^2+\beta^2)\chi_f^\epsilon \phi_f^\epsilon + \chi_m^\epsilon \phi^\epsilon \bigr) f_1, f_1  \bigr\rangle
\Bigr)
= 0 .
\end{align*}
Since $\alpha^2+\beta^2 >0$, $\phi_f,\phi^\epsilon \ge \phi_->0$,
this shows that $(\chi_f^\epsilon  f_1^\epsilon
+ \chi_m^\epsilon m_1^\epsilon )$ strongly converges to $f_1$
in $L^2(\Omega \times J)$.
A similar calculation using   (\ref{2.45}) gives the result for
$\chi_m^\epsilon(c_1^\epsilon+c_2^\epsilon)$.
The proof of the Lemma is complete.
 \end{proof}

We  now have the tools to pass to the limit in the pressure
equation. The structure of the problem satisfied by
$\theta^\epsilon$  is  similar with the one of the pressure
problem in the fractured part in the totally fractured setting
(Section 4). We thus do not detail the details of the convergence.
We claim that the limit homogenized pressure problem in $\Omega
\times J$ is the following.
\begin{gather*}
 \overline{\phi} \, \partial_t p
-\mathop{\rm div} \Bigl( \int_Y \Bigl( \chi_f(y)
\frac{k_f(y)}{\mu(f_1)} + \chi_m(y) \alpha (c_1
+c_2)\frac{k(y)}{\mu(m_1)} \Bigr) (\nabla p + \nabla_y p^1) \, dy
\Bigr) =q ,
\\
 - \mathop{\rm div}_y \Bigl( \Bigl( \chi_f(y)
\frac{k_f(y)}{\mu(f_1)} + \chi_m(y) \alpha (c_1+c_2)
\frac{k(y)}{\mu(m_1)} \Bigr) (\nabla p + \nabla_y p^1)  \Bigr) =0
\quad \text{in}\  Y ,
\\
 \Bigl( \chi_f(y) \frac{k_f(y)}{\mu(f_1)} +
\chi_m(y) \alpha (c_1+c_2) \frac{k(y)}{\mu(m_1)} \Bigr) (\nabla p
+ \nabla_y p^1) \cdot \nu =0 \quad \text{on}\  \Gamma \times J,
\end{gather*}
where $\overline{\phi}$ is defined by
\begin{equation}
\overline{\phi} = \int_{Y_f} \phi_f(y) \, dy + \int_{Y_m} \phi(y) \, dy
= \overline{\phi_f}^{Y_f} + \overline{\phi}^{Y_m} .
\label{5.1}
\end{equation}
Let $v^j(x,y,t)$ be the $Y$-periodic solution of the cell-problem
\begin{equation}
\begin{gathered}
 -\mathop{\rm div}_y \bigl( ( \chi_f(y) k_f(y) + \chi_m(y) \alpha
(c_1+c_2) k(y) ) (\nabla_y v^j+e^j)  \bigr)=0 \quad  \text{in } Y,
\\
 \int_Y v^j\, dy =0.
\end{gathered}\label{5.2}
\end{equation}
Defining an homogenized tensor of permeability
$\overline{K}=(\overline{K}_{ij})_{1\le {i,j} \le 2}$ in $\Omega$
by
\begin{equation}
\begin{aligned}
\overline{K}_{ij} (x,t)
&= \int_Y \big( \chi_f(y) k_f(y)
+ \chi_m(y) \alpha (c_1+c_2) k(y) \big)\\
&\quad\times (\nabla_y v^i(x,y,t)+e^i)(\nabla_y v^j(x,y,t)+e^j)\,dy,
\end{aligned}\label{5.3}
\end{equation}
and setting
$$ p^1(x,t,y)=\sum_{j=1}^{2} v^j(x,y,t)\, \partial_{j} p (x,t),$$
we claim  the following result.

\begin{proposition} \label{prop4}
 The homogenized pressure problem is
 \begin{gather}
 \overline{\phi}\, \partial_t p + \mathop{\rm div} \underline{v}
=q, \quad
 \underline{v} =- \frac{1}{\mu(f_1)}  \overline{K} \nabla  p  \quad
  \text{in }   \Omega \times J, \label{5.4} \\
 \underline{v} \cdot \nu \big|_{\Gamma \times J}= 0, \quad
p\big|_{t=0} =p^o,
\label{5.5}
\end{gather}
 where $\overline{\phi}$ (respectively $\overline{K}$) is the
homogenized porosity (respectively permeability tensor) given
in \eqref{5.1} (respectively in \eqref{5.3}).
\end{proposition}

As in Section 4,  we have to introduce a corrector for the Darcy velocity
in order to pass to the limit in the concentration equation.
We denote
\begin{gather*}
 \underline{v}_o(x,t,y)=- ( \chi_f(y) \frac{k_f(y)}{\mu(f_1)}
+ \chi_m(y) \alpha (c_1+c_2) \frac{k(y)}{\mu(f_1)}  )(\nabla p+\nabla_y p^1),
\\
\underline{v}_o^\varepsilon(x,t)=\underline{v}_o(x,t,x / \varepsilon),
 \quad (x,t) \in \Omega \times J .
\end{gather*}
Following the lines of Subsection 4.4, we can prove the following
corrector result.

\begin{lemma} \label{lem8}
 We have, for a subsequence,
$$
\int_{\Omega \times J} \vert \chi_f^\epsilon \underline{v}_f^\epsilon
+ \chi_m^\epsilon \underline{\mathcal{V}}^\epsilon
-\underline{v}^\varepsilon_o\vert^2 \,dx\,dt \to
0 \quad \text{as }   \varepsilon \to 0.
$$
\end{lemma}

Let us now turn to the concentration problem.
Once again, the technical difficulties are not greater than the ones
of Subsection 4.4 in the fractured part.
We thus only give the outline of the convergence study.
At the limit, we need to characterize the behavior of two global
concentrations, $f_1=m_1$ and $(c_1+c_2)$.

We begin with $f_1$. We multiply
$(\alpha^2+\beta^2)\eqref{2.15}+(\ref{2.40})$ by a test function
$\psi(x,t)+ \epsilon \psi_1(x,x/\epsilon,t)$, with $\psi \in
\mathcal{D}(\Omega \times J)$, $\psi_1 \in \mathcal{D}(\Omega
\times J;\mathcal{C}^\infty_{per}(Y))$. The terms on $\Gamma_{fm}$
are equal to zero when  we integrate by parts. Noting that
$\chi_f^\epsilon (\alpha^2+\beta^2) k_f^\epsilon \nabla
p_f^\epsilon + \chi_m^\epsilon k^\epsilon \nabla p^\epsilon =
(\alpha^2+\beta^2)(\chi_f^\epsilon k_f^\epsilon \nabla
p_f^\epsilon + \chi_m^\epsilon k^\epsilon \nabla p^\epsilon) +
(1-\alpha^2-\beta^2)\chi_m^\epsilon k^\epsilon \nabla p^\epsilon$
and passing to the limit, we get the following homogenized problem
in $\Omega \times J$.
\begin{gather}
 \begin{aligned}
&\bigl( (\alpha^2+\beta^2)\overline{\phi_f}^{Y_f}
+ \overline{\phi}^{Y_m} \bigr) \,\partial_t f_1
- \frac{(\alpha^2+\beta^2)}{\mu(f_1)}  \overline{K}  \nabla p \cdot \nabla f_1
 \\
& + (1-\alpha^2-\beta^2) \int_{Y_m} \chi_m(y)
\underline{v}_o(x,t,y) \cdot (\nabla f_1+\nabla_yf_1^1) \, dy
- \mathop{\rm div} (\overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1)\\
&=\bigl( (\alpha^2+\beta^2) \vert Y_f \vert + \vert Y_m\vert \bigr) q\,
(\hat{f_1}-f_1)   ,\label{5.6}
\end{aligned} \\
\overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1 \cdot
\nu \big|_{\Gamma \times J} =0 ,
 \quad f_1\big|_{t=0}=f_1^o .
 \label{5.7}
 \end{gather}
The function $f_1^1$ and  the  tensor of diffusion
$\overline{\mathcal{D}}=(\overline{\mathcal{D}}_{ij})_{1\le {i,j} \le 2}$
are defined by
\begin{gather}
 f_1^1(x,t,y) = \sum_{j=1}^2 w^j(x,t,y) \, \partial_j f_1,
 \label{5.8} \\
\begin{aligned}
&\overline{\mathcal{D}}_{ij}(\nabla p,\mu(f_1))\\
&=\int_{Y} \bigl( (\alpha^2+\beta^2) \chi_f(y) + \chi_m(y) \bigr)
 \mathcal{D}(\underline{v}_o )\left(\nabla_y   w^i(x,y)+e^i\right)
 \cdot \big(\nabla_y w^j(x,y)+e^j\big)dy,
\end{aligned}\label{5.9}
\end{gather}
where $w^j(x,t,y)$ is the $Y$-periodic solution of the following
cell-problem  for $(x,t) \in \Omega \times J$:
\begin{equation}
\begin{gathered}
 -\mathop{\rm div}_y \bigl(\bigl( (\alpha^2+\beta^2) \chi_f(y)
+ \chi_m(y) \bigr) \mathcal{D}(\underline{v}_o)(\nabla_y w^j+e^j )\bigr)=0
\quad \text{in }   Y,  \\
 \int_Y w^j(x,t,y)\,dy=0,\quad j=1, 2.
\end{gathered} \label{5.10}
 \end{equation}
 For the asymptotic study of $(c_1^\epsilon+c_2^\epsilon)$ we exploit
Problem (\ref{2.45})-(\ref{2.47}).
It leads to the following homogenized problem
 \begin{gather}
\begin{aligned}
& \overline{\phi}^{Y_m}  \,\partial_t (c_1+c_2)
+ \int_{Y_m} \chi_m(y) \underline{v}_o(x,t,y) \cdot ( \nabla (c_1+c_2)
+  \nabla_y (c_1^1+c_2^1) ) \, dy \\
& - \mathop{\rm div} (\overline{\mathcal{D}}_m(\nabla p,\mu(f_1))
\nabla (c_1+c_2)) \\
&= \vert Y_m\vert  q\, (1-c_1-c_2)  \quad \text{in }  \Omega \times J,
\end{aligned}\label{5.11}
\\
\overline{\mathcal{D}}_m \nabla (c_1+c_2) \cdot \nu \big|_{\Gamma \times J}=0    ,
 \quad (c_1+c_2)\big|_{t=0}=c_1^o+c_2^o .
 \label{5.12}
 \end{gather}
The functions $c_1^1$, $c_2^1$ and  the homogenized tensor of diffusion
$\overline{D}_m=(\overline{D}_{m_{ij}})_{1\le {i,j} \le 2}$ are defined by
 \begin{gather}
 (c_1^1+c_2^1)(x,t,y) = \sum_{j=1}^2 w_m^j(x,t,y) \, \partial_j (c_1+c_2),
 \label{5.13} \\
 \overline{\mathcal{D}}_{m_{ij}}(\nabla p,\mu(f_1))=\int_{Y_m}  \mathcal{D}(\underline{v}_o ) (\nabla_y   w_m^i+e^i ) \cdot (\nabla_y w_m^j+e^j ) \, dy,
\label{5.14}
\end{gather}
where $w_m^j(x,t,y)$ is the $Y$-periodic solution of the following
cell-problem  for $(x,t) \in \Omega \times J$:
\begin{equation}
\begin{gathered}
 -\mathop{\rm div}_y (\mathcal{D}(\underline{v}_o)(x,t,y)
(\nabla_y w_m^j(x,t,y)+e^j ))=0 \quad  \text{in }  Y_m,  \\
  \mathcal{D}(\underline{v}_o) (\nabla_y w_m^j+e^j )) \cdot \nu_y =0 \quad
  \text{on }  \Gamma_{fm},\quad j=1, 2.
 \end{gathered}
 \label{5.15}
 \end{equation}
We summarize  the results of the present section in the following theorem.

\begin{theorem} \label{thm3}
 The homogenized partially fractured model in $\Omega \times J$ is:
 \begin{gather*}
(\overline{\phi_f}^{Y_f} + \overline{\phi}^{Y_m}) \, \partial_t
p + \mathop{\rm div} \underline{v} =q, \quad
 \underline{v} =- \frac{1}{\mu(f_1)}  \overline{K} \nabla  p  ,
\\
 \underline{v} \cdot \nu \big|_{\Gamma \times J}= 0, \quad
p\big|_{t=0} =p^o,
\\
 p^1(x,t,y)=\sum_{j=1}^{2} v^j(x,y,t)\, \partial_{j} p (x,t),
\\
 \underline{v}_o(x,t,y)=- ( \chi_f(y) \frac{k_f(y)}{\mu(f_1)}  + \chi_m(y) \alpha (c_1+c_2) \frac{k(y)}{\mu(f_1)}  )(\nabla p+\nabla_y p^1),
\\
\begin{aligned}
& \bigl( (\alpha^2+\beta^2)\overline{\phi_f}^{Y_f} + \overline{\phi}^{Y_m} \bigr) \,\partial_t f_1
+ (\alpha^2+\beta^2)  \underline{v} \cdot \nabla f_1\\
&+ (1-\alpha^2-\beta^2) \int_{Y_m} \chi_m(y)
\underline{v}_o(x,t,y) \cdot (\nabla f_1+\nabla_yf_1^1) \, dy
 - \mathop{\rm div} (\overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1)\\
&=\bigl( (\alpha^2+\beta^2) \vert Y_f \vert + \vert Y_m\vert \bigr) q\,
(\hat{f_1}-f_1) ,
\end{aligned}\\
f_1^1(x,t,y) = \sum_{j=1}^2 w^j(x,t,y) \, \partial_j f_1,
 \\
\overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1 \cdot \nu
\big|_{\Gamma \times J} =0 ,
 \quad f_1\big|_{t=0}=f_1^o,
 \\
\begin{aligned}
&\overline{\phi}^{Y_m}  \,\partial_t (c_1+c_2)
+ \int_{Y_m} \chi_m(y) \underline{v}_o(x,t,y) \cdot ( \nabla (c_1+c_2)
+  \nabla_y (c_1^1+c_2^1) ) \, dy \\
& - \mathop{\rm div} (\overline{\mathcal{D}}_m(\nabla p,\mu(f_1))
\nabla (c_1+c_2)) \\
&= \vert Y_m\vert  q\, (1-c_1-c_2)  ,
\end{aligned}\\
  (c_1^1+c_2^1)(x,t,y) = \sum_{j=1}^2 w_m^j(x,t,y) \, \partial_j (c_1+c_2),
\\
 \overline{\mathcal{D}}_m \nabla (c_1+c_2) \cdot \nu \big|_{\Gamma \times J} =0 ,
 \quad (c_1+c_2)\big|_{t=0}=c_1^o+c_2^o .
\end{gather*}
The auxiliary functions $v^j$, $w^j$, $w^j_m$ are defined by the cell
problems \eqref{5.2}, \eqref{5.10}, \eqref{5.15}.
The homogenized permeability is given in \eqref{5.3}, while the
diffusion tensors $\overline{\mathcal{D}}$ and
$\overline{\mathcal{D}}_m$ are given in \eqref{5.9} and \eqref{5.14}.
\end{theorem}

\begin{remark} \label{rmk2} \rm
The latter model is consistent with the one corresponding to a
non-fractured porous medium.
Indeed, letting $\alpha \to 1$ and $\vert Y_m \vert \to 0$, one gets
 \begin{gather*}
 \overline{\phi}^{Y}   \partial_t p + \mathop{\rm div}
\underline{v} =q, \quad
 \underline{v} =- \frac{1}{\mu(f_1)}  \overline{K} \nabla  p  ,
\\
 \underline{v} \cdot \nu \big|_{\Gamma \times J}= 0, \quad
p\big|_{t=0} =p^o,
\\
  p^1(x,t,y)=\sum_{j=1}^{2} v^j(x,y,t)\, \partial_{j} p (x,t),
\quad
 \underline{v}_o(x,t,y)=-  \frac{k(y)}{\mu(f_1)}  (\nabla p+\nabla_y p^1),
\\
 \overline{\phi}^{Y} \partial_t f_1
+  \underline{v} \cdot \nabla f_1
 - \mathop{\rm div} (\overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1)
=  q\, (\hat{f_1}-f_1) ,
 \\
 \overline{\mathcal{D}}(\nabla p,\mu(f_1)) \nabla f_1 \cdot \nu \big|_{\Gamma \times J} =0 ,
 \quad f_1\big|_{t=0}=f_1^o.
\end{gather*}
Note that this model is the homogenized form corresponding to the
following microscopic model in a non-fractured porous medium
(see \cite{Cho1} for a rigorous derivation):
\begin{gather*}
 \phi^\epsilon \partial_t p^\epsilon +\mathop{\rm
div}(\underline{v}^\epsilon) = q, \quad
 \underline{v}^\epsilon = -
\frac{k^\epsilon}{\mu(f_1^\epsilon)} \nabla p^\epsilon \quad
\text{in}\  \Omega \times J,
\\
 \phi^\epsilon  \partial_t f_1^\epsilon + \underline{v}^\epsilon
\cdot \nabla f_1^\epsilon - \mathop{\rm
div}(\mathcal{D}(\underline{v}^\epsilon ) \nabla f_1^\epsilon ) =
q (\hat{f}_1-f_1^\epsilon) \quad \text{in}\  \Omega  \times J,
 \\
\underline{v}^\epsilon \cdot \nu \big|_{\Gamma \times J}=0, \quad
\mathcal{D}(\underline{v}^\epsilon ) \nabla f_1^\epsilon \cdot \nu \big|_{\Gamma \times J}=0 ,
\quad
p_f^\epsilon \big|_{t=0}=  p^o, \quad
f_1^\epsilon \big|_{t=0}=f_1^o .
\end{gather*}
\end{remark}


\begin{remark} \label{rmk3} \rm
Let us add some ``physical'' comments on the homogenized  model of
miscible displacement in partially fractured media described in
Theorem \ref{thm3}. Contrary to the  one-component flow considered
in \cite{Sho}, it seems that the double porosity characteristics
(that is a system of equations in the homogenized fracture coupled
with one in a microscopic matrix block) disappear as soon as a
direct flow occurs in the matrix part, that is as soon as $\alpha
>0$. This corresponds to some experimental data (\cite{And} and
references therein). However, the model of Theorem \ref{thm3} is
also not a single porosity model as the one described in Remark
\ref{rmk2}. It really contains some matrix effects. We recall in
particular that the smaller $\alpha$ is, the more important is the
storage in the matrix. And in the model of Theorem \ref{thm3}, the
homogenized permeability is concentration ($\alpha(c_1+c_2)$)
dependent, and the  convective effects in the homogenized fracture
depend on $\alpha $. One can compare this effects with some models
where the permeability is concentration dependent:  propagation in
clays (see \cite{Kac} and the references therein)  or blood flow
in micro vessels (see \cite{Sug} and the references therein) for
instance.
\end{remark}

\begin{thebibliography}{00}

\bibitem{Acer}
E.~Acerbi, V.~Chiad\`o {P}iat, G.~Dal Maso, and D.~Percivale.
\newblock An extension theorem from connected sets, and homogenization in
  general periodic domains.
\newblock {\em Nonlinear Anal.}, 18:696--496, 1992.

\bibitem{A}
G.~Allaire.
\newblock Homogenization and two-scale convergence.
\newblock {\em SIAM J. Math. Anal.}, 23:1482--1518, 1992.

\bibitem{And}
R. E. Andrews, D. R. Wunsch, and J. S. Dinger.
\newblock Evaluation of the use of fracture-flow solutions to analyze aquifer
  test data collected from wells in the eastern kentucky coal field.
\newblock In National Ground~Water Association, editor, {\em Fractured-Rock
  Aquifers 2002}, pages 119--123, 2002.

\bibitem{Dou}
T.~Arbogast, J.~Douglas {J}r., and U.~Hornung.
\newblock Derivation of the double porosity model of single phase flow via
  homogenization theory.
\newblock {\em SIAM J. Appl. Math. Anal.}, 21:823--836, 1990.

\bibitem{Bar}
G. I. Barenblatt, I. P. Zheltov, and I. N. Kochina.
\newblock Basic concepts in the theory of seepage of homogeneous liquids in
  fissured rocks.
\newblock {\em Prikl. Mat. Mech.}, 24(5):1286--1303, 1960.

\bibitem{B}
J.~Bear.
\newblock {\em Dynamics of Fluids in Porous Media}.
\newblock American Elsevier, 1972.

\bibitem{Adl}
S.~Bekri, K.~Xu, F.~Yousefian, P. M. Adler, J.-F. Thovert,
J.~Muller, K.~Iden,
  A.~Psyllos, A. K. Stubos, and M. A. Ioannidis.
\newblock Pore geometry and transport properties in {N}orth {S}ea {C}halk.
\newblock {\em J. Petroleum Sci. and Engineering}, 25(3-4):107--134, 2000.

\bibitem{BLM}
A.~Bourgeat, S.~Luckhaus, and A.~Mikeli\'c.
\newblock Convergence of the homogenization process for a double-porosity model
  of immiscible two-phase flow.
\newblock {\em SIAM J. Math. Anal.}, 27:1520--1543, 1996.

\bibitem{Cho1}
C.~Choquet.
\newblock Homogenization of a one-dimensional model for compressible miscible
  flow in porous media.
\newblock {\em Boll. Un. Mat. Ital. B}, 6:399--414, 2003.

\bibitem{Cho2}
C.~Choquet.
\newblock Derivation of the double porosity model of a compressible miscible
  displacement in naturally fractured reservoirs.
\newblock {\em Appl. Anal.}, 83:477--500, 2004.

\bibitem{Cho3}
C.~Choquet.
\newblock Existence result for a radionuclide transport model with unbounded
  viscosity.
\newblock {\em J. Math. Fluid Mech.}, 6:365--388, 2004.

\bibitem{Cio}
D.~Cioranescu, A.~Damlamian, and G.~Griso.
\newblock Periodic unfolding and homogenization.
\newblock {\em C. R. Math. Acad. Sci. Paris}, 335:99--104, 2002.

\bibitem{Sho}
G. W. Clark and R. E. Showalter.
\newblock Two-scale convergence of a model for flow in a partially fissured
  medium.
\newblock {\em Electron. J. Differential Equations}, 1999(2):1--20, 1999.

\bibitem{DS}
J.~Douglas {J}r., M.~Peszy\'nska, and R.E. Showalter.
\newblock Single phase flow in partially fissured media.
\newblock {\em Transp. Porous Media}, 28:285--306, 1997.

\bibitem{DR}
J.~Douglas {J}r. and J. E. Roberts.
\newblock Numerical methods for a model for compressible miscible displacement
  in porous media.
\newblock {\em Math. of Computation}, 41:441--459, 1983.

\bibitem{Kac}
M.~Kaczmarek, T.~Hueckel, V.~Chawla, and P.~Imperiali.
\newblock Transport through a clay barrier with the contaminant dependent
  permeability.
\newblock {\em Transp. Porous Media}, 29:159--178, 1997.

\bibitem{Ko}
E. J. Koval.
\newblock A method for predicting the performance of unstable miscible
  displacements in heterogeneous media.
\newblock {\em SPEJ trans. AIME}, 228:145--154, 1963.

\bibitem{Ng}
G.~Nguetseng.
\newblock A general convergence result for a functional related to the theory
  of homogenization.
\newblock {\em SIAM J. Math. Anal.}, 20:608--623, 1989.

\bibitem{Pe}
D. W. Peaceman.
\newblock {\em Fundamentals of Numerical Reservoir Simulation}.
\newblock Elsevier, 1977.

\bibitem{Sug}
M.~Sugihara Seki and B. M. Fu.
\newblock Blood flow and permeability in microvessels.
\newblock {\em Fluid Dynam. Res.}, 37:82--132, 2005.

\end{thebibliography}






\end{document}
