CONVERGENCE OF ADAPTIVE MIXED FEMS FOR GENERAL SECOND ORDER LINEAR ELLIPTIC PROBLEMS\thanks{This work was %supported by the Society for Industrial and Applied %Mathematics
Asha K. Dond\thanks{Department of Mathematics
Indian Institute of Technology
Bombay ({\tt asha@math.iitb.ac.in
Asha K. Dond
Neela Nataraj and Amiya K. Pani
The convergence of an adaptive mixed finite element method for general second order linear elliptic problems defined on simply connected bounded polygonal domains is analyzed in this paper. The main difficulties in the analysis are posed by the non-symmetric and indefinite form of the problem along with the lack of the orthogonality property in mixed finite element methods. The important tools in the analysis are {a~posteriori error estimators, quasi-orthogonality property and quasi-discrete reliability established using representation formula for the lowest-order Raviart-Thomas solution in terms of the Crouzeix-Raviart solution of the problem. An adaptive marking in each step for the local refinement is based on the edge residual and volume residual terms of the {a~posteriori estimator. Numerical experiments confirm the theoretical analysis.
Amfem_cgs.tex
%31th Jan 2014
%24th Jan
%18 jan 2014
%\documentclass[final,leqno,onefignum,onetabnum]{siamltex1213}
%\documentclass[final]{siamltex}
\documentclass[final,leqno]{siamltex}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{mathrsfs}
\usepackage{amssymb}
%\usepackage[colorlinks,pdfmark,dvips]{hyperref}
\usepackage{pgfplots}
\usepackage{tikz}
%\usepackage{epsfig}
%\pgfplotsset{compat=1.5}
%\usepackage{subfig}
\usepackage{pgfplots}
\pgfplotsset{compat=1.5}
\usepackage{subfig}
\usepackage{graphicx}
\usepackage{bm}
\usepackage{longtable}
\newtheorem{example}{}[section]
\newtheorem{remark}{\bf Remark}[section]
\newcommand{\norm}[1]{\lVert #1 \rVert}
\newcommand{\tnorm}[1]{\ensuremath{\left| \! \left| \! \left|} #1 \ensuremath{\right| \! \right| \! \right|}}
\newcommand{\mc}[1]{\mathcal #1}
\newcommand{\e}[1]{\epsilon #1}
\newcommand{\T}{\mathcal T}
\newcommand{\J}{\jmath}
%\def{\cT}{\mathcal{T}}
\title{CONVERGENCE OF AN ADAPTIVE MIXED FINITE ELEMENT METHOD FOR GENERAL SECOND ORDER LINEAR ELLIPTIC PROBLEMS}

%\title{CONVERGENCE OF ADAPTIVE MIXED FEMS FOR GENERAL SECOND ORDER LINEAR ELLIPTIC PROBLEMS\thanks{This work was
%supported by the Society for Industrial and Applied
%Mathematics}} 
\author{Asha K. Dond\thanks{Department of Mathematics, Indian Institute of Technology, Bombay
        ({\tt asha@math.iitb.ac.in}).}
       \and Neela Nataraj \thanks{Department of Mathematics, Indian Institute of Technology, Bombay
                ({\tt neela@math.iitb.ac.in})}
                \and Amiya Kumar Pani \thanks{Department of Mathematics, Indian Institute of Technology, Bombay
                                ({\tt akp@math.iitb.ac.in})}
                }
%\author{Asha K. Dond,\and Neela Nataraj  and Amiya K. Pani }
%\thanks{Society for Industrial and
%Applied Mathematics, Philadelphia, Pennsylvania. 
%(\email{tex@siam.org}). Questions, comments, or corrections
%to this document may be directed to that email address.}}

\begin{document}
\maketitle
\slugger{mms}{xxxx}{xx}{x}{x--x}%slugger should be set to mms, siap, sicomp, sicon, sidma, sima, simax, sinum, siopt, sisc, or sirev

\begin{abstract}
The convergence  of an adaptive mixed finite element method for general second order linear elliptic problems defined on simply connected bounded polygonal domains is analyzed in this paper. The main difficulties in the analysis are posed by the non-symmetric and  indefinite form of the problem along with the lack of the orthogonality property in mixed finite element methods. The important tools in the  analysis are {\it a~posteriori} error estimators, quasi-orthogonality property and quasi-discrete reliability established using  representation formula for the lowest-order Raviart-Thomas solution  in terms of the Crouzeix-Raviart solution of the problem. An adaptive marking in each step for 
the local refinement is based on  the edge residual  and  volume residual terms of the {\it a~posteriori} estimator. 
Numerical experiments confirm the theoretical analysis.
\end{abstract}
%
\begin{keywords}
{Adaptive mixed finite element method, {\it a~posteriori} error estimator,  contraction property, convergence and quasi-optimality.}
\end{keywords}

\begin{AMS}
65N30,65N50
\end{AMS}


\pagestyle{myheadings}
\thispagestyle{plain}
%\markboth{TEX PRODUCTION}{USING SIAM'S \LaTeX\ MACROS}
\markboth{Convergence of AMFEM}{}
\section{Introduction}
The general second-order linear elliptic PDE  
on a  simply connected bounded  polygonal Lipschitz domain  $\Omega \subset {\mathbb R}^2$ with boundary $\partial \Omega$ reads: For given right-hand side $f\in L^2(\Omega)$, seek $u$ such that
\begin{eqnarray}\label{eq1}
 \mathcal{L}u :=  -\nabla \cdot (\mathbf A \nabla u+u {\mathbf b})+ \gamma ~ u=f  \hspace{.5cm}\mbox {in}~\Omega, ~~~u=0~\hspace{.5cm}\mbox{on}~\partial\Omega.
\end{eqnarray}
The coefficients ${\bf A,b},\gamma$ are all piecewise smooth and the symmetric matrix ${\mathbf A}$ is positive definite and uniformly  bounded away from zero.\\
The flux variable ${\bf p}  =-({\mathbf A\nabla u+u {\mathbf b}})$ and ${\mathbf b}^*= \textbf A^{-1}{\bf b}$ 
 allow to  recast  (\ref{eq1}) as a  first-order system  
\begin{eqnarray}\label{eq3}
\begin{array}{lll} 
\mathbf A^{-1}{\bf p} +u {\mathbf b}^* + \nabla u =0  \; ~{\text{and}}~~
{\text{div}}~{\bf p}+\gamma \: u=f \;  {\rm in } \; \Omega. %{\rm in } \; \Omega, \\ %\bar{\mathbf b})=0\
\end{array}
\end{eqnarray}
The mixed formulation  %corresponding to (\ref{eq1}) is  to
seeks  $({\bf p},u)\in H(\text{div},\Omega) \times L^2(\Omega) $ such that 
\begin{eqnarray} \label{eq4}
\begin{array}{llll} 
(\textbf A^{-1}{ \bf p} + u {\bf b}^*,{\bf q})-({\rm div}~ {\bf q}, u)=0 \qquad \mbox {for all}\; {{\bf q}} \in H(\text{div},\Omega),\\
({\text{div}}~ {\bf p}, v )+(\gamma \: u,v)=(f,v) \qquad \qquad \mbox{for all}\; v \in L^2(\Omega).
\end{array}
\end{eqnarray}
Here and throughout the paper, $H(\text{div},\Omega)=\{{\bf q}\in L^2(\Omega; \mathbb R^2) : \: \text{div}~ {\bf q} \in L^2(\Omega) \} $ and 
$L^2(\Omega; \mathbb R^2)$ denotes the space of $\mathbb R^2$-valued  $L^2$  functions defined over the domain $\Omega$.
The existence and uniqueness of the mixed solution for elliptic problems have been proved in \cite{Brezzi,CC}. 

 The study of analysis of the adaptive finite element methods (AFEM) is an essential component  of the adaptive process. Various  {\it a~posteriori}  error estimators are  reviewed in \cite{oden}, and the references therein.
 The marking strategies, %interior node property, a concept of data oscillation,
 convergence and optimality 
 are well established for  the adaptive conforming finite element methods in literature \cite{Car_apo1,quasi, dof, morin1, morin2,  veeser, steven}. 
For the Poisson problem, the convergence and optimality   %and for  general second order linear elliptic PDE 
 have been established for the adaptive nonconforming FEM  \cite{bek_non,CH,car_non_cgt} % and   \cite{chen_hoppe} respectively. 
and for the adaptive mixed FEM   \cite{bek_mixed,car_cgt2,car_cgt1,cascon,chen_cgs}.  
The recent article `Axioms of adaptivity' {\cite{axiom}} provides a general framework  to optimality of adaptive schemes.

The non-symmetric and indefinite second order elliptic equations  with conforming, nonconforming mixed FEM have been discussed in various articles \cite{Brenner-Scott, Brezzi,CC,Chen,chen_hoppe,daglas,schatz1,schatz}. These articles discuss the existence and uniqueness of the solution with {\it a~priori}  error estimates. 
%of the non-symmetric and indefinite second order elliptic problems.
 {\it A posteriori}  error estimates and its  convergence for conforming FEM for  general second order 
 linear elliptic PDEs have been achieved using contraction of the sum of energy error plus oscillation  
 in  \cite{nochetto} and  the quasi- optimality in \cite{prac}.  %by Mekchay and Nochetto
{\it A posteriori}  error estimates  and quasi-optimal convergence of the adaptive nonconforming FEM   
have been obtained in  
  \cite{chen_hoppe}. %for  general second order linear elliptic PDE% respectively.
 % **Schatz***
%The literature regarding
To the best of our knowledge, we have not come across any work which discusses the 
convergence and optimality of the adaptive mixed finite element method (AMFEM) for non-symmetric and indefinite elliptic problems. % The lack of orthogonality in mixed finite element method and the non-symmetric form of equation are  main difficulties  and thus the subject of interest in this problem. 
The main challenges, the lack of orthogonality in MFEM and the non-symmetric form of equation  are addressed in this work.
Also as the flux variable  ${\bf p}$ involves  $u$ explicitly,  
 the analysis of variable $u$  becomes inevitable for the analysis of the flux ${\bf p}$. %This paper focuses on it, in following way
In this paper, the main contributions are summarized as:
\begin{itemize}
\item for the adaptive algorithm, the  marking strategy in each step for 
the local refinement is proposed  based on the comparison of the edge residual term and the volume residual terms of the {\it a~posteriori} estimator,
\item {\it a~posteriori} error estimator, quasi-orthogonality property and 
quasi-discrete reliability results are derived with the help of the representation formula for the lowest-order Raviart-Thomas solution  in terms of 
the Crouzeix-Raviart solution of the problem,
\item the contraction property is shown for the linear combination of %he total error that is, 
the sum of errors in  ${\bf p}$ and $u$, 
the edge residual  estimator  and the volume residual estimators,% terms. 
%\item  The contraction propety is established by the reduction of the total error, edge residual estimator and the volume residual estimator terms between two sucessive (consecutive) adaptive loops. 
\item  the convergence and the quasi-optimality results are achieved, under the assumption of  small initial mesh-size $h_0$.
%\item Under assumption of sufficiently small initial mesh-size $h$, the convergence and the quasi-optimality property are established.
\end{itemize}
An outline of the  paper is as follows.
 Section 2 introduces  notations and the adaptive algorithm for the mixed finite element method.
Section 3 describes some auxiliary results necessary for the convergence analysis. 
The contraction property
 and the quasi-optimal convergence of the adaptive mixed finite element method  are established in Section 4.
The numerical experiments  are presented in Section 5. Appendix I summarizes the constants used in the article and their interdependencies.

Here are some notations used  throughout the paper. An inequality $A\lesssim B$ abbreviates $A\leq CB$, where $C>0$ is a mesh-size independent constant that depends only on the domain and 
the shape of finite elements;
 $A\approx B$ means $A \lesssim B \lesssim _{â€¢} A$. Standard notation applies to Lebesgue and Sobolev spaces and $\norm{\cdot}$ abbreviates
 $\norm{\cdot}_{L^2(\Omega)}$ with   $L^2$ scalar product $(\cdot,\cdot)$. For a vector ${\bf q}=(q_1,q_2)\in H({\rm div},~ \Omega)$, $\norm{{\bf q}}:=(\norm{q_1}^2+\norm{q_2}^2)^{1/2}$. Let $\norm{\cdot}_T $  and $\norm{\cdot}_E $ 
 denote  $\norm{\cdot}_{L^2( T)} $ and $\norm{\cdot}_{L^2(E)} $ respectively. $H^m(\Omega)$ denotes the  Sobolev  space of order $m$ with  norm given by $\norm{\cdot}_m.$  
 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\section{AMFEM algorithm}
%This section starts with notations and discusses the  adaptive algorithm. \\
This section discusses notations and  the  adaptive algorithm. \\
Let $\mathcal{T}_h$ be a regular triangulation  the domain 
${\Omega}\subset {\mathbb R}^2 $ into triangles such that 
$\cup_{T\in \mathcal{T}_h} T =\overline{\Omega}.$
Let $\mathcal{E}_h $ be the set of all edges in $\mathcal{T}_h$ and let ${\mathcal{E}_h}({\partial \Omega})$ be the 
set of all boundary edges in $\mathcal{T}_h$. % and let $\mathcal{N}$ denote the set of vertices in $ \mathcal{T}.$
Further, let $\text{mid}(E) $ denote the midpoint of the edge $E$ and $\text{mid}(T)$ denote the centroid of the triangle $T.$ 
 The set of edges of the element $T$ is denoted by $ {\mathcal{E}} (T)$ and let 
 $h_T:=$ diam$(T)$ for $T \in \mathcal{T}_h$. Define $h_{\mathcal T} \in P_0(\mathcal T_h)$, as a piecewise constant mesh-size function such that   
$h_{\mathcal T}|_T:= h_T$ for all $T \in \mathcal T_h$. Let $h:= \max_{T\in \mathcal {T}_h} h_T$ and  
$h_E$ be the length of the edge $E\in \mathcal{E}_h $. % with unit outward normal  $\nu_E.$\\
For any edge $E $, $\nu_E $ is  the unit normal vector exterior
to $T$ and $\tau_E$ is the unit tangential vector along $E$.
Let $\Pi_h$ be the $L^2$ projection onto  $P_0(\mathcal{T}_h)$ and define $osc_h(f):=\norm{h_{\mathcal T}(1-\Pi_h)f},$ where
$$
P_r(\mathcal{T}_h)=\{v \in L^2(\Omega): \forall T\in \mathcal{T}_h, v|_T\in P_r(T) \}.
$$
%Here and throughout this paper,
Here, $P_r(T)$ denotes the algebraic polynomials  of total degree at most $r \in \mathbb {N}$ as functions 
on the triangle $T\in \mathcal{T}_h.$
%$ The $P_1$ conforming finite element space reads 
% $$V(\mathcal{T}) := P_1(\mathcal{T}) \cap H^1_0(\Omega).$$ 
%
The jump of ${\bf q}$ across $E$ is denoted by $[{\bf q}]_E$; that is, for  two  neighboring triangles
$T_{+}$ and $T_{-},$
$$
[{\bf q}]_E(x) :=({\bf q}|_{T_{+}} (x)-{\bf q}|_{T_{-}}(x))~\text{  for } ~x \in E =\partial T_{+} \cap 
\partial T_{-}.
$$
The sign of $[{\bf q}]_E$ is defined using the convention that there is a fixed orientation of $\nu_E$ pointing 
outside of $T_{+}.$ % An union of $T_+$ and $T_-$ neighboring triangles,   patch of edge $E$ is denoted by $\omega_E$.
The patch $\omega_E$ denote the union of elements that share a common edge $E$. 
The piecewise gradient $\nabla_{NC}: H^1(\mathcal {T}_h)\longrightarrow L^2(\Omega; {\mathbb R}^2)$ acts as
$ \nabla_{NC}  {v}|_T = \nabla v|_T\;\;\;\text{for all} ~T \in \mathcal{T}_h.$ 
The broken Sobolev norm $\tnorm{\cdot}_{NC}$ abbreviates $(\mathbf{A} \nabla_{NC} ~\cdot,\nabla_{NC} ~\cdot)^{1/2}$.\\% based on  the  triangulation $\mathcal{T}_h.$ \\ 
The non-conforming Crouzeix-Raviart (CR) finite element space  with respect to the triangulation $\mathcal T_h$ reads
\begin{eqnarray*}
&& CR^1(\mathcal T_h): =\{v \in P_1(\mathcal{T}_h): v~ \text{ is continuous in all midpoints mid($E$) of edges} ~E \in \mathcal E_h\},\\
&& CR^1_0(\mathcal T_h ):=\{v \in CR^1(\mathcal T_h):v(\text{mid} (E))
=0 ~~\text{for all} ~E \in \mathcal E_h ({\partial\Omega}) \}.
\end{eqnarray*}
%Let $H(\text{div},\Omega)=\{{\bf q}\in L^2(\Omega; \mathbb R^2) : \: \text{div}~ {\bf q} \in L^2(\Omega) \} $, 
%with the norm $~\norm{\cdot}_{H(\text{div},\Omega)}$.
%and
The lowest-order Raviart-Thomas space with respect to the triangulation $\mathcal T_h$  reads 
%$RT_0(\mathcal{T})$ for discretization of ${\bf p}$ as
\begin{align*}
RT_0(\mathcal{T}_h):=\{{\bf q}\in H(\text{div},\Omega): & \: \forall T \in \mathcal{T}_h 
~\exists {\bf c} \in \mathbb{R}^2~ \exists d \in \mathbb{R}~ \; \forall {\bf x} \in T, ~{\bf q}({\bf x})={\bf c}+ d ~ {\bf x}\\
&\text{ and} \; ~\forall E\in \mathcal{E}_h, [{\bf q}]_{E} \cdot \nu_E=0\} .  
\end{align*}
%
\subsection{ Algorithm}
 The standard structure of an adaptive algorithm is successive loops 
 $$\mbox{ SOLVE $\rightarrow$ ESTIMATE $\rightarrow$ MARK $\rightarrow$ REFINE}$$
 on different levels of the triangulation.\\
 %on each level $\ell$ with respect to a computed triangulation $\mathcal_\ell$.\\
 In the step { SOLVE}, the discrete mixed finite element problem (RTFEM) for \eqref{eq4} defined by: seek $ ({\bf p_h}, u_h) \in RT_0(\mathcal {T}_h)\times P_0(\mathcal T_h)$ such that
\begin{eqnarray}
 (\mathbf{A}^{-1}{\bf p_h}+ u_h {\bf b}^*,{\bf q}_h)-(\text{div}~{\bf q}_h,u_h)=0 \qquad 
 \text{for all}\, {\bf q}_h\in RT_0(\mathcal {T}_h), \label{eqna1}\\
 (\text{div}~ {\bf p_h},v_h)+ (\gamma u_h, v_h)=(f_h,v_h) \qquad  \text {for all} \, v_h \in P_0(\mathcal{T}_h), \label{eqna2}
\end{eqnarray} 
is solved.  
Recall that, $f_h:=\Pi_hf $ is the $L^2$-projection of $f$ onto $P_0(\T_h)$. \\ %The matlab implementation of the lowest order Raviart-Thomas mixed finite element method is discussed in \cite{3matlab} for  Poisson problem.\\
 %In the step {\bf Solve}, the discrete mixed finite element problem (RTFEM) for \eqref{eq4} reads as: seek $ ({\bf p_H}, u_H) \in RT_0(\mathcal {T}_H)\times P_0(\mathcal T_H)$ such that
%\begin{eqnarray}
% (\mathbf{A}^{-1}{\bf p_H}+ u_H {\bf b}^*,{\bf q} _{RT})-(\text{div}~{\bf q}_{RT},u_H)=0 \qquad 
% \text{for all}\, {\bf q}_{RT}\in RT_0(\mathcal {T_H}), \label{eqna1}\\
% (\text{div}~ {\bf p_H},v_H)+ (\gamma u_H, v_H)=(f_H,v_H) \qquad  \text {for all} \, v_H \in P_0(\mathcal{T_H}). \label{eqna2}
%\end{eqnarray} 
% Here $f_H$ is the $L^2$-projection of $f$. 
The step {ESTIMATE} consists of computation of an {\it a~posteriori} error estimator. Here, {\it a~posteriori} error estimator is a combination of the edge estimator $\eta_h$ and
the volume  estimator $\mu_h$, that is,
\begin{equation}\label{esti_fin}
\norm{{\bf A^{-{1}/{2}}(p-p_h)}}^2+\norm{u-u_h}^2\leq C_{rel} (\eta_h^2+\mu_h^2),~~~~{\rm where}
\end{equation}
 %where
 \begin{align}
   \eta_h^2&:=\sum_{E\in \mathcal E_h}\eta_h^2(E)~ {\rm with} ~ \eta_h^2(E):=\norm{h_{E}^{1/2}[{\bf A^{-1}p_h+b^*}u_h]\cdot \tau_E}^2_E ;\label{eta}\\
%   \end{align}
%   \begin{align}  
  \mu_h^2&:=\sum_{T \in \mathcal T_h }\mu_h^2(T)~~{\rm with}~\nonumber\\
  &~~~\mu_h^2(T):=osc_h^2(f)_T
  +\norm{h ~{\rm div}~{\bf p_h}}^2_T+\norm{h ({\bf A^{-1}p_h }+u_h{\bf b^*})}^2_T. \label{mu}
\end{align}
 The {\it a~posteriori} estimate (\ref{esti_fin}) is  derived in Section 3, Theorem 3.2.\\
  %\eta_h^2=h_{E_j}\norm{[{\bf A^{-1}p_h+b^*}u_h]\cdot \tau_E}^2_E
%\end{eqnarray}\\
 %\item {\bf MARK}:- %Mark minimal set $\mathcal M_\ell$\\
The step {MARK} consists of the two alternatives {\bf (A)} and {\bf (B)} which depend on the computable quantities $\eta_h$ and $\mu_h$ and a positive parameter
$\kappa$. \\
%\bigskip
  {\bf Case (A)}: if $\mu^2_h\leq \kappa \eta_h^2 $, compute the minimal set of edges $\mathcal M_h \subset \mathcal E_h $ such that %{\mathcal E}
 \begin{equation}\label{A}
  \eta_h^2 \leq \eta_h^2(\mathcal M_h) ~~{\rm with}~~0<\theta_A<1.
   \end{equation}
 {\bf  Case (B)}: if $\mu^2_h> \kappa \eta_h^2 $,
compute the minimal set of triangles $\mathcal M_h \subset \mathcal T_h$ such that 
\begin{equation}\label{B}
\theta_B \mu_h^2 \leq \mu_h^2(\mathcal M_h)~~{\rm with}~~0<\theta_B<1. \end{equation}
Here $\theta_A$, $\theta_B$ and $\kappa$ are the parameters of the marking criteria and will be chosen  appropriately.\\
%\item {\bf REFINE}:
 Newest vertex bisection (NVB) algorithm  \cite{biven,steven} is  applied for refining the marked  edges  or elements and generate a new regular triangulation in  {REFINE} step.
Note that to maintain the conformity of the triangulation, some additional edges and elements may also need refinement.
% with the mark elements to form the conforming triangulation.      
%Refine \cite{biven} %(more detail required)
\begin{remark} Instead of separate marking (\ref{A})-(\ref{B}) defined in MARK,  one could use a collective marking, that is, compute the minimal set of triangles $\mathcal M_h \subset \mathcal T_h$
\begin{equation}
\theta_A (\eta_h^2 + \mu_h^2 ) \leq (\eta_h^2+\mu_h^2)(\mathcal M_h) ~~{\rm with}~~0<\theta<1.
\end{equation}
%can be use instead of the separate marking as defined in the step MARK.
\end{remark}
\section{Auxiliary results}
This section discusses some important results required  for the convergence analysis which are the {\it a~posteriori} error estimator, error and estimator reduction properties.\\
The  nonconforming finite element method (NCFEM) for (\ref{eq1}) seeks
 $ {u}_{CR}\in CR^1_0(\mathcal T_h)$ such that 
\begin{equation}\label{non1}
(\mathbf{A}\nabla_{NC}{u}_{CR}+ {u}_{CR}  {\bf b} ,\nabla_{NC} v_{CR})+(\gamma{u}_{CR},v_{CR})= (f, v_{CR}),{\forall} ~ v_{CR} \in CR^1_0(\mathcal T_h).
\end{equation}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
{\bf Representation of RTFEM Solution via NCFEM} \cite{CC,Marini}: 
%\vspace{0.35 cm}
%\medskip
The coefficients ${\bf A,b},\gamma$ are all piecewise constants.
Now the auxiliary discrete problem is to seek  $ u_h^N\in CR^1_0(\mathcal T_h)$ such that
% CR- FE space $(\mathcal T)$
\begin{equation}\label{eqna3}
 (\mathbf{A}\nabla_{NC} u_h^N+{u}_h {\bf b},\nabla_{NC} v_{CR})+(\gamma ~{u}_h,v_{CR})=
 (f_h,v_{CR}),\text{for all}\, v_{CR}\in CR^1_0(\mathcal T_h),
 \end{equation}
where for $T \in \mathcal{T}_h$ 
\begin{eqnarray} 
&& \displaystyle  u_h({\bf x})=\left(1+\frac{S(T)}{4} \gamma \right)^{-1} \left(\Pi_0 u_h^N+\frac{S(T)}{4}f_h \right)  \;\; \text{for}~ {\bf x} \in T,\label{umt}\\
&& ~~S(T)= \displaystyle{\int_T} ({\bf x}-\text{mid}(T))\cdot \mathbf{A}^{-1}({\bf x}-\text{mid}(T)) \,d{\bf x}.\label{s}
\end{eqnarray}
%where $ u_h^N\in CR^1_0(\mathcal T)$ 
%is the non-conforming finite element approximation which 
Then, the solution ${\bf p}_h$ of the mixed finite element method formulation (\ref{eqna1})-(\ref{eqna2}) satisfies
\begin{equation} \label{pmt}
 {\bf p}_h({\bf x})=-\left(\mathbf{A}\nabla _{NC} u_h^N+ u_h{\bf b} \right)+\left(f_h-\gamma  u_h \right)
 \frac{\left({\bf x}-\text{\rm mid}(T)\right)}{2} \text{ for}~ {\bf x}\in T.
\end{equation}
The well-posedness of (\ref{non1}) and (\ref{eqna3}) and the equivalence of (\ref{eqna1})-(\ref{eqna2}) with (\ref{eqna3}) is discussed in \cite{CC}.
% Let $u_h^N (= \tilde u_{CR})$ ncfem solution (of above equation ) for $\mathcal T_h$ and let define the quantity
%  ${\bf p_h^N}:=-(\mathbf{A}\nabla _{NC} u_h^N+ u_h{\bf b})$ .\\
% The following lemma enables a refined {\it a posteriori} error analysis for $\norm{u-u_h}$  and $\norm{{\bf p}-{\bf {\bf p_h}}}.$ 
%will be useful.
\begin{lemma}\label{lm8}
Let $u_h^N$ and $({\bf p_h}, u_h)$ solve (\ref{eqna3}) and (\ref{eqna1})-(\ref{eqna2}), respectively. Then, it holds
\begin{eqnarray}
&& \norm{\nabla_{NC} u_h^N}\leq \norm{\mathbf{A}^{-1}{\bf p_h}+u_h{\bf b}^*},\label{5.16}\\% ~~~{\rm and}\\
&&\norm{{\rm div}~{\bf p_h}} = \norm{f_h-\gamma u_h} \lesssim \norm{f_h}+\norm{\mathbf{A}^{-1}{\bf p_h}+u_h{\bf b}^*}.\label{5.17}
\end{eqnarray}
\end{lemma}
% %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
{\em Proof}. 
From (\ref{pmt}), 
 $${\bf A}^{-1} {\bf p_h}+u_h{\bf b}^*=-\nabla_{NC} u_h^N+\left(f_h-\gamma u_h \right){\bf A}^{-1}
\displaystyle{ \frac{\left({\bf x}-\text{\rm mid}(T)\right)}{2}}.$$
Since % $\mathbf{A_h}\nabla_{NC} \tilde u_{CR}\perp \left(f_h-\gamma_h \tilde u_h \right)
% {\left({\bf x}-\text{\rm mid}(T)\right)}/{2}$ 
$\left((f_h-\gamma  u_h)\left({\bf x}-\text{mid}(T)\right)/2, \nabla_{NC} u_h^N \right) =0,$
the Pythagoras theorem yields
\begin{equation*}
 \norm{{\bf A}^{-1} {\bf p_h}+u_h{\bf b}^*}^2=\norm{\nabla_{NC} u_h^N}^2
 +\norm{\left(f_h-\gamma u_h \right)\mathbf{A}^{-1} 
\displaystyle{ \frac{\left({\bf x}-\text{\rm mid}(T)\right)}{2}}}^2. \hspace{1cm}%\hfill{\Box}
\end{equation*}
Hence, (\ref{5.16}) holds. A use of triangle inequality with  (\ref{umt}) and (\ref{5.16})  implies (\ref{5.17}). \hfill${\Box}$%\endproof
\\
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
The  following theorem  is on   {\it a~posteriori}  error estimates of  ${\bf e_p:= p -p_h} $~and $e_u:=u-u_h$ %in $L^2$-norm. The following Theorem on {\it a posteriori} error estimator 
the proof of which is obtained by minor  modifications in the proof of Theorem 5.5 in \cite{CC}. However, for the sake of completeness, a short proof is given below. %is presented.
% is from \cite[Theorem 5.5]{CC}, it is rewritten with some modification for the completeness of the paper.% when the exact solution $u\in H^2(\Omega) \cap H^1_0(\Omega)$ is smooth. 
\begin{theorem}\label{5.5}({\it A~posteriori} error estimate)
 ~Let $u \in H^1_0(\Omega) $ be the unique weak solution of (\ref{eq1}) and %satisfying the regularity result  for some $0< \delta \leq 1.$ Further, 
  let $({\bf p_h},u_h) $ be  the solution of (\ref{eqna1})-(\ref{eqna2}).  For  small initial mesh-size  $h_1>0$
  there holds
\begin{eqnarray}\label{estimator1}
\norm{{\bf A^{-{1}/{2}}(p-p_h)}}^2+\norm{u-u_h}^2\leq C_{rel} (\eta_h^2+\mu_h^2), 
\end{eqnarray}
where  $0<h\leq h_1$ and $\eta_h$, $\mu_h$  are as defined in (\ref{eta}) and (\ref{mu}).
 \end{theorem}\\
 {\em Proof.} 
 Consider the Helmholtz decomposition: ${\bf e_p}={\mathbf A} \nabla z+\text{Curl} ~\beta$ for $z \in H^1_0(\Omega)$
  and $\beta \in H^1(\Omega)/\mathbb R$. Then% with ${\bf e_p= p -p_h}$
 \begin{equation}\label{ep11}
  (\mathbf{A}^{-1}{\bf e_p},{\bf e_p})=({\bf e_p},\nabla z)+(\mathbf{A}^{-1}{\bf e_p}, \text{Curl}~\beta).
 \end{equation}
 For the  first term  on the right-hand side of (\ref{ep11}), an integration by parts with  (\ref{eq3}) and  the fact ${\rm div~} {\bf p}_h+\gamma u_h=f_h$ lead to
 \begin{align}\label{ep12}
  ({\bf e_p},\nabla z)&=(\text{div}~{\bf e_p},z)=(f-f_h,z)-(\gamma (u-u_h),z)\nonumber\\
 & =(f-f_h,z-\Pi_h z)-(\gamma e_u,z),\nonumber \\
 &\lesssim {osc}_h (f)\;\norm{z}_1+\norm{ e_u} \norm{z}.
 \end{align}
 %Given any $v \in H^1_0(\Omega)$, 
Define $\beta_h:=I_h \beta$, where $I_h:H^1(\Omega)\rightarrow P^1(\mathcal T_h)\cap H^1_0(\Omega)$ is the Clement's interpolation operator \cite{verf}. With
$ {\rm Curl}~ \beta_h \in RT_0(\mathcal T_h)$, ${\rm Curl}~ \beta_h \bot ~\nabla H^1_0(\Omega )~(\bot~ {\rm denotes}~ L^2(\Omega)${ orthogonality}) and (\ref{eq3}),  the  second term 
 on the right-hand side of (\ref{ep11}) can be written as
 \begin{align*}
 (\mathbf{A}^{-1}{\bf e_p},~ \text{Curl}~\beta)
& =-(\mathbf{A}^{-1}{\bf p_h}+ u_h{\bf b}^*,~\text{Curl}~(\beta-\beta_h))-(e_u {\bf b^*},~\text{Curl}~\beta).%\nonumber \\
 \end{align*}
From the integration by part formula
%  for $\psi,$
% $$ \int_\Omega (\phi \text{curl}~\psi+\psi \cdot \text{curl}~ \phi )dx =\int_{\partial \Omega}\phi(\psi\cdot t)ds~~~~t=(-\nu_2,\nu_1)$$
\begin{equation}\label{ep13}
 (\mathbf{A}^{-1}{\bf e_p},~ \text{Curl}~\beta)
%& =-(\mathbf{A}_h^{-1}{\bf p_h}+ u_h{\bf b}_h^*,~\text{Curl}~(\beta-\beta_h))+(e_u {\bf b^*},~\text{Curl}~\beta)\nonumber\\
= \sum_{E\in \mc E_h} \int_{ E} [\mathbf{A}^{-1}{\bf p_h}+ u_h{\bf b}^*]\cdot \tau_E (\beta -\beta_h) ds - (e_u {\bf b^*},~\text{Curl}~\beta).
\end{equation}
With the interpolation estimates $\norm{\beta-\beta_h}_E\leq C h_E^{1/2}\norm{\beta}_{1,\omega_E}$, 
    the bounds 
  $\norm{\nabla \beta}_{\omega_E}=\norm{{\rm Curl}~ \beta}_{\omega_E}\leq \norm{\mathbf{A}^{-{1}/{2}} e_p}_{\omega_E}$,  and $\norm{z}\lesssim \norm{z}_1\lesssim \norm{{\bf A^{-1/2} e_p}}$, (\ref{ep11})-(\ref{ep13}) result in
%  %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 \begin{align}\label{p_esti}
  \norm{\mathbf{A}^{-{1}/{2}}{\bf e_p}}&\lesssim osc_h(f)
 +\norm{h_E^{1/2}[{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*]\cdot \tau_E}_{\mathcal E_h}+ \norm{e_u}.%\nonumber\\
 \end{align}
% %
 To estimate   $\norm{e_u},$ start with the triangle inequality 
 \begin{equation}\label{esti_5}
 \norm{e_u}\leq\norm{u-u_h^N}+\norm{u_h^N-u_h}.
 \end{equation}
 For $\tilde e= u_{CR}-u_h^N,$  Lemma 4.5 of \cite{CC} with sufficiently small mesh-size $h$ shows  
 \begin{eqnarray}\label{esti_9}
  \tnorm{\tilde e}_{NC}+\norm{\tilde e} &\lesssim &  osc_h(f).%+\Big(\norm{h_{\mathcal T}}\Big)\tnorm{h_{\mathcal{T}}u_h^N}_{NC}+ \norm{h^2_{\mathcal{T}}\;f_h}.
 \end{eqnarray}
% The estimates for $\norm{u-u_h^N}$ are derived with the help of result, that is, 
For any $\epsilon >0$, from lemma 3.3 of \cite{CC}, there exists  small mesh-size $h$  such that
\begin{equation}\label{dual}
\norm{u-u_{CR}} \leq \epsilon \tnorm{u-u_{CR}}_{NC} ~~{\rm holds}.
\end{equation}
 A repeated use of the triangle inequality yields estimates for $\norm{u-u_h^N}$ in (\ref{esti_5}) as
 %This proves
\begin{eqnarray}\label{esti_1}
   \norm{u-u_h^N}&\leq& \norm{u-u_{CR}}+\norm{u_{CR}-u_h^N}\nonumber\\%\lesssim \epsilon \tnorm{u-u_{CR}}_{NC}+\norm{u_{CR}-u_h^N}\nonumber\\
&\lesssim &\epsilon (\tnorm{u-u_h^N}_{NC}+ \tnorm{u_h^N-u_{CR}}_{NC})   +\norm{u_{CR}-u_h^N}\nonumber\\
   &\lesssim & \epsilon \norm{\nabla_{NC}(u-u_h^N)}+ {osc}_h(f).%\nonumber\\
  \end{eqnarray}
%Define ${\bf \tilde p}:=-(\mathbf A\nabla_{NC}u_h^N+u_h{\bf b} )$ % and the definition of ${\bf p}$ with (\ref{pmt}) gives
Let   ${\bf p}=-({\bf A}\nabla u+u {\bf b} )$ and define ${\bf \tilde p}=-(\mathbf A\nabla_{NC}u_h^N+u_h{\bf b} )$.  Along with an 
 addition and subtraction of the term  $\mathbf A^{-1} p_h,$ %. This shows% and the definition of ${\bf \tilde p}_{CR}$ yields
\begin{align}\label{esti_6}
 \norm{{\nabla_{NC}( u-u_h^N)}}&\lesssim \norm{\mathbf A^{-1/2}{\bf e_p }} +\norm{e_u}+ \norm{{\bf p_h}-{\bf \tilde p}}.
  \end{align}
% \begin{align}\label{esti_6}
% \norm{{\nabla_{NC}( u-u_h^N)}}&\lesssim \norm{\mathbf A^{-1}{\bf e_p }} +\norm{e_u{\bf b^*}}+ \norm{\mathbf A^{-1}({\bf p_h}-{\bf \tilde p})}.
%  \end{align} 
 For the third term on the right-hand side of (\ref{esti_6}), (\ref{pmt}) leads to %and $\norm{u_h}\lesssim \norm{f}$  leads to
 \begin{equation}\label{esti_7}
\norm{{\bf p_h}-{\bf \tilde p}}\leq \norm{(f_h-\gamma u_h)({\bf x}-\text{mid}(T))}\lesssim \norm{h(f_h-\gamma\; u_h)}.
\end{equation}
 The combination of (\ref{esti_1})-(\ref{esti_7}) results in %along with stability of solutions $u $ and $u_h^N$ yields %$\norm{u}\lesssim \norm{f}, ~\tnorm{u_h^N}_{NC}\leq \norm{f}$  yields
 \begin{align}\label{esti_119}
 \norm{u-u_h^N}\lesssim &~ {osc}_h(f)+\epsilon \big(\norm{\mathbf{A^{-1/2}}{\bf e_p}}+ \norm{e_u}\big)+ \epsilon\norm{h(f_h-\gamma\;u_h)}.%+ \tnorm{h u_h^N}_{NC}.
 %+ \epsilon (\norm{(\mathbf A^{-1}-\mathbf A^{-1}_h){\bf p_h }}+\norm{u_h({\bf b}^*-{\bf b}_h^*)}).
 \end{align}
% %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
To bound $\norm{u_h^N-u_h}$ in (\ref{esti_5}), use (\ref{umt}), the triangle inequality and the fact that $S(T)\approx h^2$ to obtain
 \begin{eqnarray}\label{esti_3}
  \norm{u_h^N-u_h}&\leq &\left(1+\frac{S(\mathcal{T})}{4} \gamma \right)^{-1} \norm{u_h^N-\Pi_0 u_h^N+\frac{S(\mathcal{T})}{4}(\gamma u_h^N- f_h )}\nonumber\\
 &\lesssim & \norm{h\nabla_{NC}u_h^N} + \tnorm{h^2(u_h^N- u_h)}_{NC}+\norm{h^2\;(f_h-\gamma u_h)}.
 \end{eqnarray}
% %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
 %The combination of the previous estimates with 
 A use of (\ref{esti_9})-(\ref{esti_3}) in  (\ref{esti_5}) along with (\ref{eqna2}) leads to %and (\ref{esti_4})
 \begin{equation}\label{esti_8}
 \norm{e_u}\lesssim {osc}_h (f)+ \epsilon \Big(\norm{\mathbf{A^{-1/2}}{\bf e_p}}+ \norm{e_u}\Big)+ \epsilon\norm{h~{\rm div}~{\bf p_h}}
%  &~~~+\Big( 1 +\norm{h^{-1}_{\mathcal{T}}(\mathbf{A- A}_h)}_\infty+\norm{h^{-1}_{\mathcal{T}}(\mathbf{b- b}_h)}_\infty
%     +\norm{h^{-1}_{\mathcal{T}}(\gamma- \gamma_h)}_\infty\Big)\;\nonumber \\
 +\norm{h(\mathbf{A}^{-1}{\bf p_h}+u_h{\bf b}^*)}.
 %&+\Big(h+h^2+\norm{{\bf A-A}_h}_\infty+ \norm{{\bf b-b}_h}_\infty +\norm{\gamma- %\gamma_h}_\infty\Big)\norm{\mathbf{A}_h^{-1}{\bf p_h}+u_h{\bf b}_h^*}\nonumber \\ 
 %+ \epsilon \Big(\norm{(\mathbf A^{-1}-\mathbf A^{-1}_h){\bf p_h }}+\norm{u_h({\bf b}^*-{\bf b}_h^*)}\Big).
 \end{equation}
% %
 For  small mesh-size $h_1>0$ with $0 < h\leq h_1,$  (\ref{esti_8}) and (\ref{p_esti}) prove (\ref{estimator1}). 
 %The proof of  (\ref{estimator2}) utilizes the  additional regularity with $\epsilon =o(h^{\delta})$.% yields the required result and this concludes the proof  
% \qquad \endproof%
 \hfill{$\Box$}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\\
\begin{lemma}\label{lmeff}(Efficiency)~Let $({\bf p},u)$ be the solution of (\ref{eq4}) and $({\bf p_h},u_h)$ be the solution of (\ref{eqna1})-(\ref{eqna2})
over the triangulation $\T_h$. Then, it holds
\begin{align}\label{edeff}
C_{\rm eff}(\norm{h_E^{1/2}[{\bf A^{-1}p_h+}u_h{\bf b^*}]\cdot \tau_E}_{\mc E_h}+ &\norm{h({\bf A^{-1}p_h+b^*}u_h)})\nonumber\\
 &\leq  \norm{u-u_h}+\norm{{\bf A^{-1/2}(p-p_h)}}.%+\norm{h({\bf A^{-1}p_h+b^*}u_h)}
 \end{align}
\end{lemma}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
{\it Proof.} The proof is divided into two steps.\\
{\em Step 1}. 
Let $b_E$ denote the continuous edge bubble function satisfying $0\leq b_E\leq 1$ on $\omega_E$ and $b_E\in P_2(T)$ for each $T \subset \omega_E$. Let $\phi:= [{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*]\cdot \tau_E $ be a polynomial function along $E$. There exists an extension operator $P:\mathcal C(E)\rightarrow \mathcal C(\omega_E)$ \cite{verf}, 
where $\mathcal C(E)~({\rm resp}. ~\mathcal C (\omega_E))$ denotes the space of the continuous 
functions defined on $E$ (resp. $\omega_E$) such that the  operator $P$ satisfying $P\phi|_E=\phi$ and
\begin{equation}\label{eff_2}
h_E^{1/2}\norm{\phi}_E \lesssim \norm{b_E^{1/2}P\phi}_{\omega_E}\lesssim h_E^{1/2}\norm{\phi}_E.
\end{equation}
An equivalence of norm argument  implies
\begin{equation*}
\norm{\phi}^2_E\lesssim \norm{b_E^{1/2}\phi}_E^2= \int_E (b_EP\phi)[{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*]\cdot \tau_E~ ds.
\end{equation*}
An integration by parts and a use of  ${\rm Curl}({\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*)=0 $ result in
\begin{equation}\label{eff_42}
\norm{\phi}^2_E\lesssim \norm{b_E^{1/2}\phi}_E^2\leq -\int_{\omega_E} {\rm Curl}(b_EP\phi)\cdot ({\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*)~ ds.
\end{equation}
Note that $(\nabla u, {\rm Curl}(b_E^{1/2}P\phi))_{\omega_E}=0$. Hence, with the help of the 
Cauchy-Schwarz inequality,  (\ref{eff_42}) can be written as
\begin{equation*}
\norm{\phi}^2_E\lesssim \norm{{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*+\nabla u}_{\omega_E} \norm{{\rm Curl}(b_EP\phi)}_{\omega_E}.
\end{equation*}
The inverse inequality, (\ref{eff_2}) % and $|b_E P\phi|_{\omega_E}\lesssim h^{1/2}\norm{\phi}_E$\cite{verf}, and
and a utilization of the definition ${\bf p}=-(\mathbf A \nabla u+u{\bf b} )$, result in
\begin{align*}
%& \min_{v\in H^1_0(\Omega)}\norm{\textbf A^{-1}{ \bf p_h} +u_h{\textbf b^*} -\nabla v}\leq \norm{\textbf A^{-1}{ \bf p_h} +u_h{\textbf b^*} +\nabla u}\\
%&~~=\norm{\textbf A^{-1}({ \bf p-p_h})} +\norm{{(u-u_h)\textbf b^*}}\\
\norm{h_E^{1/2}[{\bf A^{-1}p_h}+u_h{\bf b^*}]\cdot \tau_E}_{E}~~\lesssim \norm{ {\bf p-p_h}}_{\omega_E} +\norm{u-u_h}_{\omega_E}.
\end{align*}
A summation over all the edges leads to an estimate of the first term on the left-hand side of (\ref{edeff}).\\% Second part of (\ref{edeff}) is follow from \cite[Theorem 5.7]{CC}. \\
{\em Step 2}. Define the function ${\textbf q_T}:= b_T(\mathbf A^{-1}{ \bf p_h} +u_h{\textbf b^*} )\in P_4(T)\cap W^{1,\infty}_0(T)$ and the cubic bubble function 
$
 ~b_T = 27\lambda_1 \lambda_2 \lambda_3 \in P_3(T)\cap C_0(T)
$
in terms of the barycentric coordinates $\lambda_1, \lambda_2, \lambda_3$ of $T \in \mathcal T_h$ \cite{verf}.
Since ${\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*$ is affine on $T\in \mathcal T_h$, an equivalence of norm argument shows
\begin{align*}%\label{effi2}
 \norm{\mathbf A^{-1}{ \bf p_h} +u_h{\textbf b^*} }_T^2&\lesssim  \int_T{\textbf q_T}\cdot({ \mathbf A^{-1}\bf p_h} +u_h{\textbf b^*} )~dx.
 \end{align*}
The definition of ${\bf p}$ and (\ref{eq3}) show that %the right-hand side of (\ref{effi2}) reads
\begin{align*}
\norm{{\mathbf A^{-1} \bf p_h} +u_h{\textbf b^*} }_T^2 
& \lesssim \int_T{\textbf q_T}\cdot \big ( \mathbf A^{-1}({ \bf p_h - p}) -(u-u_h){\textbf b^*} \big)~dx
%& ~~+ \int_T{\textbf q_T}\cdot\Big ( (\mathbf A^{-1}_h-\mathbf A^{-1}){\bf p_h}-u_h({ \bf b^* -b}_h^*)  \Big)~dx \\
~-\int_T{\textbf q_T}\cdot \nabla u ~dx.                      
 \end{align*}
The Cauchy-Schwarz inequality with $\norm{{\textbf q_T}}_T\lesssim \norm{\mathbf A^{-1}_h{\bf p_h}+u_h{\textbf b_h^*} }_T$ is employed in the first 
two terms. Adding the zero terms $\nabla u_h|_T$ to the right-hand side of the above equation,  an integration by parts  %with $\mathbf A $ as defined in Section 3.5 in the last term to yield
%in the last term 
 shows  that 
\begin{eqnarray*}%\label{effi3}
 h_T\norm{{\mathbf A^{-1} }{\bf p_h}+u_h
 {\bf b} }^2_T &\lesssim &h_T\norm{\mathbf A^{-1}{\bf p_h}+u_h{\bf b}^* }_T\Big ( \norm{{\bf p-p_h}}_T
+\norm{u-u_h}_T \Big )\nonumber\\
 &&+ h_T \int_T (u-u_h)\text{div}~{\bf q}_T ~dx.
 \end{eqnarray*}
Since $\textbf q_T \in P_4(T)$, an inverse estimate yields%and \cite{ }
\begin{equation*}
 h_T\norm{\text{div}~ \textbf q_T}_T\lesssim \norm{\textbf q_T}_T\lesssim \norm{\mathbf A^{-1}{\bf p_h}+u_h{\bf b}^* }_T.
\end{equation*}
 %controls the last contribution in (\ref{effi3}). 
Since $ h_T \lesssim 1$, it follows
\begin{align*}
h_T\norm{\mathbf A^{-1}{\bf p_h}+u_h{\bf b}^* }_T&\lesssim \norm{u-u_h}_T+\norm{{\bf p}-{\bf p_h}}_T.
\end{align*}
A summation  over all elements leads to an estimate of  second  term on the left-hand side of (\ref{edeff}). This
concludes the proof. \hfill{$\Box$}%\qquad \endproof%
% The sum over all triangles implies 
% \begin{align*}%\label{effi1}
%h\norm{\mathbf A^{-1}{\bf p_h}+u_h{\bf b}^* }&\lesssim \norm{u-u_h}+\norm{{\bf p}-{\bf p_h}}.\nonumber \qquad\qquad\qquad
% \end{align*}\hfill{$\Box$}
\\
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\\
 %In the following results,
 Let $\mc T_H$ and $\mc T_h$ with $H<h$ denote  nested triangulations, %be a refined triangulation of $\T_H$,
  and $({\bf p_h}, u_h)$ and $({\bf p_H}, u_H)$ denote the solutions of  (\ref{eqna1})-(\ref{eqna2}) obtained with right-hand sides $f_h$ and $f_H$ over
  $\mathcal T_h$ and $\mathcal T_H$, respectively.
 The following notations are used  in the sequel: %the purpose of simplicity.% in the paper
\begin{subequations}\label{com_1}
\begin{eqnarray}
&&E_H^2:=\norm{{\bf A}^{-1/2}{\bf (p_h-p_H)}}^2+\norm{u_h-u_H}^2,\\
&&e_{p_h}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_h)}}^2,\qquad\qquad e_{p_H}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_H)}}^2,\\
&&e_{u_h}^2:=\norm{u-u_h}^2, \qquad~~~~~~~~~~~\qquad e_{u_H}^2:=\norm{u-u_H}^2,\\
&&e_h^2~:=e_{p_h}^2+e_{u_h}^2,\qquad\qquad ~~~~~~~~~~~~~ e_H^2~:=e_{p_H}^2+e_{u_H}^2.
\end{eqnarray}
\end{subequations}
%\begin{subequations}
%Maxwell's equations:
%\begin{align}
%        B'&=-\nabla \times E,\\
%        E'&=\nabla \times B - 4\pi j,
%\end{align}
%\begin{equation}\label{com_2}
%  \begin{aligned}
%  &&E_H^2:=\norm{{\bf A}^{-1/2}{\bf (p_h-p_H)}}^2+\norm{u_h-u_H}^2,\\
%&&e_{p_h}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_h)}}^2,\qquad\qquad e_{p_H}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_H)}}^2,\\
%&&e_{u_h}^2:=\norm{u-u_h}^2, \qquad~~~~~~~~~~~\qquad e_{u_H}^2:=\norm{u-u_H}^2,\\
%&&~e_h^2:=e_{p_h}^2+e_{u_h}^2,\qquad\qquad ~~~~~~~~~~~~~~ e_H^2:=e_{p_H}^2+e_{u_H}^2.
%% &&E_H^2:=\norm{{\bf A}^{-1/2}{\bf (p_h-p_H)}}^2+\norm{u_h-u_H}^2, &&\\
%%&&e_{p_h}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_h)}}^2,\qquad\qquad\qquad &e_{p_H}^2&:=\norm{{\bf A}^{-1/2}{\bf  (p-p_H)}}^2,\\
%%&&e_{u_h}^2:=\norm{u-u_h}^2, \qquad  \qquad &&e_{u_H}^2:=\norm{u-u_H}^2,\\
%%&&e_h^2:=e_{p_h}^2+e_{u_h}^2, \qquad \qquad &&e_H^2&:=e_{p_H}^2+e_{u_H}^2.\\
%   \end{aligned}
%\end{equation}
%\begin{equation}
%  \begin{aligned}
%  w^T x_i + b &\geq 1-\xi_i &\text{ if }& y_i=1,  \\
%  w^T x_i + b &\leq -1+\xi_i & \text{ if } &y_i=-1,
%  \end{aligned}
%\end{equation}
%\end{subequations}
% Let $E_H^2=\norm{{\bf A}^{-1/2}{\bf (p_h-p_H)}}+\norm{u_h-u_H}$ and $e_h^2=e_{p_h}^2+e_{u_h}^2$, where
% $e_{p_h}^2:=\norm{{\bf A}^{-1/2}{\bf  (p-p_h)}}^2$ and $e_{u_h}^2=\norm{u-u_h}^2$. Similarly, $e_H^2$ is define over triangulation $\mathcal T_H$.

 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\ref{com_1},\ref{com_2}
\begin{lemma}\label{lm3}(Volume estimator reduction) For given $0<\theta_B \leq 1$, there exist constants $0<\delta_1,~ \rho_B < 1$ and  the positive constant $\varLambda_2$ such that $\mu_H$ and $E_H$ defined in (\ref{mu}) and (\ref{com_1}) satisfy
% Let $\mathcal T_h$ and $\mathcal T_H$ be two nested conforming triangulation, %and $\mathcal M \subset \mathcal T_H$, %, $\eta$ is edge estimator term (\ref{})
% then with  {\bf Case (B)}, for given $\theta_B>0$ there exists $ 0<\rho_B<1$ and  constant $\varLambda_1 > 0$  such that
%%If ...there exists $ 0<\rho_A, \rho_B<1$ 
\begin{eqnarray}
&&  \mu_h^2 \leq (1+\delta_1) \mu_H^2 +\varLambda_1 H^2 E_H^2 \qquad \qquad {\rm for~~the~~ {\bf Case ~(A)}},\label{et_rd_11} \\
&&  \mu_h^2 \leq \rho_B \mu_H^2 +\varLambda_1 H^2 E_H^2 \qquad \qquad \qquad {\rm for~~the ~~{\bf Case ~(B)}}.\label{et_rd11} 
  %\mu_h^2 \leq  2\mu_H^2 +\varLambda_1 H^2 E_H^2\label{et_rd1}\\
 % \eta_h^2 \leq \rho_A \eta_H^2 +\varLambda_2 E_H^2 \label{et_rd2}
  %\eta_h^2 \leq (1+\delta_2) \eta_H^2 +\varLambda_2 E_H^2 \label{et_rd3}
\end{eqnarray}
%where ).
\end{lemma}
{\it Proof}.  %First volume estimator reduction for the $\mu_h$\\
For any  triangle $K \in \mathcal T_H$, which gets refined in the level $\mc T_h$, that is, $K\in \mc T_H\setminus \mc T_h$, there exist triangles  $T_1,T_2,\cdots,T_J $
 such that $K= T_1\cup T_2 \cdots \cup T_J$. For $K \in \mathcal T_H \setminus \mathcal T_h$ %and the constant 
 \begin{eqnarray*}
\mu^2_h(K)&=&\sum_{j=1}^J\Big( \norm{h(f-f_h)}_{T_j}^2+\norm{h~{\rm div}~{\bf p_h}}^2_{T_j}+\norm{h ({\bf A^{-1}p_h }+u_h{\bf b^*})}^2_{T_j}\Big).
\end{eqnarray*}
At least one refinement of $T$ implies $h\leq H/2$,  the fact $\norm{f-f_h}_{T_j}\leq\norm{f-f_H}_{T_j}$  along with the triangle 
inequality yields
\begin{eqnarray}
\mu^2_h(K)&\leq& \frac{H^2}{4}\sum_{j=1}^J\Big( \norm{f-f_H}_{T_j}^2+\norm{~{\rm div}~{(\bf p_h- p_H)}+{\rm div}~{\bf p_H}}^2_{T_j}\nonumber \\
&&+\norm{({\bf A^{-1}p_H }+u_H{\bf b^*})+{\bf A^{-1}(p_h-p_H) }+(u_h-u_H){\bf b^*}}^2_{T_j}\Big)\nonumber
\end{eqnarray}
%The factor $\frac{1}{4}$ is result from $J\geq 2$ refinement, that is, atleast one bisection of $T\in \T_H$. Hence  
%A use of Young's inequality 
For each $K\in \T_H$, a use of $\norm{f_h-f_H}_K\leq \norm{\Pi_h(1-\Pi_H)f}_K\leq \norm {f-f_H}_K$ and  (\ref{eqna2}), with  Young's inequality yields for $\delta_1>0$
\begin{align}\label{51}
\mu^2_h(K)&\leq \frac{H^2}{4}\sum_{j=1}^J\big( \norm{(f-f_H)}_{T_j}^2+2\norm{f_h-f_H}^2_{T_j}+ 3\norm{{\rm div}~{\bf p_H}}^2_{T_j}+6\norm{\gamma (u_h-u_H)}^2\nonumber \\
&~+(1+\delta_1)\norm{{\bf A^{-1}p_H }+u_H{\bf b^*}}^2_{T_j}+(1+\frac{1}{\delta_1})\norm{{\bf A^{-1}(p_h-p_H) }+(u_h-u_H){\bf b^*}}^2_{T_j}\big)\nonumber\\
&\leq \frac{3}{4}(1+\delta_1) \mu_H^2(K) +(1+\frac{1}{\delta_1}) C_1 H^2 E_H^2(K),%{H_{(\sum_j{T_j})}}^2,
\end{align} 
where $C_1=\max\{\norm{{\bf A}^{-1/2}}_\infty^2,\norm{{\bf b^*}}_\infty^2+6 \norm{\gamma}_\infty^2\}$. Denote  $(1+\frac{1}{\delta_1}) C_1=:\varLambda_1.$\\
For $K\in \mathcal T_H \cap \mathcal T_h$,
 \begin{eqnarray}\label{52}
\mu^2_h(K)%&=&( \norm{h(f-f_h)}_{T}^2+\norm{h(f_h-\gamma u_h)}^2_{T}+\norm{h ({\bf A^{-1}p_h+b^* }u_h)}^2_{T}\nonumber\\
% &=&( \norm{H(f-f_H)}_{T}^2+\norm{H(f_H-\gamma u_H+\gamma u_H-\gamma u_h)}^2_{T}\nonumber\\
% &&+\norm{H ({\bf A^{-1}p_H+b^* }u_H+{\bf A^{-1}(p_h-p_H)+b^* }(u_h-u_H))}^2_{T}\nonumber\\
&\leq& (1+\delta_1)\mu_H^2(K) +\varLambda_1 H^2 E_{H}^2{(K)}.
 \end{eqnarray}
From (\ref{51}) and (\ref{52}), a summation over all triangles implies %$ \varLambda_1 =(1+\frac{1}{\delta_1})C_1$
%for   {\bf Case (A)} it holds %Then,  for  For any refinement (\ref{vol4}) gives%the choice of $\delta_1 =1$ in
 \begin{eqnarray}\label{et_rd1}
\mu_h^2 \leq (1+\delta_1) \mu_H^2+\varLambda_1  H^2 E_H^2,
\end{eqnarray}
for  the  {\bf Case (A)}.\\
For  the {\bf Case (B)}, a summation over all triangles, 
the marking criteria for $\mc M_H\subset \T_H$, 
\begin{equation}\label{81}
\theta_B \mu_H^2 \leq \mu_H^2(\mathcal M_H) \leq \mu_H^2 (\mathcal T_H \setminus T_h).
\end{equation}
  A use of (\ref{51}), (\ref{52}) and (\ref{81})  leads to the sharper bound 
\begin{eqnarray*}\label{vol4}
\mu_h^2 & =&\mu_h^2(\mathcal T_H \setminus T_h) +\mu_h^2(\mathcal T_h \cap \mathcal T_H)\nonumber \\
&\leq& \frac{3}{4}(1+\delta_1) \mu_H^2(\mathcal T_H \setminus T_h) +(1+\delta_1)\mu_H^2(\mathcal T_h \cap \mathcal T_H)+\varLambda_1 H^2 E_H^2.\nonumber\\
 %\end{eqnarray}
% As (\ref{vol4}) holds for any refinement of $\mathcal T_H$, that is independent of any specific marking criteria.  Hence for   {\bf Case (A)} it holds %Then,  for  For any refinement (\ref{vol4}) gives%the choice of $\delta_1 =1$ in
% \begin{eqnarray}\label{et_rd1}
%\mu_h^2 \leq (1+\delta_1) \mu_H^2+\varLambda_1  H^2 E_H^2.
%\end{eqnarray}
% {\bf Case (B)} marking criteria for $\mc M_H\subset \T_H$, $\theta_B \mu_H^2 \leq \mu_H^2(\mathcal M_H) \leq \mu_H^2 (\mathcal T_H \setminus T_h)$, leads to
%\begin{eqnarray*}
&\leq &(1+\delta_1)\mu_H^2- \frac{1}{4}(1+\delta_1) \mu_H^2(\mathcal T_H \setminus T_h)+\varLambda_1  H^2 E_H^2\nonumber\\
 & \leq &(1+\delta_1)(1-\frac{\theta_B}{4})\mu_H^2+\varLambda_1  H^2 E_H^2.
\end{eqnarray*}
%The choice of $\delta_1 =\frac{\theta_B}{8-2\theta_B}$ 
For given  $0<\theta_B<1$, the selection  $0<\delta_1< {\theta_B}/{(4-\theta_B)}$  results in  \\
$0<\rho_B=(1+\delta_1)(1-\frac{\theta_B}{4})<1$.~%AND hence (\ref{et_rd11}) holds.
 This concludes the proof.\hfill${\Box}$ %\qquad \endproof%
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\noindent 
\begin{remark}\label{rmk3} Irrespective of the marking criteria, from (\ref{et_rd_11})-(\ref{et_rd11}), it
 follows that %s (\ref{vol4}) holds for any refinement of $\mathcal T_H$, that is independent of any specific marking criteria. Then,  for any refinement (\ref{vol4}) gives%the choice of $\delta_1 =1$ in
 %\begin{equation}%\label{et_rd1}
$$\mu_h^2 \leq 2 \mu_H^2+\varLambda_1  H^2 E_H^2.$$
 %\end{equation}
\end{remark}
\begin{corollary}\label{lm13} Let  $\T_h$ be a refined triangulation of $\T_H$. Then it holds,
 \begin{eqnarray}\label{lm12}
\frac{\gamma_0^2}{2}\mu^2_H\leq \mu_h^2  +2 \norm{H(f_h-f_H)}^2+C_1 H^2 E_H^2,
 \end{eqnarray}
 where for $T\in \T_h$,  $T' \in \T_H$ and $T\subset \T'$, $\gamma_0$ is defined as
 $0<\gamma_0<1$ if $T'$ get refined, otherwise $\gamma_0=1$.
\end{corollary}
% \[ \gamma_0 = \left\{
%  \begin{array}{l l}
%    (0,1) & \quad \text{if $T \in \T_H$ is refined,}\\
%   1 & \quad \text{otherwise.}
%  \end{array} \right.\]
%  let $T\in \T_h$,  $T' \in \T_H$ and $T\subset \T'$, set $\gamma_0H=h, $ here $0<\gamma_0<1$ if $T'$ get refined, otherwise $\gamma_0=1$.
 %{\it Proof}. Let $T\in \T_h$,  $T' \in \T_H$ and $T\subset T'$. Set $\gamma_0H=h, $ here $0<\gamma_0<1$ if $T'$ get refined, otherwise $\gamma_0=1$.
% %Set $h_T=\gamma_0 h_{T'}$, here $h_T,~ h_{T'}$ are the diameter 
%%of $T $ and $T'$ respectively, and $0<\gamma_0<1$ if $T'$ get refined, otherwise $\gamma_0=1$. ~%First define $\gamma_0=\Big \{ \gamma ~~0<\gamma <1 and ~ 1 ~if T not refined$
%Recall the definition of $\mu_H$
% \begin{eqnarray*} 
%\frac{\gamma_0^2}{2}\mu^2_H(T')= \sum_{T\subset T'}\Big[ \frac{\gamma_0^2}{2} \norm{H(f-f_h+f_h-f_H)}^2+ \frac{\gamma_0^2}{2}\norm{H(f_H-f_h+f_h)}^2+\\
%\frac{\gamma_0^2}{2}\norm{H ({\bf A^{-1}(p_h-p_h+p_H) }+(u_h-u_h+u_H){\bf b^*})}^2\Big].\nonumber
%\end{eqnarray*}
%The Young's inequality, implies
%%as $\gamma_0H=h$
%\begin{eqnarray*}
%\frac{\gamma_0^2}{2}\mu^2_H(T')\leq \mu^2_h(T')+ 2\norm{H(f_h-f_H)}^2_{T'}+C_1 H^2 (\norm{{\bf A}^{-1/2}{\bf( p_h-p_H)}}^2_{T'}+\norm{u_h-u_H}^2_{T'}).%\frac{\gamma_0^2}{2}\big( 2\norm{H(f-f_h)} \big)
%\end{eqnarray*}
%A summation over all triangle results in (\ref{lm12}).
%
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

\begin{lemma}\label{lm2}( Edge estimator reduction) Given $0<\theta_A \leq 1$, there exist constants $0<\delta_2,~ \rho_A < 1$  and  a positive constant $\varLambda_2$ such  $\eta_H$ and $E_H$ defined in (\ref{eta}) and (\ref{com_1})
satisfy
% Let $\mathcal T_h$ and $\mathcal T_H$ be two nested conforming triangulation. %and $\mathcal M \subset \mathcal T_H$, %, $\eta$ is edge estimator term (\ref{})
% Then, with {\bf Case (A)} for given $0<\theta_A<1$, there exists $ 0<\rho_A<1$ and  constant $\varLambda_2 > 0$  such that
\begin{eqnarray}
 % \mu_h^2 \leq \rho_B \mu_H^2 +\varLambda_1 H^2 E_H^2\label{et_rd11} \\
  %\mu_h^2 \leq  2\mu_H^2 +\varLambda_1 H^2 E_H^2\label{et_rd1}\\
  \eta_h^2 \leq \rho_A \eta_H^2 +\varLambda_2 E_H^2 \qquad \qquad {\rm for~~the~~{\bf Case~ (A)}},\label{et_rd2} \\  
    \eta_h^2 \leq (1+\delta_2)\eta_H^2 +\varLambda_2 E_H^2 \qquad \qquad ~ {\rm for~~the ~~{\bf Case ~(B)}}.\label{et_rd_2}
     %\eta_h^2 \leq (1+\delta_2) \eta_H^2 +\varLambda_2 E_H^2 \label{et_rd3}
\end{eqnarray}
\end{lemma}
{\it Proof}. For all $F \in \mathcal E_H$, either $F\in \mathcal E_h$ or there exist $E_1, E_2, ... E_J \in \mathcal E_{h}$  with 
$F = E_1\cup E_2\cup ...\cup E_J$ for $J \geq 2$. For the case  $F \in \mathcal E_H \setminus \mathcal E_h$, a use of Young's inequality yields 
%for  $0< \delta_2 ={\theta_A}/{(4-2\theta_A)}$,
\begin{eqnarray}\label{99}
\eta_h^2(F)&=&\sum_{j=1}^J \eta^2_h (E_{j})=\sum_{j=1}^J \norm{h_{E_j}^{1/2}[{\bf A^{-1}p_h}+u_h{\bf b^*}]\cdot \tau_{E_j}}^2_{E_j} \nonumber\\%^{\frac{1}{2}}
&\leq & (1+\delta_2)\frac{H_{F}}{2}\norm{[{\bf A^{-1}p_H}+u_H{\bf b^*}]\cdot \tau_{E}}^2_F \nonumber\\
&&+(1+\frac{1}{\delta_2})\sum_{j=1}^J {h_{E_j}}\norm{[{\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}]\cdot \tau_{E_j}}^2_{E_j}.
\end{eqnarray}
Here $H_F$ denotes the lenght of edge $F$. The fact that  $J\geq 2$ implies that there exists  at least one bisection of $F \in \mathcal E_H \setminus \mathcal E_h$.\\
For $F\in \mc E_H \cap \mc E_h$,
\begin{eqnarray}\label{98}
\eta_h^2(F)%&=&\sum_{j=1}^J \eta^2_h (E_{j})=\sum_{j=1}^J \norm{h_{E_j}^{1/2}[{\bf A^{-1}p_h}+u_h{\bf b^*}]\cdot \tau_{E_j}}^2_{E_j} \nonumber\\%^{\frac{1}{2}}
&\leq & (1+\delta_2)\norm{{H_{F}}[{\bf A^{-1}p_H}+u_H{\bf b^*}]\cdot \tau_{E}}^2_E \nonumber\\
&&+(1+\frac{1}{\delta_2}) {H_{F}}\norm{[{\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}]\cdot \tau_{F}}^2_{F}.
\end{eqnarray}
For any $E\in \mathcal E_h$ with $E \nsubseteq \cup \mathcal E_H$,  $E$ is the interior edge of some element of $T\in \T_H$ and hence,  $ [{\bf A^{-1}p_H}+u_H{\bf b^*}]\cdot \tau_E =0$ along $E$ implies  $\eta_h^2(E)=\norm{h_E[{\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}]\cdot \tau_E}^2_E$. 
Now, consider this and the cases (\ref{99})-(\ref{98})   to obtain
\begin{eqnarray*}
\eta_h^2 =\sum_{E \in \mathcal E_h} \eta_h^2(E) &\leq &\frac{1+\delta_2}{2} \sum_{E\in \mathcal E_H \setminus \mathcal E_h} \eta_H^2(E)
+ (1+\delta_2) \sum_{E\in \mathcal E_H \cap \mathcal E_h} \eta_H^2(E)\nonumber \\ &&~+ (1+\frac{1}{\delta_2})\sum_{E \in \mathcal E_h} 
\norm{h_E [{\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}]\cdot \tau_E}^2_{E}.
\end{eqnarray*}
The inverse inequality $\norm{{\bf q_H}}^2_{{E}} \lesssim H^{-\frac{1}{2}}\norm{{\bf q_H}}^2_{{\omega_E}}$ results in
\begin{eqnarray*}
\norm{h_E [{\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}]\cdot \tau_E}^2_{E}%&\leq& \norm{h_E \big({\bf A^{-1}(p_h-p_H)+b^*}(u_h-u_H)\big)|_{T+}}^2_{L^2(E)}\nonumber\\
%&&+\norm{h_E \big({\bf A^{-1}(p_h-p_H)+b^*}(u_h-u_H)\big)|_{T-}}^2_{L^2(E)}\nonumber\\
&\lesssim &  \norm{ \big({\bf A^{-1}(p_h-p_H)}+(u_h-u_H){\bf b^*}\big)}^2_{\omega_E}
\end{eqnarray*}
for the edge patch $\omega_E$ of $E$ in $\mathcal T_h$. Since there is only a finite overlap of all edge patches, there holds
\begin{eqnarray*}%\label{edge_4}
 \eta_h^2\leq \frac{1+\delta_2}{2}\sum_{E \in \mathcal E_H \setminus \mathcal E_h } \eta_H^2+(1+\delta_2)\sum_{E \in \mathcal E_H \cap \mathcal E_h } \eta^2_H
+C_2 (1+\frac{1}{\delta_2}) E_H^2.
\end{eqnarray*}
%Where $C_2=\max \{ \}$.\\
%{\bf Case (A)} marking  criteria  for $\mc M_h\in \mc E_H$, $\theta_A \eta_H^2 \leq \eta_H (\mathcal M_H)\leq \eta (\mathcal E_H \setminus \mathcal E_h)$ and 
%  denoting the constant  $C_2 (1+\frac{1}{\delta_2})$ by $\varLambda_2$, lead to
Denote $C_2 (1+\frac{1}{\delta_2})=:\varLambda_2$ to obtain 
\begin{eqnarray}\label{edge_4}
 \eta_h^2\leq \frac{1+\delta_2}{2}\sum_{E \in \mathcal E_H \setminus \mathcal E_h } \eta_H^2+(1+\delta_2)\sum_{E \in \mathcal E_H \cap \mathcal E_h } \eta^2_H
+\varLambda_2 E_H^2.
\end{eqnarray}
 For the {\bf Case (A)}, the marking criteria  for $\mc M_H\in \mc E_H$,   $$\theta_A \eta_H^2 \leq \eta_H (\mathcal M_H)\leq \eta (\mathcal E_H \setminus \mathcal E_h)$$
 leads to 
\begin{eqnarray*}
\eta_h^2 &\leq & (1+\delta_2) \eta_H^2 -\frac{(1+\delta_2)}{2} \sum_{E \in \mathcal E_H \setminus \mathcal E_h}\eta_H^2+\varLambda_2E_H^2 \nonumber\\
& \leq &(1+\delta_2)(1-\frac{\theta_A}{2}) \eta_H^2 + \varLambda_2 E_H^2.
\end{eqnarray*}
%The choice of $\delta_2 ={\theta_A}/{4-\theta_A}$
For any given  $\theta_A$, the choice of $\delta_2={\theta_A}/{(4-2\theta_A)}$ implies $0<\rho_A=(1+\delta_2)(1-\frac{\theta_A}{4})<1$.  Hence, (\ref{et_rd2}) holds. 
For  the {\bf Case (B)},  a summation over the all edges  implies 
\begin{equation*}\label{et_rd3}
 \eta_h^2 \leq (1+\delta_2) \eta_H^2 + \varLambda_2 E_H^2.
\end{equation*}
This completes the rest of the proof. \hfill${\Box}$ %\qquad \endproof % $\hfill{\Box}$
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\begin{remark}~\label{rmk2} Note that Equation (\ref{edge_4}) irrespective of the marking criteria for the refinement.
% Hence,  for these cases selection of  $\delta_2=1$  gives following result
% %hold for any refinement, which is independent of the marking criteria. Hence,  for these cases following result holds
%\begin{eqnarray}\label{et_rd3}
%\eta_h^2 \leq 2 \eta_H^2 + \varLambda_2 E_H^2.
%\end{eqnarray}
%\end{remark}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\begin{lm}\label{lm5} Let $({\bf p_h}, u_h)$ and $({\bf p_H}, u_H)$ be the solution of  (\ref{eqna1})-(\ref{eqna2}) with 
%the nested triangulation $\mathcal T_h$ and $\mathcal T_H$, then for sufficiently small mesh-size $h_1$,  it holds
% \begin{eqnarray}\label{ortho_5}
% 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|&+&2|(u-u_h,u_h-u_H)|\leq \alpha_1 e_h^2+\alpha_2 E_H^2+\varLambda_3 \mu^2_H,
%\end{eqnarray}
%%  \begin{eqnarray}
%%  2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|&+&2|(u-u_h,u_h-u_H)|\leq (\delta_3+\frac{2C_1^2\epsilon^2}{\delta_4})(e_{p_h}^2+e_{u_h}^2) \nonumber\\
%% &&~~~+(\delta_4+\frac{4C^2}{\delta_4}\varLambda_1 H^2) E_H^2+(\frac{8C^2}{\delta_4}+\delta_4+\frac{C^2}{2\delta_4})\mu^2_H
%% \end{eqnarray}
%where the constants $0<\alpha_1,\alpha_2<1$ and $0<\varLambda_3$ shall be specified in the following proof.
%\end{lm}
%{\it Proof}.  Introduce an intermediate term ~${\bf \tilde p}= -({\bf A}\nabla_{NC}u_h^N+u_h{\bf b})$ ~in ~${\bf  (p- p_h)}$, for the estimates of the first term on left-hand side of (\ref{ortho_5})% $({{\bf A^{-1} (p- p_h), p_h-p_H}})$ 
%% Lets  find first bound for $({{\bf A^{-1} (p- p_h), p_h-p_H}})$ by %and $({{ u- u_h, u_h-u_H}})$
%% introducing an intermediate term ${\bf \tilde p}$ in ${\bf  (p- p_h)}$ 
% \begin{eqnarray}\label{ortho1}
%({{\bf  A^{-1} (p- p_h), p_h-p_H}}) &=& ({{\bf A^{-1} (p- \tilde p), p_h-p_H}})+({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}})\nonumber\\
%&=& -(\nabla_{NC} (u-u_h^N), {\bf p_h- p_H})-((u-u_h){\bf b^*},{\bf p_h-p_H})\nonumber\\
%&&+({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}}).
% \end{eqnarray}
%% For the bound of first term, integration by part elementwise, result in 
%An elementwise integration by part in the first term of right-hand side yields
%\begin{eqnarray*}\label{ortho2}
%-(\nabla_{NC} (u-u_h^N), {\bf p_h- p_H})= ((u-u_h^N), \text{div}~ ({\bf p_h- p_H}))-\sum_{E\in \mathcal E_h} \int_E [u-u_h^N]({\bf p_h -p_H})\cdot \nu_E ds.
%\end{eqnarray*}
%As $({\bf p_h -p_H})\in RT^1_0(\mathcal T_h),$ it follows that $({\bf p_h -p_H})\cdot \nu_E$ is continuous along edge $E$ and constant on $E$ and $\int_E [u-u_h] ds =0$.
%The second equation of weak formulation (\ref{eqna2}) and the definition of $L^2$-projection gives
%\begin{eqnarray*}\label{ortho3}
% (u-u_h^N, \text{div}~ ({\bf p_h- p_H})) = ((u-u_h^N)-\Pi_H (u-u_h^N), f_h-f_H)-(u-u_h^N,\gamma (u_h-u_H)).
%\end{eqnarray*}
%Substitution of the combination of above identities in (\ref{ortho1}),   the Cauchy-inequality and the Poincar$\acute{e}$ inequality result in 
%\begin{eqnarray*}
%({{\bf  A^{-1} (p- p_h), p_h-p_H}}) &\lesssim &   \tnorm{u-u_h^N}_{NC} \norm{H(f_h-f_H)}+\norm{u-u_h^N} \norm{u_h-u_H}  \nonumber\\
%&&+\norm{{\bf A^{-1} (\tilde p- p_h)}} \norm{{\bf p_h- p_H}}+\norm{u-u_h} \norm{{\bf A^{-1} (p_h-p_H)}}.
%\end{eqnarray*}
%% \begin{eqnarray}
%% ({\bf b^*}(u-u_h),{\bf p_h-p_H})&\lesssim& \norm{u-u_h} \norm{{\bf A^{-1} (p_h-p_H)}}\nonumber\\
%% &\lesssim& \big( {osc}_h(f)+\norm{[{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*]\cdot t_E}_{\mathcal E}\nonumber\\  
%% &&+\norm{h_{\mathcal{T}}(\mathbf{A}^{-1}{\bf p_h}+u_h{\bf b}^*)} + \norm{h_{\mathcal{T}}\;(f_h-\gamma\;u_{M})}\big) \norm{{\bf A^{-1}( p_h-p_H)}} \nonumber\\
%% (u-u_h^N, f_h-f_H)&\lesssim& ((u-u_h^N)-\Pi_H (u-u_h^N), f_h-f_H)\nonumber\\
%% &\lesssim& \tnorm{u-u_h^N} \norm{H(f_h-f_H)}\nonumber\\
%% &\lesssim& (e_{p_h}+e_{u_h} +\norm{h(f_h-\gamma u_h)})\norm{H(f_h-f_H)}\nonumber\\
%% (u-u_h^N,\gamma (u_h-u_H))&\lesssim &\norm{u-u_h^N} \norm{u_h-u_H} \nonumber\\
%% &\lesssim & \big(osc_h (f)+\epsilon(e_{p_h}+e_{u_h} +\norm{h(f_h-\gamma u_h)})\big) \norm{u_h-u_H} \nonumber\\
%% ({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}})&\lesssim &\norm{{\bf A^{-1} (\tilde p- p_h)}} \norm{{\bf (p_h- p_H)}}\nonumber\\
%% &\lesssim &\norm{h(f_h-\gamma u_h)} \norm{{\bf (p_h- p_H)}}
%% \end{eqnarray}
%%\vspace{4cm}
%Using the estimates (\ref{esti_1})-(\ref{esti_8}) and the addition of the term $(u-u_h,u_h-u_H)$ yields% combining all equation term
%%
%\begin{eqnarray*}
%&& 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq C_3 (e_{p_h}+e_{u_h})osc_H^2 (f)\nonumber\\
%&&~~+C_3(\norm{h f_h}+\norm{h( {\bf A^{-1} p_h}+ u_h{\bf b^*})})osc_H^2 (f) +C_3\epsilon (e_{p_h}+e_{u_h})(\norm{u_h-u_H}+\norm{\bf A^{-1/2}(  p_h-p_H)})\nonumber\\
%&&~~+C_3(osc_h (f)+\norm{hf_h}+\norm{h{\bf A^{-1} p_h}+ u_h{\bf b^*}})(\norm{u_h-u_H}+\norm{\bf A^{-1/2}(  p_h-p_H)}).
%\end{eqnarray*}
%Then, by  the Young's inequality and the simple calculation  result in 
%%for some $0<\delta_3$  rearrengement of the terms results in
%% \begin{eqnarray}
%% && 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq (\delta_3+{8C_1^2\epsilon^2})(e_{p_h}^2+e_{u_h}^2)+\frac{1}{2} E_H^2 \nonumber\\
%% &&~~+\frac{C^2}{\delta_4}\norm{h(f_h-\gamma u_h)}^2+(\delta_4+\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2 \nonumber\\
%% &&~~+\frac{3C^2}{\delta_4}(osc_h^2 f+\norm{h(f_h-\gamma u_h)}^2+\norm{{\bf A^{-1} p_h+b^*} u_h}^2.)%+(\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2
%%  \end{eqnarray}
%\begin{eqnarray*}
% &&2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|  \nonumber\\
%&&\qquad\qquad\leq (\delta_3+{4C_3^2\epsilon^2})(e_{p_h}^2+e_{u_h}^2)+\frac{1}{2} E_H^2
%+7C_3^2 \mu_h^2+(\frac{1}{2}+\frac{C_3^2}{2\delta_3})osc_H ^2(f).
%%&&~~+\frac{3C^2}{\delta_4}(osc_h^2 f+\norm{h(f_h-\gamma u_h)}^2+\norm{{\bf A^{-1} p_h+b^*} u_h}^2.)%+(\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2
% \end{eqnarray*}
%% use (\ref{l{et_rd1}} )  to result for
%% % \begin{eqnarray}
%% % && 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq (\delta_3+\frac{2C_1^2\epsilon^2}{\delta_4})(e_{p_h}^2+e_{u_h}^2)+\delta_4 E_H^2 \nonumber\\
%% % &&+C\norm{h(f_h-\gamma u_h)}\norm{h(f_h-f_H)} \nonumber\\
%% % &&+\frac{3C^2}{\delta_4}(osc_h^2 f+\norm{h(f_h-\gamma u_h)}^2+\norm{{\bf A^{-1} p_h+b^*} u_h})^2+(\frac{C^2}{2\delta_2})\norm{h(f_h-f_H)}
%% %  \end{eqnarray}
%% \begin{eqnarray}
%% && 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq (\delta_3+\frac{2C_1^2\epsilon^2}{\delta_4})(e_{p_h}^2+e_{u_h}^2)+\delta_4 E_H^2 \nonumber\\
%% &&~~~+(\frac{4C^2}{\delta_4})\mu_h^2+(\delta_4+\frac{C^2}{2\delta_4})\mu^2_H
%%  \end{eqnarray}
%As $\mathcal T_h$ and $\mathcal T_H$ are the nested triangulations, Remark \ref{rmk3}  implies that $\mu_h^2 \leq  2\mu_H^2 +\varLambda_1 H^2 E_H^2$. % (\ref{et_rd1}),
%And hence, with notations  $0<\alpha_1=(\delta_3+{4C_3^2\epsilon^2}), ~\alpha_2= (\frac{1}{2}+{7C_3^2}\varLambda_1 H^2)$  and $\varLambda_3 =({14C_3^2}+\frac{1}{2}+\frac{C_3^2}{2\delta_3})$ above equation reduce to
%\begin{eqnarray*}
%&& 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq \alpha_1(e_{p_h}^2+e_{u_h}^2)+ \alpha_2 E_H^2+\varLambda_3\mu^2_H.~~~~~~ \hfill{\Box}
%\end{eqnarray*}
%%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\begin{lemma}\label{lm6}(Quasi-orthogonality)
Let $\mathcal T_h$  be  a refined triangulation of $\mathcal T_H$. Then for  small  initial mesh-size $h\leq h_2$, there exist constants $0<\alpha_1 ,\alpha_3<1$ such that
\begin{eqnarray}\label{er_re}
(1-\alpha_1){ e_{h}}^2\leq { e_H}^2-\alpha_3 E_H^2+ \varLambda_3 \mu_H^2.
\end{eqnarray}
%where the constants $0<\alpha_1,\alpha_3<1$ and $0<\varLambda_3$ shall be specified in the following proof.
\end{lemma}
{\it Proof}. The following hold %Let start with the term $\norm{{\bf A^{-1/2} (p- p_h)}}$ and $\norm{{u- u_h}}^2$, 
\begin{align}
~\norm{{\bf A^{-1/2} (p- p_h)}}^2&=\norm{{\bf A^{-1/2} (p- p_H)}}^2-\norm{{\bf A^{-1/2} (p_h- p_H)}}^2\nonumber \\
&~~-2\big({{\bf A^{-1} (p- p_h), p_h-p_H}}\big),\label{nov2}\\
\qquad\qquad \norm{{u- u_h}}^2&=\norm{{u- u_H}}^2-\norm{{ u_h- u_H}}^2-2({{ u- u_h, u_h-u_H}}).\label{nov3}
 \end{align}
%Similarly
%{\it Proof}. 
The term $\big({{\bf A^{-1} (p- p_h), p_h-p_H}})$ is estimated  by introducing the intermediate term ~${\bf \tilde p}= -({\bf A}\nabla_{NC}u_h^N+u_h{\bf b})$ as  %~in ~${\bf  (p- p_h)}$,  %the first term on left-hand side of (\ref{ortho_5})% $({{\bf A^{-1} (p- p_h), p_h-p_H}})$ 
% Lets  find first bound for $({{\bf A^{-1} (p- p_h), p_h-p_H}})$ by %and $({{ u- u_h, u_h-u_H}})$
% introducing an intermediate term ${\bf \tilde p}$ in ${\bf  (p- p_h)}$ 
 \begin{eqnarray}\label{ortho1}
({{\bf  A^{-1} (p- p_h), p_h-p_H}}) &=& ({{\bf A^{-1} (p- \tilde p), p_h-p_H}})+({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}})\nonumber\\
&=& -(\nabla_{NC} (u-u_h^N), {\bf p_h- p_H})-((u-u_h){\bf b^*},{\bf p_h-p_H})\nonumber\\
&&+({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}}).
 \end{eqnarray}
% For the bound of first term, integration by part elementwise, result in 
An elementwise integration by parts of the first term on the  right-hand side of (\ref{ortho1})  yields
\begin{align}\label{ortho2}
-(\nabla_{NC}(u-u_h^N), {\bf p_h- p_H})&=(u-u_h^N, \text{div}~ ({\bf p_h- p_H}))\nonumber\\
&~~~-\sum_{E\in \mathcal E_h} \int_E [u-u_h^N]({\bf p_h -p_H})\cdot \nu_E ~ds.
\end{align}
%As $({\bf p_h -p_H})\in RT^1_0(\mathcal T_h),$ it follows that $({\bf p_h -p_H})\cdot \nu_E$ is continuous along edge $E$ and constant on $E$ and $\int_E [u-u_h^N] ds =0$.
The second term on right-hand side of (\ref{ortho2}) is zero, as $({\bf p_h -p_H})\cdot \nu_E$ is continuous along the edge $E$ and constant on $E$, and  $u_n^N\in CR_0^1 (\T)$, $\int_E [u-u_h^N] ds =0$. 
The second equation of weak formulation (\ref{eqna2}) and the definition of $L^2$-projection show
\begin{align}\label{ortho3}
 (u-u_h^N, & \text{div}~ ({\bf p_h- p_H})) \nonumber \\
 &= ((u-u_h^N)-\Pi_H (u-u_h^N), f_h-f_H)-(u-u_h^N,\gamma (u_h-u_H)).
\end{align}
A substitution of (\ref{ortho2})-(\ref{ortho3}) in (\ref{ortho1}) with a use of   the Cauchy-Schwarz inequality and the Poincar$\acute{e}$ inequality results in 
\begin{align*}
&({{\bf  A^{-1} (p- p_h), p_h-p_H}}) \lesssim   \norm{\nabla_{NC} (u-u_h^N)} \norm{H(f_h-f_H)}+\norm{u-u_h^N} \norm{u_h-u_H}  \nonumber\\
&\hspace{10mm}+\norm{u-u_h} \norm{{\bf A^{-1/2} (p_h-p_H)}}+\norm{{\bf \tilde p- p_h}} \norm{{\bf A^{-1/2} (p_h-p_H)}}.
\end{align*}
% \begin{eqnarray}
% ({\bf b^*}(u-u_h),{\bf p_h-p_H})&\lesssim& \norm{u-u_h} \norm{{\bf A^{-1} (p_h-p_H)}}\nonumber\\
% &\lesssim& \big( {osc}_h(f)+\norm{[{\mathbf A}^{-1}{\bf p_h}+u_h{\bf b}^*]\cdot t_E}_{\mathcal E}\nonumber\\  
% &&+\norm{h_{\mathcal{T}}(\mathbf{A}^{-1}{\bf p_h}+u_h{\bf b}^*)} + \norm{h_{\mathcal{T}}\;(f_h-\gamma\;u_{M})}\big) \norm{{\bf A^{-1}( p_h-p_H)}} \nonumber\\
% (u-u_h^N, f_h-f_H)&\lesssim& ((u-u_h^N)-\Pi_H (u-u_h^N), f_h-f_H)\nonumber\\
% &\lesssim& \tnorm{u-u_h^N} \norm{H(f_h-f_H)}\nonumber\\
% &\lesssim& (e_{p_h}+e_{u_h} +\norm{h(f_h-\gamma u_h)})\norm{H(f_h-f_H)}\nonumber\\
% (u-u_h^N,\gamma (u_h-u_H))&\lesssim &\norm{u-u_h^N} \norm{u_h-u_H} \nonumber\\
% &\lesssim & \big(osc_h (f)+\epsilon(e_{p_h}+e_{u_h} +\norm{h(f_h-\gamma u_h)})\big) \norm{u_h-u_H} \nonumber\\
% ({{\bf A^{-1} (\tilde p- p_h), p_h-p_H}})&\lesssim &\norm{{\bf A^{-1} (\tilde p- p_h)}} \norm{{\bf (p_h- p_H)}}\nonumber\\
% &\lesssim &\norm{h(f_h-\gamma u_h)} \norm{{\bf (p_h- p_H)}}
% \end{eqnarray}
%\vspace{4cm}
Using the estimates (\ref{esti_1})-(\ref{esti_8}), Lemma \ref{lm8} and  the addition of the term
\\ $(u-u_h,u_h-u_H)$ yield% combining all equation term
%
\begin{eqnarray*}
&& 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq C_3 (e_{p_h}+e_{u_h})osc_H (f)\nonumber\\
&&~~+C_3(\norm{h f_h}+\norm{h( {\bf A^{-1} p_h}+ u_h{\bf b^*})})osc_H (f) +C_3\epsilon (e_{p_h}+e_{u_h})(\norm{u_h-u_H}\nonumber\\ 
&&~~+\norm{{\bf A^{-1/2}(  p_h-p_H)}})+C_3(osc_h (f)+\norm{hf_h}+\norm{h{\bf A^{-1} p_h}+ u_h{\bf b^*}})\nonumber \\
&&~~(\norm{u_h-u_H}+\norm{\bf A^{-1/2}(  p_h-p_H)}).
\end{eqnarray*}
The Young's inequality and rearrangement of terms result in 
%for some $0<\delta_3$  rearrengement of the terms results in
% \begin{eqnarray}
% && 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|\leq (\delta_3+{8C_1^2\epsilon^2})(e_{p_h}^2+e_{u_h}^2)+\frac{1}{2} E_H^2 \nonumber\\
% &&~~+\frac{C^2}{\delta_4}\norm{h(f_h-\gamma u_h)}^2+(\delta_4+\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2 \nonumber\\
% &&~~+\frac{3C^2}{\delta_4}(osc_h^2 f+\norm{h(f_h-\gamma u_h)}^2+\norm{{\bf A^{-1} p_h+b^*} u_h}^2.)%+(\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2
%  \end{eqnarray}
\begin{eqnarray}\label{pp1}
 &&2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+2|(u-u_h,u_h-u_H)|  \nonumber\\
&&\qquad\qquad\leq (\delta_3+{4C_3^2\epsilon^2})(e_{p_h}^2+e_{u_h}^2)+\frac{1}{2} E_H^2
+7C_3^2 \mu_h^2+\big(\frac{1}{2}+\frac{C_3^2}{2\delta_3}\big)osc_H ^2(f).
%&&~~+\frac{3C^2}{\delta_4}(osc_h^2 f+\norm{h(f_h-\gamma u_h)}^2+\norm{{\bf A^{-1} p_h+b^*} u_h}^2.)%+(\frac{C^2}{2\delta_4})\norm{h(f_h-f_H)}^2
 \end{eqnarray}
As $\mathcal T_h$ and $\mathcal T_H$ are the nested triangulations, Remark \ref{rmk3}  implies that \\ $\mu_h^2 \leq  2\mu_H^2 +\varLambda_1 H^2 E_H^2$. % (\ref{et_rd1}).
Hence, with notations  $0<\alpha_1=(\delta_3+{4C_3^2\epsilon^2}),$ $\alpha_2= (\frac{1}{2}+{7C_3^2}\varLambda_1 H^2)$  and $\varLambda_3 =({14C_3^2}+\frac{1}{2}+\frac{C_3^2}{2\delta_3})$, (\ref{pp1})  reduces to
\begin{align}\label{nov1}
 2|({{\bf A^{-1} (p- p_h), p_h-p_H}})|+& 2|(u-u_h,u_h-u_H)| \nonumber \\
 & \leq \alpha_1(e_{p_h}^2+e_{u_h}^2)+ \alpha_2 E_H^2+\varLambda_3\mu^2_H.% \hfill{\Box}
\end{align}
A combination of (\ref{nov2})-(\ref{nov3}) with (\ref{nov1}) yields
\begin{eqnarray}
({ e_{p_h}}^2+e_{u_h}^2)\leq ({ e_{p_H}}^2+e_{u_H}^2)- E_H^2+\alpha_1(e_{p_h}^2+e_{u_h}^2)+ \alpha_2 E_H^2+\varLambda_3\mu^2_H.\nonumber
 \end{eqnarray}
 %$A_0=(\delta_3+\frac{2C_1^2\epsilon^2}{\delta_4})<1$ and  -$(1-)$,  hence using defined notation of $e_h$ above equation reduced to
%Combining first two equation with above
Hence $ (1-\alpha_1)e_{h}^2\leq e_{H}^2-(1-\alpha_2) E_H^2+ \varLambda_3 \mu_H^2.$\\
%\begin{equation*}
% (1-\alpha_1)e_{h}^2\leq e_{H}^2-(1-\alpha_2) E_H^2+ \varLambda_3 \mu_H^2.
%\end{equation*}
 For any $\epsilon >0$,  
it always possible to find  a  small initial mesh-size  $h_2$ and $0<\delta_3<1$, such that
$0<\alpha_1=\delta_3+4C_3^2\epsilon^2<1$ and $0<\alpha_3=1-\alpha_2=1/2-7C_3^2\varLambda_1 H^2$. % where $\alpha_4 $ will be chosen later.
Hence, (\ref{er_re}) holds true and this completes the rest of the proof.\hfill${\Box}$ % \qquad \endproof %\hfill${\Box}$
%
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
\begin{remark} \label{rmk2} For proving the contraction property, a more sharper bound for \\
$\alpha_1 <\alpha_4:={1}/{2}\min\{1, {\alpha_3(1-\rho_A)}/({(\varLambda_1+\varLambda_2)C_{rel}})\}<1$ will be selected. 
% For any $\epsilon >0$ and given $0<\alpha_4<1$,  
%it always possible to find  a sufficiently small initial mesh-size  $h_1$ and $0<\delta_3<1$, such that
%$0<\alpha_1=\delta_3+4C_3^2\epsilon^2<\alpha_4<1$ and $0<\alpha_3=1-\alpha_2=1/2-7C_3^2\varLambda_1 H^2$. % where $\alpha_4 $ will be chosen later.
%Hence (\ref{er_re}) holds true. \hfill${\Box}$
\end{remark}
%For any selected $0<\delta_3<1$ and $\epsilon^2$ there exist suffientily small mesh-size $h_1$ such that 
% $\alpha_1<1$ and $\alpha_3=1-\alpha_2<1$, and hence (\ref{er_re}) holds. %(=(\frac{1}{2}-{8C_3^2}\varLambda_1 H^2))
%here $\varLambda_3 =  1-(\delta_4+\frac{4C^2}{\delta_4}\varLambda_1 H^2)~~      \varLambda_4 =$
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~``
 \begin{corollary}\label{coro3} Under the assumption that  small initial mesh-size $h_2 >0$, and   the constants  defined in Lemma \ref{lm6}, the following result holds for $0<h\leq h_2$
\begin{eqnarray}
 {\bf e_{H}}^2&\leq & (1+\alpha_1){\bf e_h}^2+(1+\alpha_2) E_H^2+ \varLambda_3 \mu_H^2\nonumber \\
 &\leq&  2{\bf e_h}^2+ 2 E_H^2+ \varLambda_3 \mu_H^2.
\end{eqnarray}
%An assumtion of initial mesh-size $h_0$ and $\delta_3$ implies $\alpha_1, \alpha_2 <1 $, and hence
%\begin{equation}
% {\bf e_{H}}^2\leq  2{\bf e_h}^2+ 2 E_H^2+ \varLambda_3 \mu_H^2.
%\end{equation}
 \end{corollary}
 % Proof-
%  begin{eqnarray}
% \norm{{\bf A^{-1/2} (p- p_H)}}^2=\norm{{\bf A^{-1/2} (p- p_h)}}^2+\norm{{\bf A^{-1/2} (p_h- p_H)}}^2+2({{\bf A^{-1} (p- p_h), p_h-p_H}})
%  \end{eqnarray}
% And above procedure
%This concludes the proof. 
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\begin{lemma}({Quasi-discrete reliability})\label{coro2}% Let $\mathcal  T_H$  and  $\mathcal T_h$ are the nested triangulation. Then,
% with constant  $C_5>0$ 
~Let $({\bf p}_h, u_h )$ and $({\bf p}_H, u_H)$ be the MFEM solutions of (\ref{eqna1})-(\ref{eqna2}) over the triangulations $\T_h$ and $\T_H$, respectively.  There exists a constant $C_4 >0$ such that for any $\epsilon >0$ 
%from (\ref{dual})  %Using quasi-discrete reliability result
 \begin{align*}
\norm{{\bf A^{-1/2} (p_h -p_H)}}^2+\norm{u_h-u_H}^2 &\leq C_4\big( \eta_H^2(\mathcal E_H \setminus \mathcal E_h) + \norm{H (f_h- f_H)}^2\\
&~~~ +  \epsilon^2 ( e_h^2+  e_H^2) +\mu_h^2 +\mu_H^2\big),
\end{align*}
where $E_H, e_h, e_H$ are defined  in (\ref{com_1}) and $\mu_h$ in (\ref{mu}).
%it holds
%
%\begin{eqnarray}
%\norm{{\bf p_h -p_H}}^2 \lesssim \eta_H^2(\mathcal E_H \setminus \mathcal E_h) + \norm{H (f_h- f_H)}^2  + \epsilon^2 ( e_h^2+  e_H^2) +\mu_h^2 +\mu_H^2.
%\end{eqnarray}
\end{lemma}\\
{\it Proof}.
%  Let $({\bf p_H},u_H)$ and $u_H^N$ be the solution of  (\ref{eqna1})-(\ref{eqna2}) and (\ref{eqna3})  respectively with the  right-hand side $f_H$ over the  triangulation $\mathcal T_H$.
%Similarly,  $({\bf p_h},u_h)$ and $u_h^N$ are the solutions of  (\ref{eqna1})-(\ref{eqna2}) and (\ref{eqna3}) respectively  with the right-hand side $f_h$ over the triangulation  $\mathcal T_h$.
Introduce the discrete mixed finite element problem: seek
 $ (\tilde{\bf p}_h, \tilde u_h) \in RT_0(\mathcal {T}_h)\times P_0(\mathcal T_h)$ such that
  \begin{eqnarray}
  &&(\mathbf{A}^{-1}\tilde {\bf p}_h,{\bf q_h})-(\text{div}~{\bf q_h},\tilde u_h)=-( u_H {\bf b}^*,{\bf q_h}) \qquad 
  \text{for all}\, {\bf q_h}\in RT_0(\mathcal {T}_h), \label{eqna8}\\
&&  (\text{div}~ \tilde{\bf p}_h,v_h)=(f_H,v_h)- (\gamma u_H, v_h) \qquad  \text {for all} \, v_h \in P_0(\mathcal{T}_h). \label{eqna9}
 \end{eqnarray}
  %where $ u_h^N\in CR^1_0(\mathcal T)$ 
 %is the non-conforming finite element approximation which 
 Define  the nonconforming discrete  problem  corresponding to (\ref{eqna8})-(\ref{eqna9}): seek  $\tilde u_h^N\in CR_0^1(\mathcal T_h)$ as the solution of
 % CR- FE space $(\mathcal T)$
 \begin{equation}
  (\mathbf{A}\nabla_{NC} \tilde u_h^N+{u}_H {\bf b},\nabla_{NC} v_{CR})+(\gamma {u}_H,v_{CR})=
  (f_H,v_{CR}), {\forall} v_{CR}\in CR^1_0(\mathcal T_h).\label{eqna10}
 \end{equation}
 %\begin{eqnarray} 
% \rm{where} ~~~\displaystyle  \tilde u_h({\bf x})=\left(1+\frac{S(T)}{4} \gamma \right)^{-1} \left(\Pi_0 u_H^N+\frac{S(T)}{4}f_H \right)  \;\; \text{for}~ x \in T\in \mathcal T_h\label{umt1}.
 %&& ~~S(T)= \displaystyle{\int_T} ({\bf x}-\text{mid}(T))\cdot \mathbf{A}^{-1}({\bf x}-\text{mid}%(T)) \,d{\bf x}\;\; \text{for} ~ T\in \mathcal T_H,\label{s}
% \end{eqnarray}
  %\begin{eqnarray} 
% && \displaystyle  \tilde u_h({\bf x})=\left(1+\frac{S(T)}{4} \gamma \right)^{-1} \left(\Pi_0 u_H^N+\frac{S(T)}{4}f_H \right)  \;\; \text{for}~ x \in T\in \mathcal T,\label{umt1}
% %&& ~~S(T)= \displaystyle{\int_T} ({\bf x}-\text{mid}(T))\cdot \mathbf{A}^{-1}({\bf x}-\text{mid}%(T)) \,d{\bf x}\;\; \text{for} ~ T\in \mathcal T_H,\label{s}
% \end{eqnarray}
The solution ${\tilde {\bf p}_h}$ of (\ref{eqna8})-(\ref{eqna9}) can be written in the terms of $\tilde u_h^N$ as
 \begin{equation} \label{pmt1}
  \tilde {\bf p}_h({\bf x})=-\left(\mathbf{A}\nabla _{NC} \tilde u_h^N+u_H {\bf b} \right)+\left(f_H-\gamma  u_H\right)
  \frac{\left({\bf x}-\text{\rm mid}(T)\right)}{2} \text{ for}~~ {\bf x}\in T.% \in \mathcal{T}_h.
 \end{equation}
 Now for the estimates of $\norm{{\bf p_h- p_H }}$,  use ${\tilde {\bf p}_h}$ as an intermediate term to split  ${\bf p_h}-{\bf p_H}:= ({\bf p_h}-{\bf {\tilde p_h}})+({\bf {\tilde p_h}}-{\bf p_H})$ and  the triangle 
 inequality. 
  From the representation formula (\ref{pmt}) and (\ref{pmt1}) of ${\bf p_h}$ and ${\bf \tilde p_h }$, respectively, it follows that
\begin{eqnarray*}
{{\bf p_h- \tilde p_h }}= - ({\bf A} \nabla (u_{h}^N-\tilde u_h^N)+ (u_h -u_H){\bf b})+\frac{1}{2} (f_h-f_H-
\gamma (u_h-u_H))(x- {\rm mid} (T)).
 \end{eqnarray*}
 The triangle inequality shows 
\begin{align}\label{rea_4}
\norm{  {\bf p_h- \tilde p_h } } \lesssim \norm{{\bf A} \nabla (u_{h}^N-\tilde u_{h}^N)}+\norm{ u_h -u_H}+\norm{h(f_h-f_H)}. %+ \norm{h (u_h-u_H)}
\end{align}
%The difference of the nonconforming systems 
Subtracting (\ref{eqna3}) from (\ref{eqna10}) leads to
 \begin{align}\label{rea_6}
({\bf A} \nabla (u_{h}^N-\tilde u_{h}^N), \nabla v_{CR})&= (f_h-f_H,v_{CR})+(\gamma(u_h-u_H),v_{CR})\nonumber \\
&~~-((u_h-u_H){\bf b},\nabla v_{CR}).
 \end{align}
A substitution $v_{CR}={\bf A} (u_{h}^N-\tilde u_{h}^N)$ in (\ref{rea_6}) with the $L^2$-projection property and the Cauchy-Schwarz inequality results in%(\ref{rea_6})
$$\norm{{\bf A}\nabla (u_{h}^N-\tilde u_{h}^N)}\lesssim \norm{H(f_h-f_H)}+\norm{u_h-u_H}.$$ 
With this estimate, (\ref{rea_4}) reduces to 
 \begin{eqnarray}\label{rea_9}
\norm{  {\bf p_h- \tilde p_h } } \lesssim\norm{H(f_h-f_H)}+ \norm{ u_h -u_H}.\label{rea_7}
\end{eqnarray}
%\norm{{\bf p_h- p_H}} &\leq & \norm{{\bf p_h -\tilde p_H}}+\norm{{\bf \tilde p_h -p_H}} \nonumber \\
%&\lesssim & \eta_H (\mc E_H \setminus \mc E_h) +\norm {h (f_h-f_H)} +\norm{u_h-u_H}
% \end{eqnarray}
For a bound of  the term $ (u_h-u_H),$ a use of (\ref{esti_8}) shows
 \begin{eqnarray}\label{rea_8}
\norm{u_h-u_H}^2 &\lesssim& \norm{u_h-u}^2+ \norm {u-u_H}^2 \nonumber \\
&\lesssim & \epsilon^2 e_h^2 +\mu_h^2 + \epsilon^2 e_H^2+\mu_H^2.
 \end{eqnarray}
For the estimate of $\norm{{\bf p_h- p_H }}$, note that
 %Now $\norm{{\bf p_h- p_H }}\leq \norm{ {\bf p_h -\tilde p_h}} +\norm{{\bf \tilde p_h- p_H}}$
${\rm div}~ {\bf \tilde p_h} =f_H-\gamma u_H={\rm div}~ {\bf p_H}$. It implies $\rm {div} ~(\tilde {\bf p}_h-{\bf p}_H)=0$, and hence,  
$ {\bf \tilde p_h- p_H}$ is a piecewise constant vector function over $\mathcal T_h$.
%The discrete Helmholtz decomposition \cite{arnold} states for ${\bf A( p_h- p_H) }$ equals  a discrete gradient ${\mathbf A}\nabla_{NC} \alpha_{CR}$ of some nonconforming function $\alpha_{CR}\in CR_0^1(\mathcal T_h)$ plus the ${\rm Curl}~ \beta_h$ of some piecewise affine conforming function $\beta_h \in P_1(\mathcal T_h)\cap C(\bar \Omega)$. That is,
% ${\bf A (\tilde p_h-p_H)= A}\nabla \alpha_{CR}+ \rm {Curl}~\beta_h$, and hence %for $\alpha_{CR}\in CR^1_0(\mathcal T_h)$
  %and $\beta_c \in P_1(\mc T_h) $ 
  The discrete Helmholtz decomposition  states for % \cite{arnold} ${\bf A( p_h- p_H) }$ equals  a discrete gradient ${\mathbf A}\nabla_{NC} \alpha_{CR}$ of some nonconforming function $\alpha_{CR}\in CR_0^1(\mathcal T_h)$ plus the ${\rm Curl}~ \beta_h$ of some piecewise affine conforming function $\beta_h \in P_1(\mathcal T_h)\cap C(\bar \Omega)$. That is,
 ${\bf A (\tilde p_h-p_H)= A}\nabla \alpha_{CR}+ \rm {Curl}~\beta_h$, where $\alpha_{CR}\in CR_0^1(\mathcal T_h)$ and $\beta_h \in P_1(\mathcal T_h)\cap C(\bar \Omega)$, and hence, %for $\alpha_{CR}\in CR^1_0(\mathcal T_h)$
   \begin{eqnarray*}
 \norm{{\bf \tilde p_h -p_H}}^2=({\bf \tilde p_h -p_H}, \nabla \alpha_{CR}+ A^{-1} \rm {Curl}~ \beta_h ).
 \end{eqnarray*}
 Let $\beta_H := I_H \beta_h$ be the Scott-Zhang quasi-interpolation operator  for  $E \in \mc E_H$, %\cite{scottzhang} 
 where $\mathcal E_H$ the set of  edges on the triangulation $\mathcal T_H$ and its neighbourhood $\omega_E$ with
  \begin{equation}\label{drel22}
  \norm{\beta_h- \beta_H}_E\leq C h_E^{{1}/{2}} \norm{\beta_h}_{H^1{(\omega_E})}.
  \end{equation}
 Note that $\norm{\beta_h- \beta_H}_E=0 $ if $E \in \mc E_h\cap\mc E_H.$\\
 %The $L^2$-orthogonality
 The fact $(({\bf \tilde p_h -p_H}), \nabla_{NC} \alpha_{CR})=0$ shows~$ \norm{{\bf \tilde p_h -p_H}}^2=({\bf \tilde p_h -p_H, A^{-1}}\rm {Curl}~ \beta_h )$. % First term of right-hand side
  %\begin{eqnarray}
%({\bf \tilde p_h -p_H}, \nabla \alpha_{CR})=-(div({\bf  p_h -p_H}), \nabla \alpha_{CR})+\int_E [\alpha_{CR}]({\bf \tilde p_h -p_H}\cdot \mu)=0
% \end{eqnarray}
 %for the second term
The weak formulation (\ref{eqna1}) with  ${\bf q_H}= {\rm Curl}~ \beta_H \in RT_0(\mc T_H)\subset RT_0(\mc T_h)$ over $\mc T_H$ and (\ref{eqna8})  with ${\bf q_h}= {\rm Curl}~ \beta_H \in RT_0(\mc T_h)$ and integration by parts lead to
% \begin{eqnarray}
% %\norm{{\bf \tilde p_h -p_H}}^2=
% ({\bf \tilde p_h -p_H, A^{-1}}\rm {Curl}~ \beta_h )&=&({\bf A^{-1} p_h}, \rm {Curl}~ \beta_h)-({\bf A^{-1} p_H}, \rm {Curl}~ \beta_h)\nonumber \\
%&\leq & -(b^* u_H, \rm {Curl}~ \beta_h)-({\bf A^{-1} p_H, \rm {Curl}~ \beta_h})
% \end{eqnarray}
\begin{eqnarray}\label{rea_2}
 \norm{{\bf \tilde p_h -p_H}}^2& = & -( u_H{\bf b^*}, \rm {Curl}~ \beta_h) - ({\bf A^{-1}p_H}, \rm {Curl}~\beta_h)\nonumber \\
 &= &({\bf A^{-1} p_H}+ u_H {\bf b^*}, {\rm Curl}~ (\beta_H-\beta_h)) \nonumber\\%[int by part]
& = & \sum_{E \in \mc E_h} \int_E [{\bf A^{-1}p_H+b^*}u_H]\cdot \tau_E (\beta_H-\beta_h) ~ds \nonumber\\
&&~-\sum_{T \in \mc T_h}\int_T  {\rm Curl}~ ({\bf A^{-1}p_H}+u_H{\bf b^*}) (\beta_H-\beta_h)~ ds. \nonumber
\end{eqnarray}
Since $\rm {Curl}~ ({\bf A^{-1}p_H}+u_H{\bf b^*})=0$
\begin{eqnarray}
\norm{{\bf \tilde p_h -p_H}}^2& \lesssim & \sum_{E \in \mc E_h \setminus \mc E_H } \norm{[{\bf A^{-1}p_H}+u_H {\bf b^*}]\cdot \tau_E}_{E} \norm{\beta_h-\beta_H}_{E}.
\end{eqnarray}
The estimates of $\norm{\beta_h-\beta_H}_{E}$ from (\ref{drel22}) and the bound  $\norm{\nabla \beta_h}=\norm{{\rm Curl}~ \beta_h} \lesssim \norm{{\bf \tilde p_h -p_H}}$ together with (\ref{rea_2}) result in
\begin{equation}\label{rea_1}
 \norm{{\bf \tilde p_h -p_H}} \lesssim \eta_H (\mc E_H \setminus \mc E_h).
\end{equation}
%\begin{eqnarray}\label{rea_1}
% \norm{{\bf \tilde p_h -p_H}}^2%& \lesssim & \sum_{E \in \mc E_H \setminus \mc E_h } \norm{[{\bf A^{-1}p_H}+u_H{\bf b^*}]\cdot \tau_E}_{ } h_E^{\frac{1}{2}} ~|{\beta_h}|_{H^1(\omega_E)}\nonumber\\
% &\lesssim& \eta_H (\mc E_H \setminus \mc E_h)~\norm{{\bf \tilde p_h -p_H}}.\label{rea_3}
% \end{eqnarray}
% Hence, %a division by $\norm{{\bf \tilde p_h -p_H}}$ shows
 %\begin{eqnarray}
 %\end{eqnarray}
 A combination of (\ref{rea_1}), (\ref{rea_9}) and (\ref{rea_8}) leads to
 \begin{eqnarray}\label{53}
 \norm{  {\bf p_h-  p_H } } ^2\lesssim \eta_H^2(\mathcal E_H \setminus \mathcal E_h) + \norm{H (f_h- f_H)}^2 + \epsilon^2 ( e_h^2+  e_H^2) +\mu_h^2 +\mu_H^2.
 \end{eqnarray} 
 The combination of (\ref{rea_8}) and (\ref{53}) completes the proof.\hfill${\Box}$
 %\qquad \endproof %\hfill${\Box}$
 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% \begin{corollary}\label{coro2} For any $\epsilon >0$ and there exists constant $C_4 >0$ such that %Using quasi-discrete reliability result
% \begin{eqnarray*}
%\norm{{\bf A^{-1/2} (p_h -p_H)}}^2+\norm{u_h-u_H}^2 \leq C_4\big( \eta_H^2(\mathcal E_H \setminus \mathcal E_h) + \norm{H (f_h- f_H)}^2 +  \epsilon^2 ( e_h^2+  e_H^2) +\mu_h^2 +\mu_H^2\big).
%\end{eqnarray*}
% \end{corollary}
 %Proof- 
 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\section{Convergence Analysis}
This section is devoted to the convergence analysis of the adaptive mixed finite element method.
 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
\subsection{Contraction property}
Denote two consecutive adaptive loop levels as $\ell$ and $\ell +1$. Let $e_\ell$ , $\eta_\ell$ and $\mu_\ell$  denote
the error and the estimator terms on the level  $\ell$ with triangulation $\mathcal T_\ell$. 
 Based on the reduction properties of the error, the error estimators and quasi-orthogonal property developed in last section, the  contraction property is proved  for 
 the weighted term $\xi^2_\ell$, which is a linear combination of error $e^2_\ell$ and the  estimator terms $\eta^2_\ell$ and $\mu^2_\ell$, between two consecutive adaptive loops.  
 \begin{theorem}\label{thm1}(Contraction Property) Let $\mathcal T_{\ell+1}$ be a refinement of $\mathcal T_{\ell}$ using  AMFEM algorithm. Given $0<\theta_A, \theta_B <1$, there
exist positive parameters $\alpha,\beta, \kappa$ and $0<\rho <1$ depending on constants $\alpha_1,\alpha_2, \varLambda_1, \varLambda_2, \varLambda_3$ from 
Lemmas \ref{lm2}, \ref{lm3} and \ref{lm6} such that on any level $\ell\geq 0$,  the weighted term $\xi_\ell^2$ satisfies the following contraction property:
\begin{eqnarray}
 \xi_{\ell+1}^2 \leq \rho ~\xi_\ell^2, \label{contra}~~
 where ~~\xi_\ell^2:= \eta_\ell^2+\alpha e_{\ell}^2 +\beta \mu_\ell^2,\nonumber
\end{eqnarray}
whenever the initial mesh is chosen with small mesh-size. 
\end{theorem}\\
{\it Proof}.  %Consider  the first {\bf Case (A)} $\mu_\ell^2 \leq \kappa \eta_\ell^2$ of algorithm.
  For the {\bf Case (A)},  $\mu_\ell^2 \leq \kappa \eta_\ell^2$. The combination of (\ref{et_rd_11}) and (\ref{et_rd2})
 with a positive parameter $\beta$, to be chosen later, yields 
 % $0<\beta =2 \max\big{\{}C_5 \varLambda_3+C_5 \alpha_4 C_{rel}, \frac{C_5\varLambda_3+C_5\alpha_4C_{rel}+D}{1-\rho_B}\big{\}}, D$ will be specified later
\begin{eqnarray}\label{con_2}
\eta_{\ell+1}^2+\beta \mu_{\ell+1}^2 \leq  \rho_A \eta_\ell^2  + (1+\delta_1) \beta \mu_\ell^2 +(\varLambda_2+\beta \varLambda_1 h_\ell^2) E_\ell^2,
\end{eqnarray}
where $h_\ell$ denotes the mesh-size at  the level $\ell$ of the triangulation.
With a choice of the initial mesh-size $h_3\leq\min \{ h_1,h_2, \frac{1}{\sqrt{\beta}}\}$,  multiply (\ref{er_re})
with the constant  $C_5=({\varLambda_1+ \varLambda_2 })/{\alpha_3}$, and then add %to (\ref{er_re})
with (\ref{con_2}) to obtain  
%Let $C_4=\frac{\varLambda_2+ \varLambda_1 }{\varLambda_3}$ multiply to (\ref{er_re}) and for $H_3^2=H^2< \frac{1}{\beta}$ , gives
%\begin{eqnarray}
%C_4(1-\alpha_1)e_{\ell+1}^2\leq C_4 e_{\ell}^2-(\varLambda_2+ \varLambda_1) E_\ell^2+C_4 \varLambda_4 \mu_\ell^2.
%\end{eqnarray}
%Combining above two equation 
\begin{eqnarray*}
\eta_{\ell+1}^2+C_5(1-\alpha_1)e_{{\ell+1}}^2+\beta \mu_{\ell+1}^2 \leq  C_5 e_{{\ell}}^2+ \rho_A \eta_{\ell}^2  +( (1+\delta_1) \beta +C_5 \varLambda_3) \mu_\ell^2.
%+(\varLambda_2+\beta \varLambda_1 H^2) E_H^2.
 %-(\varLambda_2+ \varLambda_1) E_H^2+C_4 \varLambda_4 \mu_H^2.
\end{eqnarray*}
%As we define $\xi_h^2:= osc_h^2 (f)+\alpha \eta_h^2+\beta \mu_h^2$ ,
Define $\alpha:=C_5(1-\alpha_1)$; $0<\alpha_4:=\frac{1}{2}\min\{1, \frac{1-\rho_A}{C_5C_{rel}} \}$,  %such that $0<A_0<B_0<1$ 
 and use the reliability result (\ref{estimator1}) that is, $e_{\ell}^2\leq C_{rel}(\eta_\ell^2+\mu_\ell^2)$ to obtain
\begin{eqnarray}\label{con_4}
\xi_{\ell+1}^2&&\leq C_5(1-\alpha_4)e_{\ell}^2+C_5\alpha_4 C_{rel}(\eta_\ell^2+\mu_\ell^2) +\rho_A \eta_\ell^2  +( (1+\delta_1) \beta +C_5 \varLambda_3) \mu_\ell^2 \nonumber \\
&&\leq C_5(1-\alpha_4)e_{\ell}^2+(\rho_A+C_5\alpha_4 C_{rel})\eta_\ell^2 +( (1+\delta_1) \beta +C_5 \varLambda_3+C_5\alpha_4 C_{rel}) \mu_\ell^2.
\end{eqnarray}
Since the marking criteria  implies  $E  \kappa \eta_\ell^2-E \mu_\ell^2>0$, where $E=2(C_5 \varLambda_3+C_5\alpha_4 C_{rel})$, an addition of this term on the right-hand side of (\ref{con_4}) yields
\begin{eqnarray}
\xi_{\ell+1}^2&&\leq C_5(1-\alpha_4)e_{\ell}^2+(\rho_A+C_5\alpha_4C_{rel}+E \kappa  )\eta_\ell^2\nonumber \\
&&~~~ +((1+\delta_1)\beta-(C_5 \varLambda_3+C_5\alpha_4 C_{rel}) ) \mu_\ell^2.
\end{eqnarray}
A use of Remark \ref{rmk2} and the definition   of $\alpha_4$ with the choice of the parameters   
$$\kappa <\kappa_0:= \frac{({1-\rho_A-C_5\alpha_4C_{rel}})}{E}, ~~\beta\geq 2(C_5 \varLambda_4+C_5\alpha_4C_{rel}), %
 ~~{\rm and}~~~ \delta_1<\min \Big \{ \frac{E}{2\beta}, \frac{\theta_B}{4-\theta_B} \Big { \} },$$
   yield on any level $\T_{\ell+1}$, %with {\bf Case (A)},
 a contraction   $\xi_{\ell+1}^2 \leq  \rho_1 \xi_\ell^2 $, where % with $\rho_1$ defined by
\begin{eqnarray*}
0< \rho_1=\max \Big{\{ } \frac{C_5(1-\alpha_4)}{\alpha}, \rho_A+C_5\alpha_4C_{rel}+E \kappa , \frac{((1+\delta_1)\beta-E/2)}{\beta} \Big \}<1.
\end{eqnarray*}
%\begin{eqnarray*}
%0< \rho_1=\max \big\{  \frac{C_4(1-\alpha_4)}{\alpha}, \rho_A+C_4\alpha_4C_{rel}+2\beta \kappa , \frac{C_4 \varLambda_4+C_4\alpha_4C_{rel}}{\beta} \big \}<1.
%\end{eqnarray*}
%Selecting   and $C_0=2\beta$ shows that the third term less than one.\\
%The chioce of $B_0=\min\{1, \frac{1-\rho_A}{C_4C_{rel}} \}$ and $K=\frac{1-\rho_A-C_4B_0C_{rel}}{C_0	}$ implies second term less than one.
%And As $A_0=\delta_3+\frac{2C^2\epsilon^2}{\delta_4}$ finally selecting $\delta_3<B_0$ and sufficiently small $H_3$ such that
 %$\epsilon^2 <\frac{(B_0-\delta_3)\delta_4}{2C^2}$ implies $A_0<B_0$, hence
%first term is less than one. \\
%Let us now consider  {\bf Case (B)}  $\mu_\ell^2 > \kappa \eta_\ell^2$  of the algorithm.\\
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
For the {\bf Case (B)}: $\mu_\ell^2 > \kappa \eta_\ell^2$. 
Similar to proof of the {\bf Case (A)},  the  equation  corresponding to (\ref{con_4}) for the {\bf Case (B)} is 
%\begin{eqnarray}
%\eta_h^2+\beta \mu_h^2 \leq  (1+ \delta_B) \eta_H^2  + \rho_B \beta \mu_H^2 +(\varLambda_2+\beta \varLambda_1 H^2) E_H^2.
%\end{eqnarray}
%Let $C_4=\frac{(\varLambda_2+ \varLambda_1 }{\varLambda_3})* $ (\ref{er_re}) and for $H^2< \frac{1}{\beta}$ , gives
%\begin{eqnarray}
%C_4(1-A_0)({\bf e_{p_h}}^2+e_{u_h}^2)\leq C_4({\bf e_{p_H}}^2+e_{u_H}^2)-(\varLambda_2+ \varLambda_1) E_H^2+C_4 \varLambda_4 \mu_H^2.
%\end{eqnarray}
%Combining above two equation 
%\begin{eqnarray}
%\eta_h^2+C_4(1-A_0)({\bf e_{p_h}}^2+e_{u_h}^2)+\beta \mu_h^2 \leq  C_4({\bf e_{p_H}}^2+e_{u_H}^2)+ (1+\delta_2) \eta_H^2  +(\rho_B \beta   +C_4 \varLambda_4) \mu_H^2.
%\end{eqnarray}
%As we define $\xi_h^2:= osc_h^2 (f)+\alpha \eta_h^2+\beta \mu_h^2$ , choose $B_0$  such that $0<A_0<B_0<1$ with reliability result
%${\bf e_{p_H}}^2+e_{u_H}^2\leq C_{rel}(\eta_H^2+\mu_H^2)$ result in
\begin{eqnarray}\label{gg}
\xi_{\ell+1}^2&&\leq C_5(1-\alpha_4)e_{\ell}^2+C_5\alpha_4C_{rel}(\eta_\ell^2+\mu_\ell^2) +(1+\delta_2) \eta_\ell^2  +(\rho_B \beta +C_5 \varLambda_3) \mu_\ell^2\nonumber \\
&&\leq C_5(1-\alpha_4)e_{\ell}^2+(1+\delta_2+C_5\alpha_4C_{rel})\eta_\ell^2   +(\rho_B \beta +C_5 \varLambda_3+C_5\alpha_4C_{rel}) \mu_\ell^2.
\end{eqnarray}
The  marking  criteria  implies $  D \mu_\ell^2-D \kappa \eta_\ell^2>0$, where $0<D=\frac{1.5+C_5\alpha_4C_{rel}}{\kappa}$.
Add this term on the right-hand side of (\ref{gg}) to obtain 
\begin{eqnarray*}
\xi_{\ell+1}^2\leq C_5(1-\alpha_4)e_{\ell}^2+(2+C_5\alpha_4C_{rel}-D \kappa )\eta_\ell^2
+(\rho_B\beta  +C_5 \varLambda_3+C_5\alpha_4 C_{rel}+D ) \mu_\ell^2.
\end{eqnarray*}
A use of Remark \ref{rmk2}, the parameters choice and  $\beta >\frac{C_5 \varLambda_4+C_5\alpha_4C_{rel}+D}{1-\rho_B}$ %and $ \kappa <\kappa_0:= \frac{1-\rho_A-C_4\alpha_4C_{rel}}{2 \beta_1}$,
 yields that, on any level $\T_\ell$  the contraction   $\xi_{\ell+1}^2 \leq  \rho_2 \xi_\ell^2 $  holds true, where $\rho_2$ is defined by %with {\bf Case (B)},
%Proof of $\xi_h^2 \leq  \rho \xi_H^2 $ will complete by showing
\begin{eqnarray*}
0< \rho_2=\max \Big\{  \frac{C_5(1-\alpha_4)}{\alpha}, 2+C_5\alpha_4C_{rel}-D\kappa, \frac{\rho_B\beta+C_5 \varLambda_3+C_5\alpha_4C_{rel}+D }{\beta} \Big \}<1.
\end{eqnarray*}
%Selection of $ D_0=\frac{0.5+\delta+C_4B_0C_{rel}}{K}<\frac{1}{2}$ and $\beta_2 =2\frac{C_4 \varLambda_4+C_4\alpha_4C_{rel}+D}{1-\rho_B}$ and remaining as in case 
%A show $0<\rho_2<1$.\\
 Finally, the combination of  both cases with %$\beta:= \max \{2(C_4 \varLambda_4+C_4B_0C_{rel}) , 2\frac{C_4 \varLambda_4+C_4B_0C_{rel}+D_0 }{1-\rho_B} \}$
 %$\beta=\max\{\beta_1, \beta_2 \}$ and
    $0<\beta =2 \max\big{\{}C_5 \varLambda_3+C_5 \alpha_4 C_{rel},$ $ \frac{C_5\varLambda_3+C_5\alpha_4C_{rel}+D}{1-\rho_B}\big{\}},   $  $\rho=\max\{\rho_1,\rho_2\}$
   % will be specified later
     and the initial mesh-size $h_3$ implies that (\ref{contra}) holds. \hfill${\Box}$ %\qquad \endproof %\hfill${\Box}$
%~~~~~~~~~~~~~~~~~~~
\begin{theorem}
 (Convergence) Under the assumptions of Theorem \ref{thm1}, there exist a constant $\rho \in (0, 1)$ and $C_0>0$ 
   depending only on the given data and the initial triangulation such that %and $\alpha, \beta, C>0 $ 
\begin{equation*}
 \eta_\ell^2 +\alpha e_{\ell}^2 +\beta \mu_\ell^2\leq C_0 \rho^\ell.
\end{equation*}
\end{theorem}
{\it Proof}.  The proof is a consequence of the contraction property in Theorem \ref{thm1}.\hfill${\Box}$%\qquad \endproof %\hfill${\Box}$\\
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\begin{remark}\label{rmk4} Using the reliability  result (\ref{estimator1}) and  the relation of  the error estimator $\eta_\ell^2$ and $\mu_\ell^2$ in the marking stratergy, the weighted term $\xi_\ell$, the  error  term $e_\ell$ and the error estimator terms $\eta_\ell$ and $\mu_\ell$ over triangulation $\mc T_\ell$ are obtained. 
Now $\xi_\ell$  satisfies
%\begin{eqnarray}
\[ \xi^2_\ell \leq \left\{
  \begin{array}{l l}
   C_6 \eta^2_\ell & \quad \text{for~~the~~  {\bf Case (A)}},\\
    \frac{C_6}{\kappa} \mu^2_\ell & \quad \text{for ~~the~~{\bf Case (B)}},
  \end{array} \right.\]
%\end{eqnarray}
where $C_6=1+\alpha C_{rel}+ (\alpha C_{rel}+\beta)\kappa$. % it can be proved .  
The reliability (\ref{estimator1}), Lemma \ref{lm8} with the  efficiency result (\ref{edeff}) implies
\begin{equation*}
\xi_\ell ^2\approx e_\ell^2 +\norm{h_{\mathcal T} f}^2.\end{equation*}
%\begin{\}
%\end{} 
\end{remark}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\begin{lemma}\label{lm14} Let $\T_h$ and $\T_{H}$ be two nested triangulations. 
Then for small mesh-size $h_4>0$, it holds for  $0<h\leq h_4$
%there exists the positive constant $C_C$ such that
\begin{eqnarray}\label{1}
\xi_h^2\lesssim \xi_{H}^2. 
\end{eqnarray}
 %$$e_{\ell+\e_1}^2 \leq C_6 \tau_1^2 \eta_\ell^2  $$
\end{lemma}
{\it Proof}.  %$\alpha e_\T^2 \leq \xi_\ell^2$\\ The definition of $\xi$ with
 From Lemma \ref{lm6}, Remark \ref{rmk3} and the efficiency result (\ref{lmeff}) implies
\begin{eqnarray}
\xi_h^2&&=\eta^2_h+\alpha e^2_h+\beta\mu_h^2\leq (C^{-1}_{\rm eff}+\alpha)e^2_{h} +2\beta \mu^2_{H} + \beta \varLambda_1 H ^2 E^2_{H}  \nonumber\\
&&\leq \frac{(C^{-1}_{\rm eff}+\alpha)}{1-\alpha_1}(e^2_{H}-(\frac{1}{2}-8C_3\varLambda_1 H^2) E^2_{H}+\varLambda_3 \mu^2_H )+2\beta\mu^2_{H} +\beta\varLambda_1 H^2 E_{H}^2.\label{2}%\nonumber\\
%&&\leq \frac{(C^{-1}_{\rm eff}+\alpha)}{1-\alpha_1}e^2_{\e_1}-\big({(C^{-1}_{\rm eff}+\alpha)}(\frac{1}{2}-8C_3\varLambda_1 h^2_{\e_1})  -\beta \varLambda_1 h_{\e_1}^2\Big)E^2_{\e_1}+\frac{(C^{-1}_{\rm eff}+\alpha)\varLambda_3}{1-\alpha_1} \mu^2_{\e_1} )\label{2}%+\beta(2\mu^2_{\e_1} E_{\e_1}^2 )
\end{eqnarray}
Select the initial mesh-size $h_4>0$ such that the coefficient  of $E_H^2$ is non-positive for $0<h\leq h_4.$
 Then (\ref{2}) implies (\ref{1}). Note that the inequality constant in (\ref{1}) is independent of AMFEM marking parameters. This concludes the proof. \hfill${\Box}$%\qquad \endproof %\hfill{$\Box$}
%~~~~~~~~~~~~~~~~%%%%%%%%%~~~~~~~~~~~~~~~~%%%%%%%%%%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~~~~~~###############~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~%%%%%%%%~~~~~~~~~~~~~~~~~%%%%%%%%%%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\subsection{Quasi-optimality}
In this subsection, the quasi-optimal convergence \cite{steven} of the adaptive algorithm MFEM is discussed with the help of the quasi-discrete reliability and the contraction property. 
\begin{definition} \cite{Car_apo, steven}(Approximation class) Given an initial triangulation $\mathcal T_0$ of $\Omega$ and $s>0$, the approximation class is defined as
\begin{eqnarray*}
\mc A_s :=\{ ({\bf p},u,f)\in H({\rm div},~ \Omega) \times L^2(\Omega) \times L^2(\Omega)~|~ \norm{({\bf p},u,f)}_{\mc A_s} <\infty \} ~~~{\rm with}
\end{eqnarray*}
\begin{eqnarray*}
 && \norm{({\bf p},u,f)}_{\mc A_s} := \sup_{N\in \mathbb N}\Big( N^s \inf_{|\mc T|-|\mc T_0|\leq N} (e^2(\mc T)+\norm{h_{\mathcal T} f}^2)^{{1}/{2}}\Big) {\rm and} ~ e^2(\T)=e^2_p(\T)+e^2_u(\T).
\end{eqnarray*}
Here, the infimum is  over all regular and NVB-generated refinements $\mc T$ of $\mc T_0$% called admissible triangulations,
 ~with the number of element domains $|\mc T|\leq N+|\mc T_0|$ and the exact error $e(\mc T)$ in MFEM solution.
\end{definition}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\begin{remark} \cite{steven,Car_apo} Characterisation of the  approximation class: $({\bf p}, u, f)\in \mathcal A_s$ implies for all $\epsilon_1 >0$, 
%there exists an admissible triangulation $\mc T_{\epsilon_1}$ suct that the associated MFEM solution $({\bf p_{\epsilon_1}},u_{\epsilon_1})$
%in $RT_0(\mc T_{\epsilon_1}) \times P_0 (\mc T_{\epsilon_1})$  satisfies
%\begin{eqnarray*}
%e^2_{\epsilon_1}+\norm{h f}\leq \epsilon_1^2~~~{\rm and}
% ~~~~~~|\mathcal T_\ell|-|\mathcal T_0| \lesssim \epsilon_1^{-1/s}\norm{({\bf p},u,f)}_{\mc A_s}^{1/s}.
%\end{eqnarray*}
%\end{remark}
%
An adaptive mixed finite element method is quasi-optimal convergent in the sense that 
given $({\bf p}, u, f)\in \mathcal A_s,~ {\rm and} ~\J \in \mathbb N$, the  AMFEM  algorithm  generates a triangulation $\mathcal T_\J$ with discrete  solution $({\bf p_\J},u_\J)\in RT_0(\mathcal T _\J)\times P_0(\mathcal T_\J)$ such that 
$$|\mathcal T_\J|-|\mathcal T_0| \leq \xi_\J ^{-\frac{1}{s}}\approx (e_\J^2 +\norm{h_{\mathcal T}f}^2)^{-\frac{1}{2s}}.$$
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~``
%~~~~~~~~~~~Quasi-optimality ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\begin{theorem}(Quasi-optimality) Assume $({\bf p}, u, f)\in \mathcal A_s.$ Let $\{\mc T_\J \}_{\J\geq 0}$ be the sequence of the meshes generated  by AMFEM algorithm and $\{ ({\bf p_\J},u_\J)\in {RT_0(\mc T_\J)}\times P_0(\mc T_\J)\}_{\J\geq 0}$ be a corresponding sequence of approximate solutions.
Then, for   small initial mesh-size $h_0$, the following estimates holds true for $0<h\leq h_0$
\begin{equation}\label{obtim}
|\mathcal T_\J|-|\mathcal T_0| \lesssim \xi_\J^{-\frac{1}{s}}. %\qquad \text{forall}~~~ \ell\geq 0.
\end{equation}
\end{theorem}
{\em Proof}.
Consider
 \begin{eqnarray}
\qquad \qquad |\mathcal T_\J|-|\mc T_0| \leq \sum_{\ell=0}^{\J-1} \Big(|\mc T_{\ell+1}|-|\mc T_\ell| \Big) 
\end{eqnarray}
%Now by the combining both the cases from (\ref{opti_3}) and (\ref{opti_5}),
Now $|\T_{\ell+1}|-|\mc T_\ell|\lesssim |\mc M_\ell|$ for $\ell\geq 0$ \cite{biven,steven},%\cite{%}
~ where $\T_{\ell+1}$ is a refinement of $\T_\ell$ and $\mc M_\ell$ denotes the set of the marked edges or elements in the triangulation level $\T_\ell$. Then
\begin{eqnarray}\label{qo1}
|\mathcal T_\J|-|\mc T_0| %&\leq& \sum_{j=0}^{\ell-1} |\mc T_{j+1}|-|\mc T_j|\nonumber \\
&\lesssim& \sum_{\ell=0}^{\J-1} |\mc M_\ell|. %\qquad (\mc T_{j+1}-{\mc T_j}\leq C_{11} {\mc M_j})\nonumber\\
%\lesssim \sum_{j=0}^{\ell-1} \xi_j^{-\frac{1}{s}}.\nonumber %\\
%&\lesssim & \xi^{-\frac{1}{s}_0+\xi^{-\frac{1}{s}_1+\cdots+\xi^{-\frac{1}{s}_{\ell-1}\nonuber\\ (contraction property)
\end{eqnarray}
To estimate of $|\mc M_\ell|$, use the characterization of the approximate class and the overlay.\\
If $({\bf p}, u, f)\in \mathcal A_s$, then for $\e_1:=\tau_1 \xi_\ell$, %=\min \big{\{} , \big{\}}$,  
there exists some admissible triangulation $\mc T_{\e_1}$ obtained as refinement of the initial triangulation $\T_0$  such that %\cite{steven,car_cgt1}
\begin{eqnarray}\label{opt1}
\xi^2_{\e_1}=\eta^2_{\e_1}+\alpha e^2_{\e_1}+\beta \mu^2_{\e_1}\leq \e_1^2 ~~~\rm{and}~~~ |\T_{\e_1}|-|\T_0|\lesssim \e_1^{-\frac{1}{s}}.
\end{eqnarray}
Here, $0<\tau_1=\min\{ \tau_2,\tau_3 \}$, and $\tau_2 ~{\rm and} ~\tau_3$ will be specified later.\\
Let $\mc T_{\ell+{\e_1}}:=\mc T_{\e_1}\oplus \T_\ell $ be the overlay of $\T_{\e_1}$ and $\T_\ell$. As in  \cite{quasi}, the  number of elements of the overlay $\mc T_{\ell+{\e_1}}$  can be bounded by%and (\ref{opt1})
\begin{eqnarray}\label{opt2}
|\mc T_{\ell+{\e_1}}|-|\T_\ell|\leq|\T_{\e_1}|-|\T_0|.%\lesssim\e_1^{-\frac{1}{s}}.
\end{eqnarray}
To estimate  $|\mc M_\ell|$ with the first case of Mark algorithm, that is, $\kappa \eta^2_\ell \geq \mu^2_\ell$,
first define $\mc M^\dagger :=\mc E_\ell\setminus \mc E_{\ell+\e_1} \subset \mc E_\ell$ as the set of edges of $\T_\ell$ being refined in $\T_{\ell+\e_1}$. Note that, if  $\mc M^\dagger$ satisfies the marking criteria of the {\bf Case (A)}, that is, 
\begin{eqnarray}\label{opti_2}
\theta_A \eta_\ell^2 \leq \eta^2_\ell(\mc M^\dagger),
\end{eqnarray} 
then,  $|\mc M_\ell|\leq|\mc M^\dagger|$, where $|\mc M_\ell|$ be the set of marked edges at level $\ell$.\\
The quasi-discrete reliability result Lemma \ref{coro2} over the triangulations $\T_\ell$ and $\T_{\ell+\e_1}$  with  $\epsilon^2=1/4C_4 $ and the initial mesh-size $h_5$ shows
\begin{eqnarray*}
C_4 \eta_\ell^2 (\mc M^\dagger)\geq  E_\ell^2 -C_4 osc^2_\ell (f)-\frac{1}{4}(e^2_\ell+e^2_{{\ell+\e_1}})-C_4(\mu^2_\ell+\mu^2_{\ell+\e_1}).
\end{eqnarray*}
Remark \ref{rmk3} yields $\mu^2_{{\ell+\e_1}} \leq 2 \mu_\ell^2+\Lambda_1 h_\ell^2 E^2_\ell$, where $h_\ell$ denotes the mesh-size over the  triangulation $\T_\ell$ and note that $ osc^2_\ell (f)\leq \mu_\ell^2$. Hence,  %and $osc^2 f \leq \mu_\ell^2$ implies 
\begin{eqnarray}\label{qo3}
C_4 \eta_\ell^2 (\mc M^\dagger)\geq(1-C_4\varLambda_1 h_\ell^2)  E_\ell^2 -4 C_4 \mu_\ell^2 -\frac{1}{4}(e^2_\ell+e^2_{{\ell+\e_1}}).% -C_4\varLambda_1^2 h_\ell^2E_\ell^2.
\end{eqnarray}
For  small initial mesh-size $h_4$, Lemma \ref{lm14}, (\ref{opt1}), the choice of $\epsilon_1$ and Remark \ref{rmk4} imply
 $$e_{{\ell+\e_1}}^2 \leq \xi_{\ell+\e_1}^2\lesssim  \xi^2_{\e_1}\lesssim  \tau_1^2 \xi^2_\ell \lesssim \tau_1^2  \eta_\ell^2 \leq  \tau_2^2 \eta_\ell^2,$$ 
 that is, $e_{{\ell+\e_1}}^2 \leq C_7 \tau_2^2 \eta_\ell^2$.  Thus, Corollary \ref{coro3} supplies the lower bound for $E_\ell$. With these results, (\ref{qo3}) leads to
\begin{eqnarray}
C_4 \eta_\ell^2 (\mc M^\dagger)&\geq &\frac{1}{2} (\frac{1}{2}-C_4\varLambda_1^2 h_\ell^2) e_\ell^2-\frac{5}{4}C_7 \tau^2 \eta_\ell^2-\big( 4C_4+ \frac{\varLambda_3}{2}\big)\mu_\ell^2.%-4 C_1 \mu_\ell^2%-\frac{1}{4}(e^2_\ell+e_\T^2)
\end{eqnarray}
%
Note that some positive terms are neglected from the right-hand side. 
Choose the initial mesh-size $h_6 :=\min\{h_2,h_5,\frac{1}{(4C_4\varLambda_1)^{1/2}} \}$ and  the marking parameter $\theta_A=\min \{1 , \frac{c_{\rm eff}}{9C_4}\}$.  % and $e_\T^2 \leq C_6 \tau^2  \eta_\ell^2$ and 
 The {\bf Case (A)} relation $\kappa \eta_\ell^2\geq \mu_\ell^2$ along-with  the efficiency result (\ref{edeff}) leads to
%Using bounds $e_\T^2$ and $C_{\rm eff_1} \eta^2\leq e^2_\ell$ and for $\epsilon_2<\frac{1}{C_1}$ there exist $H_4$ such that
\begin{eqnarray}
C_4 \eta_\ell^2 (\mc M^\dagger)&\geq &\Big(\frac{ C_{\rm eff}}{8} -\frac{5}{4} C_7\tau_2^2 - \kappa\big( 4C_4+ \frac{\varLambda_3}{2}\big)-C_4 \theta_A\Big)\eta_\ell^2+ C_4 \theta_A \eta_\ell^2. \nonumber
\end{eqnarray}
%\begin{eqnarray}
%C_1 \eta_\ell^2 (\mc E_\ell \setminus \mc E_{\T})&\geq & [ C_{\rm eff}^{-1}\frac{1}{2}(\frac{1}{2}-C_1 \Lambda_1 \tau^2)-(\frac{\Lambda^4}{4}+4C_1)+ 
%C_1\theta_A] \eta_\ell^2 +C_1 \theta_A \eta_\ell^2 .
%\end{eqnarray}
The  selection of $ \kappa < \min \{ \kappa_0,  \frac{{C_{\rm eff}}/{8}-C_4\theta_A }{4C_4+ {\varLambda_3}/{2}}\big\}$, and $\tau_2^2:=\big( \frac{(c_{\rm eff}/8)-\kappa ( 4C_4+ {\varLambda_3}/{2})-C_4 \theta_A}{(5/4)C_7} \big) $ imply that (\ref{opti_2}) holds.\\
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Since  $\mc M_\ell$ is chosen to be  the minimal cardinality set satisfying (\ref{opti_2}),   Lemma 4.4 of \cite{Car_apo},  %\cite[Lemma 4.4]{Car_apo},  
(\ref{opt2}) and (\ref{opt1}) altogether imply
\begin{eqnarray}\label{opti_3}
|\mc M_\ell|\leq |\mc M^\dagger_\ell| \lesssim |\T_{\ell+\e_1}| -|\mc T_\ell |\leq|\T_{\e_1}|-|\T_0|\lesssim\e_1^{-\frac{1}{s}}\lesssim\xi^{-\frac{1}{s}}_\ell.
\end{eqnarray}
Consider the {\bf  Case (B)} of Mark algorithm, that is, $\kappa \eta^2_\ell \leq \mu^2_\ell $.  Let $\mc M^* = \T_\ell\setminus \T_{\ell+\e_1}\subset \T_\ell$ be the set of elements of $\T_\ell$ refined in $\T_{\ell+\e_1}$. The proof of $\mc M^*$ satisfies the marking criteria of the {\bf Case (B)}, that is,
\begin{eqnarray}\label{opti_4}
\theta_B \eta^2_\ell\leq \eta^2_\ell(\mc M^*),
\end{eqnarray} 
and this will imply $|\mc M_\ell|\leq|\mc M^*|$, where $\mc M_\ell$ is the set of the marked 
elements at the level $\ell$.\\
Corollary \ref{lm13} and   the quasi-discrete reliability result in Lemma \ref{coro2}, with  $\epsilon=1 $ and initial mesh-size $h_7$ over the nested triangulations $\T_{\ell+\e_1}$ and $\T_\ell$ imply %. In Lemma \ref{coro2}, for $\epsilon=1 $, there exists sufficiently small initial mesh-size $h_7>0$ such that
%The lemma \ref{lm13}, and the discrete reliability (\ref{}) with $\epsilon =1$ sufficiently small initial mesh-size $h_2$, leads to
\begin{align*}
\mu^2_{{\ell+\e_1}} &\geq \frac{\gamma_0^2}{2} \mu^2_\ell -2\norm{h_\ell(f_\ell -f_{{\ell+\e_1}})}^2-C_{1} h_\ell^2 E_\ell^2 \nonumber\\
 &\geq \big(\frac{\gamma_0^2}{2}-C_4C_{1} h_\ell^2\big) \mu^2_\ell -(2+C_4C_{1}h_\ell^2)\norm{h_\ell(f_\ell -f_{{\ell+\e_1}})}^2\nonumber\\
&~~~ -C_{1}C_4 h_\ell^2\big(\eta^2_\ell +{}(e_\ell^2+e^2_{{\ell+\e_1}})
+\mu^2_{\ell+\e_1}\big).
\end{align*}
The reliability result, the  relation $\kappa \eta_\ell^2 <\mu_\ell^2$ for the  {\bf Case (B)} and the rearrangement  of  terms   imply %(for simplicity,  by $2\beta$ LHS and $\beta$ to RHS as $\beta>0$
\begin{eqnarray}\label{op2}
~~\mu^2_{{\ell+\e_1}}\geq\frac{(\frac{\gamma_0^2}{2}-C_{8} h_\ell^2 )}{(2+C_4C_{1} h_\ell^2)}\mu^2_\ell -\norm{h_\ell(f_\ell -f_{{\ell+\e_1}})}^2- {C_{1}C_4 h_\ell^2}e^2_{{\ell+\e_1}},
\end{eqnarray}
where $C_{8}:= C_4C_{1}(Crel+1)(1+1/\kappa)$.
 Lemma \ref{lm14}, the choice of $\epsilon_1$, that is, $\epsilon_1=\tau_1\xi_\ell$ and   the  marking criteria    $\kappa \eta_\ell^2< \mu_\ell^2$ in Remark {\ref{rmk4}} result in
\begin{eqnarray}\label{op1}
\alpha e_{{\ell+\e_1}}^2+\beta \mu^2_{{\ell+\e_1}}\leq\xi^2_{{\ell+\e_1}} \lesssim \xi_{\e_1}^2\lesssim \tau_1^2 \xi^2_\ell \lesssim \tau_3^2 \mu_\ell^2,
\end{eqnarray}
%
that is, $\alpha e_{{\ell+\e_1}}^2+\beta \mu^2_{{\ell+\e_1}}\leq C_9\tau_3^2 \mu_\ell^2 $. The combination of (\ref{op2})-(\ref{op1}), for  small  initial mesh-size 
$h_8=\min\big{\{ }h_7, h_4, ({\alpha}/{\beta C_{1}C_4})^{1/2},({\gamma_0^2}/{4C_{8})})^{1/2}\big{\}}$ and some simplifications show
%substitute in (\ref{mu1}) and for $h_3$ such that $\alpha-\beta\frac{C_{10} h_\ell^2}{4} >0$, implies
\begin{equation*}
\beta{\mu^2(\mc M^*)}\geq \beta\norm{h_\ell(f_\ell -f_{{\ell+\e_1}})}^2\geq \Big[\frac{\beta{\gamma_0^2}}{4(2+C_4C_{1} h_\ell^2)}-2C_9 \tau_3^2  % - C_9 \tau ^2(\frac{1+\alpha C_{rel}}{\kappa}+\alpha  C_{rel}
- \beta \theta_B\Big] \mu^2_\ell+ \beta \theta_B \mu^2_\ell.
\end{equation*}
%&&\beta(1+) \mu_\ell^2(\T_\ell \setminus \T) \geq \beta \norm{H(f_\ell -f_\T)}^2 &\geq &\big(\frac{}{} -\tau_1^2-\beta %\theta_B  \big)  \mu^2_\ell+\beta \theta_B \mu^2_\ell\\
% &\geq & \beta \theta_B \mu^2_\ell
%\end{eqnarray}
The selections $ \tau^2_3: =\frac{\beta{\gamma_0^2}}{16(2+C_4C_{1} h_\ell^2)C_9}$ and $ \theta_B =\frac{{\gamma_0^2}}{8(2+C_4C_{1} h_\ell^2)}$ lead to (\ref{opti_4}).
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Since $\mc M_\ell$ is chosen to be the minimal cardinality set satisfying (\ref{opti_4}),  (\ref{opt2}) and (\ref{opt1}) yield
\begin{eqnarray}\label{opti_5}
|\mc M_\ell|\leq|\mc M^*|\leq |\T_{\ell+\e_1}|-|\T_\ell|\leq|\T_{\e_1}|-|\T_0|\lesssim\e_1^{-\frac{1}{s}}\lesssim\xi^{-\frac{1}{s}}_\ell.
\end{eqnarray}
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~!!!!!!!!!!!!!!!!!~~~~~~~~~~~~~~~~~````
Now a combination of both the cases, that is, (\ref{opti_3}) and (\ref{opti_5}) with (\ref{qo1}) leads to % $|\T_{j+1}|-|\mc T_j|\lesssim |\mc M_\ell|$ for $j\geq 0$ \cite{steven}%\cite{%}
\begin{eqnarray}
|\mathcal T_\J|-|\mc T_0| %&\leq& \sum_{j=0}^{\J-1} |\mc T_{j+1}|-|\mc T_j|\nonumber \\
&\lesssim& \sum_{\ell=0}^{\J-1} |\mc M_\ell| %\qquad (\mc T_{j+1}-{\mc T_j}\leq C_{11} {\mc M_j})\nonumber\\
\lesssim \sum_{\ell=0}^{\J-1} \xi_\ell^{-\frac{1}{s}}.\nonumber %\\
%&\lesssim & \xi^{-\frac{1}{s}_0+\xi^{-\frac{1}{s}_1+\cdots+\xi^{-\frac{1}{s}_{\ell-1}\nonuber\\ (contraction property)
\end{eqnarray}
A use of  the contraction property in Theorem \ref{thm1} shows
\begin{align}
|\mathcal T_\J|-|\mc T_0| %& \lesssim & \xi^{-\frac{1}{s}}_0+\xi^{-\frac{1}{s}}_1+ \cdots +\xi^{-\frac{1}{s}}_{\ell-1}(contraction ~property) \nonumber\\ 
&\leq ( \rho^{\J(\frac{1}{2s})}+ \rho^{(\J-1)(\frac{1}{2s})}+\cdots+\rho^{\frac{1}{2s}} \big) \xi^{-\frac{1}{s}}_{\J}\nonumber\\
 &\leq  \xi^{-\frac{1}{s}}_{\J} \sum_{j=1}^{\J} \rho^{\frac{j}{2s}}
%\leq  \xi^{-\frac{1}{s}_{\ell}}\Big[ \frac{1-\rho^\frac{\ell}{2s}}{1-\rho^\frac{1}{2s}}\Big]\nonumber\\
\lesssim \xi_{\J}^{-\frac{1}{s}}\Big[ \frac{1}{1-\rho^\frac{1}{2s}}\Big] \lesssim  \xi^{-\frac{1}{s}}_{\J},
\end{align}
 and this concludes the proof. \hfill${\Box}$%\qquad \endproof %\hfill${\Box}$
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
\section{Numerical Experiments}
This section shows the performance of the adaptive algorithm (AMFEM) on some benchmark problems. 
\begin{example} \label{example1} For the problem (\ref{eq1}), set $\Omega=(0,1)\times (0, 1)$ and the coefficients ${\bf A}=I,{\bf b}=(1,1),\gamma=2$. Choose the right-hand side $f$ and $u|_{\partial \Omega}$  such that $u=exp(-100\norm{x-x_0}^2)$.
%and  in the equation (\ref{eq1}).
%
%Consider the domain $\Omega=(0,1)\times (0, 1)$, with the right-hand side $f$ and $u|_{\partial \Omega}$  such that $u=exp(-100\norm{x-x_0}^2)$
%and the coefficients ${\bf A}=I,{\bf b}=(1,1),\gamma=2$ in the equation (\ref{eq1}).
\end{example}\\
 The initial uniform criss-cross triangulation  $\T_0$ has mesh-size $h=0.25$ as shown in Fig. 5.1(a). The adaptive algorithm  is performed  with parameters  $\theta_A=0.5, \theta_B=0.5$ and  $\kappa=0.8$.
\begin{figure}[!ht]%\nonumber
\centering
\subfloat[] 
{\includegraphics[width=6.0cm]{initia.png}} \subfloat[ ]{\includegraphics[width=6.0cm]{Exam_2_refinement.png}}
 %\subfloat[ ] {\includegraphics[width=10.0cm]{ex_2_plots.png}}
\caption{(a) The initial triangulation (b) The adaptive refined mesh at the level 10.}%.(c) The convergence rates of errors and the estimators with respective Ndof  }%$C_{rel}$ for different $\gamma$ with adaptive refinement}
\end{figure}
\begin{table}[h!]
%\tiny
\footnotesize
\begin{center}
  \begin{tabular}{ c | c | c  |c|c |c}%\footnotesize%\label{tb1}
    \hline
   Ndof &$\norm{u-u_\ell}$ &$\norm{{p-p}_\ell}$&$\eta_\ell$&$\mu_\ell$& ~Marking Case\\
 \hline \hline
 168    &  0.0571&    1.3649 &   1.5499 &   4.8186  & B\\
 203   &   0.0435 &   0.9670  &  2.7266  &  4.1602  &  B\\
 233    &  0.0336  &  0.8209   & 2.9629   & 3.3392   & B\\
368   & 0.0254&    0.6160&    2.2798&    2.3559&   B\\
523    &0.0197&    0.4471 &   1.5828 &   1.6907 &  B\\ \hline
948   & 0.0152 &   0.3287  &  1.1607  &  1.1389  &  B\\
1783   & 0.0103 &   0.2451  &  0.8943  &  0.8198  &  B\\
3373    &0.0074  &  0.1647   & 0.6013   & 0.5797   & B\\
 6138   &0.0058   & 0.1334&    0.4908&    0.4182&    A\\
6918    &0.0048    &0.1027 &   0.3692 &   0.3997 &   B\\ \hline
12053   & 0.0041&    0.0866 &   0.3120 &   0.2941 &   B\\
21298  &  0.0033 &   0.0698  &  0.2527  &  0.2177  &  A\\
 23198&   0.0026  &  0.0547   & 0.1991   & 0.2137   & B\\
 40658   &0.0022   & 0.0461&    0.1685&    0.1599&    B\\
  %70453  &0.0018    &0.0378 &   0.1378 &   0.1198 &   A\\
 %78663  & 0.0014    &0.0296  &  0.1078  &  0.1155  &  B\\
  \hline
  \end{tabular}
\end{center}
    \caption{Numerical results of AMFEM for Example \ref{example1}}
\end{table}
\begin{figure}[!ht]%\nonumber
\centering
{\includegraphics[width=6.2cm]{ex_2_solution.png}}
{\includegraphics[width=6.2cm]{cht_plot3.pdf}}

%\begin{tikzpicture}[scale=.65]
%\begin{loglogaxis}[
%%title=Convergence Plot,
%xlabel={Number of degrees of freedom},
%ylabel={ Errors and estimators},%}
%font=\sffamily\scriptsize,
%]
%%\addplot table {ex1/datH1.dat};
% %\addlegendentry{ $\|\! {\bf p} {-} {\bf p}_h\!\|$ uniform}
%%\addplot table {ex1/datL2.dat};
%\addplot table {ex1/datH1L2.dat};
%%\addplot table {ex1/edge_eta.dat};
%%\addplot table {ex1/vol_eta.dat};
%\addplot table {ex1/totaleta.dat};
%\addplot table {ex1/e+hf.dat};
%\addplot table {ex1/ndof.dat};
%%\legend{$\|\! {\bf p} {-} {\bf p}_\ell\!\|$,$\|\! u-u_\ell\!\|$,$e_\ell$,$\eta_\ell$,$\mu_\ell$,$\eta_\ell+\mu_\ell$,$e_\ell+\|\!hf\|\!$, ${\rm ndof}^{-\frac{1}{2}}$}tiny
%\legend{$e_\ell$,$\eta_\ell+\mu_\ell$,$e_\ell+\|\!hf\|\!$, ${\rm Ndof}^{-\frac{1}{2}}$}
%\end{loglogaxis}
%\end{tikzpicture}
\caption{Discrete solution and convergence plot}
\end{figure}
Table 5.1 displays the experimental results for  errors $\norm{u-u_\ell}$ and  $\norm{{\bf p-p}_\ell}$, the edge and volume estimators $\eta_\ell$ and $\mu_\ell$,  from (\ref{eta}) and (\ref{mu})
for several consecutive levels $\ell$  of AMFEM with the number of degrees of freedom and the marking {\bf Case (A) or (B)}. % for the Adaptive refinement
For this example, the load function $f$ exhibits a relatively large
variation in the domain and hence, within the elements. % which forces adaptive algorithm to use Case (B) of the marking 
It is observed that  the volume estimator $\mu_\ell$ is more than the edge-estimator $\eta_\ell$ in  several levels due to large  data oscillations, which results in the use of the  marking {\bf Case (B)} for most of the refinement levels. The reduction of $\mu_\ell$ influences the reduction in the edge-estimator and both lead to the optimal convergence rate ${\rm Ndof}^{-1/2}$ for errors, where Ndof denotes the number of degrees of freedom.
Fig. 5.1(b) shows the adaptive refined mesh by the AMFEM algorithm at the level 10, where the number of degrees of freedom is 6918. Fig. 5.2 depicts the approximate solution $u_h$ and summarises the convergence rates for variables ${\bf p}, u$, and sum of the estimators with respect to 
the number of degrees of freedom. % These convergence rates are optimal by adaptive MFEM algorithm.
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%\newpage
\begin{example}\label{example2}
For the {\it Crack problem} \cite{3matlab}, consider the PDE (\ref{eq1}) with coefficients $\mathbf A=I, ~{\bf b}=(x-1,y+1)$ and $\gamma=4$ on%$\gamma=10^{-2}$ on
$~\Omega =\{(x,y)\in \mathbb{R}^2:|{{\bf x}}|\leq 1\setminus[0, ~1]\times \{0\}\}$ with Dirichlet boundary condition and
  exact solution $~u(r,\theta)=r^{{1}/{2}} \sin{\theta}/{2} -{r^2}/{2} \sin ^2(\theta)$. % (in polar coordinates).
\end{example}   
The initial uniform  triangulation  $\T_0$ has  mesh-size $h=0.25$. The adaptive algorithm  is performed with parameters  $\theta_A=0.3, \theta_B=0.3$ and  $\kappa=1$.    
%\begin{figure}[!h]
%\centering
%\subfloat[ ] {\includegraphics[width=12.0cm]{crack_err.png}}\\
%  %{\includegraphics[width=10.0cm]{crack_refinement.png}}
%\caption {Convergence of error and estimator with respective to NDof}
%\end{figure}
%\newpage
\begin{table}
%\tiny
\footnotesize
\begin{center}
  \begin{tabular}{ c | c | c  |c|c |c}
     \hline
   Ndof &$\norm{u-u_\ell}$ &$\norm{{p-p}_\ell}$&$\eta_\ell$&$\mu_\ell$& Marking Case\\
 \hline \hline%    4 & 5 & 6 \\ \hline
%    7 & 8 & 9 \\
%    \hline
%  \end{tabular}
%\end{center}
%==========================================================\\
%$Ndof~~~H1~~~L2~~H1+L2~~~H1+L2+sqrt(Thf)~~\eta+\mu~~~~~markcase$\\
%==========================================================\\
46 &0.0880&    0.4021 &   2.0875 &   2.5734   & A\\
 56  & 0.0899 &   0.3979 &   3.4950 &   1.9191   & B\\
98 & 0.0771  &  0.3557  &  2.3290   & 1.8010  &  B\\
161 &   0.0633   & 0.3085    &1.3106    &1.4866 &   A\\
207  &  0.0556    &0.2643   & 1.1180    &1.0874  & B\\ 
383   & 0.0429  &  0.2401   & 0.6512    &0.8751   & A\\ \hline
469    &0.0371   & 0.2042  &  0.5813    &0.6579    &A\\
628  &  0.0308   & 0.1695 &   0.4662    &0.4896    &A\\
901  &  0.0269&    0.1400&    0.3579 &   0.3431 &   B\\
1436  &  0.0207&    0.1225&    0.1971 &   0.2704 &   A\\  
1821   & 0.0183 &   0.1023 &   0.1762  &  0.1949  &  A\\
2618   & 0.0150  &  0.0844  &  0.1433   & 0.1340   & B\\ \hline
3957   & 0.0123   & 0.0746   & 0.0795    &0.1086    &A\\
5162    &0.0111&    0.0603    &0.0669&    0.0698    &A\\
7482    &0.0092 &   0.0486&    0.0528 &   0.0453    &B\\  
12145    &0.0074 &   0.0434&    0.0282 &   0.0384&    A\\
15726&    0.0065  &  0.0351 &   0.0231  &  0.0246 &   A\\
22520 &   0.0054   & 0.0280  &  0.0169   & 0.0158  & B\\ \hline 
%36407  &  0.0042&    0.0247   & 0.0094    &0.0130   & A\\
%47656   & 0.0036 &   0.0199    &0.0080    &0.0082    &A  \\     
           \hline
  \end{tabular}
\end{center}
    \caption{Numerical results of AMFEM for Example \ref{example2}}
\end{table}

\begin{figure}[!ht]%\nonumber
\centering
{\includegraphics[width=6.2cm]{cgs_plot1.pdf}}
{\includegraphics[width=6.2cm]{cgs_plot2.pdf}}

\caption{(a) Convergence plots for $\theta_A,\theta_B=0.3$ (b) Error convergence for $\theta=\theta_A=\theta_B$ equal to 0.3, 0.5, 0.8, 0.9, 1}
\end{figure}     
       Table 5.2 displays the experimental results for the errors  $\norm{u-u_\ell}$ and $\norm{{\bf p-p}_\ell}$, the estimators $\eta_\ell$ and $\mu_\ell$ 
for several consecutive levels $\ell$  of AMFEM with the number of degrees of freedom and the marking {\bf Case (A) or {\bf (B)} } for this problem. % for the Adaptive refinement
The right-hand function $f$ of (\ref{eq1}) in this example  is not smooth and the solution has singularity at the origin, and hence  the adaptive algorithm utilizes both the marking cases to achieve the optimal convergence rate. Figure 5.3(a) shows  optimal convergence rate ${\rm Ndof}^{-1/2}$ for adaptive algorithm while uniform refinement yields a suboptimal convergence  rate ${\rm Ndof}^{-1/4}$ for the above mentioned parameters. Figure 5.3(b) shows the convergence rate of the error  for  different values of $\theta_A$ and $ \theta_B$. This supports the theoretical prediction of the choice of marking parameters $\theta_A$ and $\theta_B$ in the proof of Theorem 4.5, that is, the quasi-optimality can be 
achieved only for  small values of  $\theta_A~ {\rm and}~ \theta_B$.
\begin{example}\label{example3} % \cite{3matlab}
Consider the PDE (\ref{eq1}) with  coefficients  $~ \mathbf A=I,~~\mathbf {b} =(0, 0),$~ %(r \cos\theta,r \sin\theta)$ 
$\gamma=-8.9$,  Dirichlet boundary condition on the L-shaped domain
$\Omega =(-1,1)\times (-1,1)\setminus [0,1]\times [-1,0]$ and the exact solution  given in polar coordinates as $u(r,\theta)=r^{2/3} \sin\big( 2\theta/3 \big).$
 \end{example}    
%\begin{figure}[!h]%\nonumber
%\centering
%\subfloat  []{\includegraphics[width=8.0cm]{Ex_1_refinement.png}} 
%%\subfloat{\includegraphics[width=8.0cm]{ex_1_plot.png}}
%\caption{(a) Adaptive refinement (b) Convergence  }%$C_{rel}$ for different $\gamma$ with adaptive refinement}
%\end{figure}
\begin{figure}[!ht]%\nonumber
\centering
{\includegraphics[width=6.0cm]{Ex_1_refinement.png}} 
{\includegraphics[width=6.2cm]{cgs_plot4.pdf}} 
%
%\begin{tikzpicture}[scale=0.7]
%\begin{loglogaxis}[
%%title=Convergence Plot,
%xlabel={Number of degrees of freedom},
%ylabel={ Errors and estimators},%}
%font=\sffamily\scriptsize,
%]
%%\addplot table {ex1/datH1.dat};
% %\addlegendentry{ $\|\! {\bf p} {-} {\bf p}_h\!\|$ uniform}
%\addplot table {gamma123/H1.dat};
%\addplot table {gamma123/L2H1.dat};
%%\addplot table {ex1/edge_eta.dat};
%%\addplot table {ex1/vol_eta.dat};
%\addplot table {gamma123/te.dat};
%%\addplot table {ex1/e+hf.dat};
%\addplot table {gamma123/sqndof.dat};
%legend style={at={(-0.1,.5)}, anchor=south west},
%\legend{$\|\! {\bf p} {-} {\bf p}_\ell\!\|$,$e_\ell$,$\eta_\ell+\mu_\ell$ , ${\rm ndof}^{-\frac{1}{2}}$}
%%\legend{$\|\! {\bf p} {-} {\bf p}_\ell\!\|$,$\|\! u-u_\ell\!\|$,$e_\ell$,$\eta_\ell$,$\mu_\ell$,$\eta_\ell+\mu_\ell$,$e_\ell+\|\!hf\|\!$, ${\rm ndof}^{-\frac{1}{2}}$}tiny
%\end{loglogaxis}
%\end{tikzpicture}
\caption{Adaptive mesh-refinement and the convergence plot}
\end{figure}     
The numerical experiment  is performed with initial mesh-size $h=0.25$ and uniform triangulation. The adaptive algorithm  is performed with parameter choice $\theta_A=0.5, \theta_B=0.5$ and  $\kappa=2$.
Since  the  forcing function $f$ is smooth, the oscillation term has higher order convergence than the remaining terms in the estimator. In this example, the marking is based on the edge estimator  to obtain the optimal convergence. %, rather than the  higher order  estimator terms. %(less convergence rate) 
 The refinement due to the data oscillations and the volume estimators  can be avoided since they play a minor role in the convergence. % and hence  to these terms can be avoided.
        Figure 5.4 displays the adaptive mesh-refinement and the optimal convergence rate ${\rm Ndof}^{-1/2}$  of the errors and the estimators for a sufficiently small mesh-size $h$.
% If  we consider all estimator terms together, it may mark extra number of elements
%at every level, while due to separate marking help to obtain 
%convergence in less NDOF.
\subsection*{5.4.  Observations}
This subsection deals with a few observations.
\begin{itemize}
 \item If the given function $ f$ has large variation in the domain and within the elements, 
 then adaptive algorithm chooses 
    {\bf Case (B)} marking criteria and achieves the convergence (see 5.1).%which forces Algorithm C to refine additional
%elements due to data oscillation (step 5 of Algorithm C), not yet marked for refinement\\
%3) Figure shows the optimal convergence rate $N^{-1/2}$ for adaptive algorithm while uniform refiment gives the suboptimal $N^{-1/4}$ .\\
\item When $f$ is not smooth and the solution also has singularity, then the adaptive algorithm utilizes both the marking cases to achieve the optimal convergence rate (see 5.2).
\item For the  smooth function $f$, the oscillation has higher order convergence than the remaining  estimator terms.
 In such cases, the data oscillations and the  volume estimators 
 have a very minor role in the convergence. Thus, the marking would be based on the edge-based error estimator to capture singularity of the solution (see 5.3).  %Hence it will be beneficial  to obtain the optimal convergence.% and convergence may be attain in less Ndof.\\ 
%Downloaded
%On the other hand, it explains why data oscillation
%seems to play a minor role for (piecewise) smooth data f, and hints at the underlying
%reasons why most adaptive strategies, although neglecting data oscillation, converge
%in practice
%\end{document}
\end{itemize}
\section{Conclusions} In this work, the convergence and the quasi-optimality  of adaptive mixed finite element method
is analysed for the non-symmetric and indefinite second order elliptic equations using the lowest order Raviart-Thomas elements. The adaptive algorithm in Subsection 2.1 is designed as a combination of  the edge and the volume error estimators. The numerical experiments confirm the efficiency of this algorithm and  support the theoretical findings.
\subsection*{Acknowledgements} The first author acknowledges the financial support of  Council of Scientific and Industrial Research (CSIR), Government of India. The authors sincerely thank Professor Carsten Carstensen,  Humboldt
University, Berlin  for his constructive comments and suggestions.
%\end{acknowledgements}
\begin{thebibliography}{1}
%\bibitem{GoMiSa} {\sc M. Goossens, F. Mittelbach, and A. Samarin},
%{\em The} \LaTeX\ {\em Companion}, Addison-Wesley, Reading, MA, 1994.
%
%\bibitem{Higham} {\sc N.~J. Higham}, {\em Handbook of Writing for
%the Mathematical Sciences}, Society for Industrial and Applied
%Mathematics, Philadelphia, PA, 1993.
%
%\bibitem{Lamport} {\sc L. Lamport}, \LaTeX: {\em A Document
%Preparation System}, Addison-Wesley, Reading, MA, 1986.
%
%\bibitem{SerLev} {\sc R. Seroul and S. Levy}, {\em A
%Beginner's Book of} \TeX, Springer-Verlag, Berlin, New
%York, 1991.

\bibitem{oden}
{\sc M. Ainsworth and J. T. Oden}, {\em A posteriori error estimation in finite element analysis},
Wiley, New York, 2000.

%\bibitem{arnold}
% {\sc Arnold, N. Douglas, and R. S. Falk},
%  {\em A uniformly accurate finite element method for the {R}eissner-{M}indlin plate},
% SIAM J. Numer. Anal., {26} ({1989}), pp. {1276-1290}.
 \bibitem{3matlab}
{\sc C. Bahriawati and C. Carstensen},
{\em Three matlab implementation of the lowest order Raviart-Thomas MFEM with a posteriori error control},
Comput. Methods Appl. Math., 5 (2005), pp. 333-1361.
% \bibitem{Braess}
% D. Braess,
% {\it Finite elements- theory, fast solvers and applications in solid mechanics},
% Cambridge university press, 1992.

\bibitem{Brenner-Scott}
{\sc S. C. Brenner and L. R. Scott},
{\em The mathematical theory of finite element methods},
third ed., Texts in Applied Mathematics, 15, Springer, New York, 2008.
% \bibitem{Brennerenrichment}
% S.~C. Brenner.
% \newblock Two-level additive {S}chwarz preconditioners for  
% nonconforming finite
%   element methods.
% \newblock {\em Math. Comp.}, 65(215):897--921, 1996.
\bibitem{bek_mixed}
 {\sc R. Becker  and S. Mao},
   {\em An optimally convergent adaptive mixed finite element method},
  {Numer. Math.},  111 (2008), pp. 35-54.
  %FJOURNAL = {Numerische Mathematik},
%    VOLUME = {111},
%      YEAR = {2008},
%    NUMBER = {1},
%     PAGES = {35--54},
\bibitem{bek_non}
 {\sc R. Becker  and S. Mao},
     {\em A convergent nonconforming adaptive finite element method with
              quasi-optimal complexity},
   {SIAM J. Numer. Anal.},
  %FJOURNAL = {SIAM Journal on Numerical Analysis},
  {47} ({2010}), pp.   {4639-4659}.
 \bibitem{biven}
 {\sc P. Binev, W. Dahmen, and  R. DeVore},
    {\em Adaptive finite element methods with convergence rates},
    {Numer. Math.},
  %FJOURNAL = {Numerische Mathematik},
     {97} ({2004}), pp. {219-268}.
\bibitem{Brezzi}
{\sc F. Brezzi and  M. Fortin},
{\em Mixed and hybrid finite element methods},
Springer Verlag, 1991.

\bibitem{Car_apo}
{\sc C. Carstensen},
{\em A~posteriori error estimate for the mixed finite element method},
 Math. Comp., 66 (1997), pp. 465-476.

\bibitem{Car1}
{\sc C. Carstensen},
 {\em A unifying theory of a posteriori finite element error control},
Numer. Math., 100 (2005), pp. 617-637.

\bibitem{Car_apo1}
{\sc C. Carstensen},
{\em Convergence of adaptive finite element methods in computational mechanics},
Appl. Numer. Math., 59 (2009), pp.  2119-2130.
\bibitem{CH}
{\sc C. Carstensen and J. Hu},
{\em A unifying theory of  a~posteriori  error control for nonconforming  finite element methods},
Numer. Math., 107 (2007), pp. 473-502. 
 
% \bibitem{Carhop}
% C. Carstensen, M. Eigel, C. L{\"o}bhard and R.H.W. Hoppe,
% {\it A review of unified a posteriori finite element error control},
% Numer. Math. Theory Methods Appl.,  5, 2012, 509-558. 
% 
\bibitem{CC}
{\sc C. Carstensen, Asha. K. Dond, N. Nataraj, and Amiya. K. Pani},
{\em Error analysis of nonconforming and mixed FEMs for
second-order linear non-selfadjoint and indefinite elliptic
problems}, {\em http://arxiv.org/abs/1401.4810}.

 
% 
% \bibitem{CHZ}
% P. Ciarlet, Jr., J. Huang and J. Zou,
% {\it Some observations on generalized saddle-point problems},
% SIAM J Matrix Anal. Appl. 25, 2003, 224-236.
% 
% 
% \bibitem{Doug}
% J. Douglas, Jr. and J. E. Roberts,
% {\it Global estimates for mixed methods for second order elliptic equations},
% Math. Comp. 44, 1985, 39-52.
% 
% \bibitem{dauge}
% M. Dauge, {Elliptic boundary problems on corner domains},
% Lecture Notes in Math. Vol. 1341, Springer verlag, Berlin 1988.
% 
% \bibitem{evans}
% L.C. Evans,
% {\it Partial Differential Equations},
% Graduate Studies in Mathematics Vol. 19, AMS, Providence, Rhode Island, 2002.
% 
% 
% \bibitem{GT}
% D. Gilberg, N. S. Trudinger,
% {\it Elliptic partial differential equations of second order},
% Springer Verlag, 1982.
% 
% \bibitem{G}
% P. Grisvard,
% {\it Elliptic Problems in Nonsmooth Domains}, Pitman, Boston, 1985.
\bibitem{car_cgt2}
 {\sc C. Carstensen and R. Hella},
     {\em An optimal adaptive mixed finite element method},
   {Math. Comp.},
  %FJOURNAL = {Mathematics of Computation},
     {80} ({2011}), pp. {649-667}.
\bibitem{car_non_cgt}
{\sc C. Carstensen and R.H.W. Hoppe},
{\em Convergence analysis of adaptive nonconforming finite element methods},
Numer. Math., 103 (2006), pp. 251-266.

\bibitem{car_cgt1}
    {\sc C. Carstensen and R. H. W.  Hoppe},
    {\em Error reduction and convergence for an adaptive mixed finite
              element method},
  {Math. Comp.}, {75} ({2006}), pp.  {1033-1042}.
      
\bibitem{axiom}
{\sc C. Carstensen, M. Feischl, M.Page, and  D. Praetorius}, {\em Axioms of Adaptivity}, ASC Report 38/2013, Institute for Analysis and Scientific Computing, Vienna University of Technology.
    
    \bibitem{quasi}
 {\sc J. Cascon, C.  Kreuzer, R. H. Nochetto,  and  K. G. Siebert},
     {\em Quasi-optimal convergence rate for an adaptive finite element
              method},
  {SIAM J. Numer. Anal.},
 % FJOURNAL = {SIAM Journal on Numerical Analysis},
 {46}  (2008),  pp. {2524-2550}.
 \bibitem{cascon}
 {\sc J. Cascon, R. H. Nochetto, and  K. G. Siebert},
     {\em Design and convergence of AFEM in $H({\rm div})$},
  {Math. Models Methods Appl. Sci.},
 % FJOURNAL = {SIAM Journal on Numerical Analysis},
 {17} ({2007}), pp. {1849-1881}.
    
\bibitem{Chen}
{\sc J. Chen and L. Li},
 {\em Convergence and domain decomposition algorithm for
 nonconforming and mixed methods for nonselfadjoint and indefinite
 problems},
 Comput. Methods Appl. Engrg., 173 (1999), pp. 1-20.
%
\bibitem{chen_cgs}
{\sc L. Chen, M. Holst, and J. Xu}, 
{\em Convergence and optimality of adaptive mixed finite element methods},
 Math. Comp., 78 (2009), pp. 35-53. 

\bibitem{chen_hoppe}
{\sc H. Chen, X. Xu, and R.H.W. Hoppe},
{\em Convergence and Quasi-optimality of adaptive nonconforming finite element methods for some nonsymmetric and indefinite problems},
Numer. Math., 116  (2010), pp. 383-419.

\bibitem{l2}
     {\sc A. Demlow  and R Stevenson},
     {\em Convergence and quasi-optimality of an adaptive finite element  method for controlling {$L_2$} errors},
    {Numer. Math.},
  %FJOURNAL = {Numerische Mathematik},
    {117} ({2011}), pp. {185-218}.
\bibitem{daglas}
{\sc J. Jr. Douglas and J.E. Robert}, {\em Global estimates for mixed methods for second order elliptic equations}
Math. Comput. 44(1985), pp. 39-51.
%\bibitem{Du_xie}
%    {\sc S. Du and X. Xie},
%   {\em Error reduction, convergence and optimality of an adaptive
%              mixed finite element method},
%   {J. Syst. Sci. Complex.},
%  %FJOURNAL = {Journal of Systems Science \& Complexity},
%   {25} (2012), pp.  {195-208}.
%   
\bibitem{prac}
{\sc M. Feischl, T. F$\ddot{u}$hrer, and D. Praetorius}, {\em Adaptive FEM with optimal convergence rates for a certain class of non-symmetric and possibly non-linear problems}, ASC Report 43/2012, Institute for Analysis and Scientific Computing, Vienna University of Technology.

\bibitem{nochetto}
{\sc K. Mekchay and R.H. Nochetto},
{\em Convergence of adaptive finite element methods for general second order elliptic PDE},
SIAM J. Numer. Anal., 43 (2005), pp. 1803-1827.

\bibitem{Marini}
 {\sc   L. D.  Marini},
    {\em An inexpensive method for the evaluation of the solution of
              the lowest order {R}aviart-{T}homas mixed method},
{SIAM J. Numer. Anal.},
  %FJOURNAL = {SIAM Journal on Numerical Analysis},
     {22} ({1985}), pp. {493-496}.
      
\bibitem{dof} 
{\sc W. D$\ddot{o}$rfler},
{\em A convergent adaptive algorithm for Poisson's equation},
SIAM J. Numer. Anal. 33 (1996), pp. 1106-1124. 
 
%\bibitem{mao}
%  {\sc S. Mao, X. Zhao, and Z. Shi},
%     {\em Convergence of a standard adaptive nonconforming finite
%              element method with optimal complexity},
% {Appl. Numer. Math.},
 % FJOURNAL = {Applied Numerical Mathematics. An IMACS Journal},
 %{60} ({2010}), pp. {673-688}.
\bibitem{morin1}
{\sc P. Morin, R.H. Nochetto, and K.G. Siebert}, 
{\em Data oscillation and convergence of adaptive fem},
SIAM J. Numer. Anal., 38 (2000), pp. 466-488. 

\bibitem{morin2}
{\sc P. Morin, R.H. Nochetto, and K.G. Siebert}, 
{\em Convergence of adaptive finite element methods},
SIAM Rev.,  44 (2002), pp. 631-658. 

\bibitem{veeser}
  {\sc P. Morin, K. G. Siebert, and A. Veeser},
    {\em A basic convergence result for conforming adaptive finite
              elements},
{Math. Models Methods Appl. Sci.},
 % FJOURNAL = {Mathematical Models \& Methods in Applied Sciences},
  {18} ({2008}),  pp.  {707-737}.
    

\bibitem{schatz1}
{\sc A. H. Schatz}, {\em  An observation concerning Ritz-Galerkin methods with indefinite bilinear forms}, Math.
Comput., 28 (1974), pp. 959-962.

\bibitem{schatz}
{\sc A. H. Schatz and J. Wang},
{\em Some new error estimates for Ritz-Galerkin methods with minimal regularity assumptions}, Math. Comp., 65 (1996), pp. 19-27.

%\bibitem{scottzhang}
%{\sc L. R. Scott and S.Zhang}, {\em Finite element interpolation of nonsmooth function satisfying boundary conditions}, Math. Comp., 54 (1990), pp. 483-493.

\bibitem{steven}
{\sc R. Stevenson}, 
{\em   Optimality of a standard adaptive finite element method},
Found. Comput. Math., 7 (2007), pp.245-269.

\bibitem{verf}
{\sc R. Verf\"{u}rth},
{\em A posteriori error estimation and adaptive mesh-refinement techniques},
J. Comput. Appl. Math., 50 (1994), pp. 67-83.


\end{thebibliography}
%\newpage
{\footnotesize
  \begin{center}
  {\bf Appendix I}\\
{List of the constants}
\begin{longtable}{| c | p{6cm} |  p{2.8cm}|p{1.1cm}|}

%\caption[Feasible triples for a highly variable Grid]{Feasible triples for 
%highly variable Grid, MLMMH.} \label{grid_mlmmh} \\
%Constant & Dependency on & Appears first in &  Value \\
\hline \multicolumn{1}{|c|}{\textbf{Constant}} & \multicolumn{1}{p{6cm}|}{\textbf{Dependency on }} & \multicolumn{1}{p{2.8cm}|}{\textbf{Appears first in}} & \multicolumn{1}{p{1.1cm}|}{\textbf{Value }} \\ \hline 
\endfirsthead

\multicolumn{4}{c}%
{{\bfseries \tablename\ \thetable{} -- continued from previous page}} \\
\hline \multicolumn{1}{|c|}{\textbf{Constant}} &
\multicolumn{1}{p{6cm}|}{\textbf{Dependency on }} &
\multicolumn{1}{p{2.8cm}|}{\textbf{Appears first in}}  &
\multicolumn{1}{p{1.1cm}|}{\textbf{Value }} \\ \hline 
\endhead

\hline \multicolumn{4}{|r|}{{Continued on next page}} \\ \hline
\endfoot

\hline \hline
\endlastfoot
&{\bf I. Natural Constants~~~~} && \\ \hline
$ C_{\rm rel}$   & Coefficients of (1.1) and interpolation  &    &  Positive \\
&constants&Theorem 3.2 &\\ \cline{1-2} \cline{4-4}
$\epsilon$   & Coefficients of (1.1) and interpolation  & ({\it A~posteriori }   &  Positive \\
&constants& {\it estimates}) &\\  \cline{1-2} \cline{4-4}
$h_1$& $\epsilon$& &  Positive \\ \hline
$  C_{\rm eff} $ & Coefficients of (1.1),  inverse inequality  &   Lemma 3.3 &  Positive \\  %\cline{1-4}%, finite overlap
&and finite overlap&({\it efficiency})  &\\ \cline{1-3} \cline{3-4}
$ C_{1}$   & Coefficients of (1.1)&   Lemma 3.4 ({\it Volume} &  Positive \\
&&{\it  estimator reduction})& \\ \cline{1-3}\cline{4-4}
$  C_{2} $ & Coefficients of (1.1),  inverse inequality  &   Lemma 3.6 ({\it edge}&  Positive \\ 
&and finite overlap&({\it  estimator reduction})& \\\cline{1-4}%, finite overlap
$ C_{3}$   & Coefficients of (1.1)  &  Lemma 3.7 ({\it Quasi-} &  Positive \\ 
 &&{\it orthogonality})&\\ \cline{1-4}
$  C_{4} $ & Coefficients of (1.1), &  Lemma 3.9 &  Positive \\
& interpolation constants &({\it Quasi-reliability})&\\ \cline{1-4}%, finite overlap
$ \gamma_0$      &$ \left \{
  \begin{array}{l l}
    0<\gamma_0<1 & ~~ \text{if $T\in \mathcal T_H,$ refined}\\
    1 & \quad \text{otherwise}
  \end{array} \right.$& Corollary  3.5 &  Positive \\ \hline \hline 
  & {\bf II. Contraction property}& & \\ \hline 
 $ \theta_B$     & $0<\theta_B <1$, to be chosen          &  &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \delta_1$    &$<\frac{\theta_B}{4-\theta_B}$     &  Lemma 3.4  &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \rho_B$    &$(1+\delta_1)(1-\frac{\theta_B}{4})$     &  ({\it volume-estimator} &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \varLambda_1$ & $C_1(1+\frac{1}{\delta_1})$     & {\it reduction }) &  (0,1)\\ \hline%  \cline{1-2} \cline{4-4}
%
$ \theta_A$     &  $0<\theta_A <1$, to be chosen                                &   & (0,1) \\  \cline{1-2} \cline{4-4}
$ \delta_2$     &$\frac{ \theta_A}{4-2 \theta_A}$ & Lemma  3.6  &  (0,1)\\   \cline{1-2} \cline{4-4}
$ \varLambda_2$ & $C_2(1+\frac{1}{\delta_2})$     & ({\it edge-estimator}&  (0,1)\\   \cline{1-2} \cline{4-4}
$ \rho_A$       &$(1+\frac{1}{\delta_2})(1-\frac{\theta_A}{2})$&   {\it reduction }) &  (0,1) \\  \hline
%
$C_3$            & Coefficient of (1.1) &                         &Positive \\   \cline{1-2} \cline{4-4}
$\alpha_2$       & $0.5+7C_3^2\varLambda_1 H^2$ &  Lemma 3.6 &  (0,1) \\   \cline{1-2} \cline{4-4}
$ \alpha_3$      & $0.5-7C_3^2\varLambda_1 H^2$ &  ({\it Quasi-} &  (0,1) \\   \cline{1-2} \cline{4-4}
$\delta_3$       & $<$1                              & {\it orthogonality}) & (0,1) \\  \cline{1-2} \cline{4-4}
$\alpha_1$       &$\delta_3+4C_3^2\epsilon^2$                             &&   (0,1) \\  \hline
%
$\alpha_4$       & $\frac{1}{2}\min \{1, \frac{(1-\rho_A) \alpha_3}{({\varLambda_1+ \varLambda_2 }) C_{\rm rel}}\}$&  &  (0,1) \\   \cline{1-2} \cline{4-4}
$\delta_3$        &$<\alpha_4$&Remark 3.2 &(0,1)\\ \cline{1-2} \cline{4-4}
 $\epsilon$       & $(\alpha_4-\delta_3)/(4 C_3^2 )$ &&Positive\\ \cline{1-2} \cline{4-4}
         $h_2$    & $\epsilon ~{\rm and}~ < \frac{1}{\sqrt{14C_3^2\varLambda_1}} $ &&Positive\\ \cline{1-2} \cline{4-4}
$\alpha_1$      &$<\alpha_4$&& (0,1) \\ \hline
%&&&&
$C_5$& $({\varLambda_1+ \varLambda_2 })/{\alpha_3}$&&Positive\\  \cline{1-2} \cline{4-4}
$\alpha$  & $C_5(1-\alpha_1)$&                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$h_3$  & $\min \{ h_1, h_2, \frac{1}{\sqrt{\beta}} \}$&                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$E$       &$2(C_5 \varLambda_3+C_5\alpha_4 C_{rel})$ &  Theorem 4.1  &  Positive \\ \cline{1-2} \cline{4-4}
$\kappa$  & $<\kappa_0:= \frac{({1-\rho_A-C_5\alpha_4C_{rel}})}{E}$&   &  Positive \\  \cline{1-2} \cline{4-4}
$D$&$({1.5+C_5\alpha_4C_{rel}})({\kappa})$ &( {\it Contraction} &  Positive \\  \cline{1-2} \cline{4-4}
 $ \delta_1$   &$\min\{ \frac{E}{2\beta} ,\frac{\theta_B}{4-\theta_B}\}$ &   {\it property} )&  (0,1) \\ \cline{1-2} \cline{4-4}
$\beta$& $2 \max\big{\{}C_5 \varLambda_3 + C_5 \alpha_4 C_{rel}, \frac{C_5\varLambda_3+C_5\alpha_4C_{rel}+D}{1-\rho_B}\big{\}}$&                                 
                                                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$\rho_1$& $\max \Big{\{ } \frac{C_5(1-\alpha_4)}{\alpha}, \rho_A+C_5\alpha_4C_{rel}+E \kappa,$&    &  (0,1)\\ %\cline{2-2} \cline{4-4}
&\qquad  \hspace{2cm}$ \frac{((1+\delta_1)\beta-E/2)}{\beta} \Big \}$&&\\  \cline{1-2} \cline{4-4}
$\rho_2$&$\max \Big\{  \frac{C_5(1-\alpha_4)}{\alpha}, 2+C_5\alpha_4C_{rel}-D\kappa,$  & &  (0,1) \\ % \cline{4-4}
&\qquad  \hspace{2cm}$ \frac{\rho_B\beta+C_5 \varLambda_3+C_5\alpha_4C_{rel}+D }{\beta} \Big \}$&&\\ \cline{1-2} \cline{4-4}
$\rho$& $\max\{\rho_1,\rho_2\}$&  & (0,1) \\  \hline  \hline
  & {\bf III. Quasi-Optimality~}& & \\ \hline 
$h_4$&  $ \frac{1}{2(7+C_3 \varLambda_1+{\beta \varLambda_1})/({C_{\rm eff}^{-1}+\alpha})}$  & Lemma 4.3  &  Positive \\ \cline{1-3} \cline{4-4}
$ \theta_A$     & $\min\{1,\frac{C_{\rm eff}}{9C_4}\}$&   & (0,1) \\ \cline{1-2} \cline{4-4}
$ \theta_B$     & $({\gamma_0^2})/({8(2+C_4C_1h^2)})$& Theorem  4.5 &  (0,1) \\  \cline{1-2} \cline{4-4} %$\frac{\gamma_0^2}{24} <\frac{\gamma_0^2}{8(2+C_4C_1h^2)}$
$\kappa$& $   \min \{ \kappa_0,  \frac{{C_{\rm eff}}/{8}-C_4\theta_A }{4C_4+ {\varLambda_3}/{2}}\big\}$&  ({\it Quasi-} &  Positive \\  \cline{1-2} \cline{4-4}
%$\tau_1$& $\min\{ \tau_2, \tau_3 \}$& {optimality}) &  Positive \\ \cline{1-2} \cline{4-4}
$\epsilon$& $\min\{ \epsilon, 1, \frac{1}{4C_4} \}$&{optimality}) &  Positive \\ \cline{1-2} \cline{4-4}
$h_5,h_7$& $ \epsilon $& &  Positive \\ \cline{1-2} \cline{4-4}
$\tau_2^2$& $ \frac{(c_{\rm eff}/8)-\kappa ( 4C_4+ {\varLambda_3}/{2})-C_5 \theta_A}{(5/4)C_7} $&  &  Positive \\ \cline{1-2} \cline{4-4}
$h_6$&  $\epsilon $ and $\min\{ h_2, h_4,\frac{1}{\sqrt{4C_4 \varLambda_1}} \}$&  &  Positive \\ \cline{1-2} \cline{4-4}
$\tau_3^2$& $\frac{\beta \gamma_0^2}{8(2+C_4C_1h^2)C_8},C_8=C_4C_1(C_{\rm rel +1})(1+\frac{1}{\kappa}) $ &  &  Positive \\ \cline{1-2} \cline{4-4}
%& $C_8=C_4C_1(C_{\rm rel +1})(1+\frac{1}{\kappa})$&&\\  \cline{1-2} \cline{4-4}
$h_8$&  $\epsilon $ and $\min\{ h_2,h_4, \sqrt{\frac{\alpha}{\beta C_1C_4 }},\sqrt{\frac{\gamma_0^2}{4C_8 }} \}$&  &  Positive \\ \cline{1-2} \cline{4-4}
$h_0$& $\min\{ h_7, h_8 \}$& &  Positive \\ \cline{1-2} \cline{4-4}
$\epsilon_1$& $\tau_1 \xi_\ell $, where $\tau_1:=\min\{ \tau_2, \tau_3 \}$&   &  Positive 
\end{longtable}
\end{center} }

\end{document}




\begin{table}[!ht]
%\tiny
\caption*{\bf Appendix }
\caption*{List of the constants}
\begin{center}
{\bf I. Natural Constants~~~~} \\
{\footnotesize
  \begin{tabular}{| c | p{5cm} | c  |p{1.2cm}|}%\footnotesize%\label{tb1}
     \hline
 Constant & Dependency on & Appears first in &  Value \\
 \hline \hline
$ C_{\rm rel}$   & Coefficients of (1.1) and interpolation  &    &  Positive \\
&constants& &\\ \cline{1-2} \cline{4-4}
$\epsilon$   & Coefficients of (1.1) and interpolation  &   Theorem 3.2 &  Positive \\
&constants&({\it A~posteriori estimates}) &\\  \cline{1-2} \cline{4-4}
$h_1$& $\epsilon$& &  Positive \\ \hline
$  C_{\rm eff} $ & Coefficients of (1.1),  inverse inequality  &   Lemma 3.3 &  Positive \\  %\cline{1-4}%, finite overlap
&and finite overlap&({\it efficiency})  &\\ \cline{1-3} \cline{3-4}
$ C_{1}$   & Coefficients of (1.1)&   Lemma 3.4 &  Positive \\
&&({\it Volume estimator reduction})& \\ \cline{1-3}\cline{4-4}
$  C_{2} $ & Coefficients of (1.1),  inverse inequality  &   Lemma 3.6 &  Positive \\ 
&and finite overlap&({\it edge estimator reduction})& \\\cline{1-4}%, finite overlap
$ C_{3}$   & Coefficients of (1.1)  &  Lemma 3.7 &  Positive \\ 
 &&({\it Quasi-orthogonality})&\\ \cline{1-4}
$  C_{4} $ & Coefficients of (1.1), &  Lemma 3.9 &  Positive \\
& interpolation constants &({\it Quasi-reliability})&\\ \cline{1-4}%, finite overlap
$ \gamma_0$      &$ \left \{
  \begin{array}{l l}
    0<\gamma_0<1 & ~~ \text{if $T\in \mathcal T_H,$ refined}\\
    1 & \quad \text{otherwise}
  \end{array} \right.$& Corollary  3.5 &  Positive \\ \hline 
%&& $1 ~~~~{\rm otherwise}$ && \\ \hline  \hline
\end{tabular} }
\end{center}
 \end{table}
\begin{table}
\begin{center}
{\footnotesize 
{\bf II. Contraction property}\\
\begin{tabular}{|p{1cm}| p{6cm} | p{2.7 cm} |p{1.5cm}|}%\footnotesize%\label{tb1}
         \hline
 Constant & Dependency on & Appears first in &  Value \\
 \hline \hline
 %Constant & Dependency on & Appears first in &  Value \\
 %\hline %
 $ \theta_B$     & $0<\theta_B <1$, to be chosen          &  &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \delta_1$    &$<\frac{\theta_B}{4-\theta_B}$     &  Lemma 3.4  &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \rho_B$    &$(1+\delta_1)(1-\frac{\theta_B}{4})$     &  ({\it volume-estimator} &  (0,1) \\  \cline{1-2} \cline{4-4}
$ \varLambda_1$ & $C_1(1+\frac{1}{\delta_1})$     & {\it reduction }) &  (0,1)\\ \hline%  \cline{1-2} \cline{4-4}
%
$ \theta_A$     &  $0<\theta_A <1$, to be chosen                                &   & (0, 1) \\  \cline{1-2} \cline{4-4}
$ \delta_2$     &$\frac{ \theta_A}{4-2 \theta_A}$ & Lemma  3.6  &  (0,1)\\   \cline{1-2} \cline{4-4}
$ \varLambda_2$ & $C_2(1+\frac{1}{\delta_2})$     & ({\it edge-estimator}&  (0,1)\\   \cline{1-2} \cline{4-4}
$ \rho_A$       &$(1+\frac{1}{\delta_2})(1-\frac{\theta_A}{2})$&   {\it reduction }) &  (0, 1) \\  \hline
%
$C_3$            & Coefficient of (1.1) &                         &Positive \\   \cline{1-2} \cline{4-4}
$\alpha_2$       & $0.5+7C_3^2\varLambda_1 H^2$ &  Lemma 3.6 &  (0,1) \\   \cline{1-2} \cline{4-4}
$ \alpha_3$      & $0.5-7C_3^2\varLambda_1 H^2$ &  ({\it Quasi-} &  (0,1) \\   \cline{1-2} \cline{4-4}
$\delta_3$       & $<$1                              & {\it orthogonality}) & (0,1) \\  \cline{1-2} \cline{4-4}
$\alpha_1$       &$\delta_3+4C_3^2\epsilon^2$                             &&   (0,1) \\  \hline
%
$\alpha_4$       & $\frac{1}{2}\min \{1, \frac{(1-\rho_A) \alpha_3}{({\varLambda_1+ \varLambda_2 }) C_{\rm rel}}\}$&  &  (0,1) \\   \cline{1-2} \cline{4-4}
$\delta_3$        &$<\alpha_4$&Remark 3.2 &(0,1)\\ \cline{1-2} \cline{4-4}
 $\epsilon$       & $(\alpha_4-\delta_3)/(4 C_3^2 )$ &&\\ \cline{1-2} \cline{4-4}
         $h_2$    & $\epsilon ~{\rm and}~ < \frac{1}{\sqrt{14C_3^2\varLambda_1}} $ &&\\ \cline{1-2} \cline{4-4}
$\alpha_1$      &$<\alpha_4$&& (0,1) \\ \hline
%&&&&
$C_5$& $({\varLambda_1+ \varLambda_2 })/{\alpha_3}$&&Positive\\  \cline{1-2} \cline{4-4}
$\alpha$  & $C_5(1-\alpha_1)$&                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$h_3$  & $\min \{ h_1, h_2, \frac{1}{\sqrt{\beta}} \}$&                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$E$       &$2(C_5 \varLambda_3+C_5\alpha_4 C_{rel})$ &  Theorem 4.1  &  Positive \\ \cline{1-2} \cline{4-4}
$\kappa$  & $<\kappa_0:= \frac{({1-\rho_A-C_5\alpha_4C_{rel}})}{E}$&   &  Positive \\  \cline{1-2} \cline{4-4}
$D$&$\frac{1.5+C_5\alpha_4C_{rel}}{\kappa}$ &( {\it Contraction} &  Positive \\  \cline{1-2} \cline{4-4}
 $ \delta_1$   &$\min\{ \frac{E}{2\beta} ,\frac{\theta_B}{4-\theta_B}\}$ &   {\it property} )&  (0,1) \\ \cline{1-2} \cline{4-4}
$\beta$& $2 \max\big{\{}C_5 \varLambda_3 + C_5 \alpha_4 C_{rel}, \frac{C_5\varLambda_3+C_5\alpha_4C_{rel}+D}{1-\rho_B}\big{\}}$&                                 
                                                                                      &  Positive \\  \cline{1-2} \cline{4-4}
$\rho_1$& $\max \Big{\{ } \frac{C_5(1-\alpha_4)}{\alpha}, \rho_A+C_5\alpha_4C_{rel}+E \kappa,$&    &  (0,1)\\ %\cline{2-2} \cline{4-4}
&\qquad  \hspace{2cm}$ \frac{((1+\delta_1)\beta-E/2)}{\beta} \Big \}$&&\\  \cline{1-2} \cline{4-4}
$\rho_2$&$\max \Big\{  \frac{C_5(1-\alpha_4)}{\alpha}, 2+C_5\alpha_4C_{rel}-D\kappa,$  & &  (0,1) \\ % \cline{4-4}
&\qquad  \hspace{2cm}$ \frac{\rho_B\beta+C_5 \varLambda_3+C_5\alpha_4C_{rel}+D }{\beta} \Big \}$&&\\ \cline{1-2} \cline{4-4}
$\rho$& $\max\{\rho_1,\rho_2\}$&  & (0,1) \\  \hline
%& & Lemma &  Positive \\  \hline
%& &  Lemma &  Positive \\  \hline
%& &  Lemma &  Positive \\  \hline
  \end{tabular} }
  \end{center}
 \end{table}
\begin{table}
\begin{center}
    {\bf III. Quasi-Optimality~}\\
    {\footnotesize
\begin{tabular}{|p{1cm}| p{6cm} | p{2.3 cm}  |p{1.5cm}|}%\footnotesize%\label{tb1}
  %   \hline
 %Constant & Dependency on & Appears first in &  Value \\
 \hline 
 $h_4$&  $ \frac{1}{2(7+C_3 \varLambda_1+{\beta \varLambda_1})/({C_{\rm eff}^{-1}+\alpha})}$  & Lemma 4.3  &  Positive \\ \cline{1-3} \cline{4-4}
$ \theta_A$     & $\min\{1,\frac{C_{\rm eff}}{9C_4}\}$&   & (0, 1] \\ \cline{1-2} \cline{4-4}
$ \theta_B$     & $\frac{\gamma_0^2}{8(2+C_4C_1h^2)}$& Theorem  4.5 &  (0,1] \\  \cline{1-2} \cline{4-4} %$\frac{\gamma_0^2}{24} <\frac{\gamma_0^2}{8(2+C_4C_1h^2)}$
$\kappa$& $   \min \{ \kappa_0,  \frac{{C_{\rm eff}}/{8}-C_4\theta_A }{4C_4+ {\varLambda_3}/{2}}\big\}$&  ({\it Quasi-} &  Positive \\  \cline{1-2} \cline{4-4}
%$\tau_1$& $\min\{ \tau_2, \tau_3 \}$& {optimality}) &  Positive \\ \cline{1-2} \cline{4-4}
$\epsilon$& $\min\{ \epsilon, 1, \frac{1}{4C_4} \}$&{optimality}) &  Positive \\ \cline{1-2} \cline{4-4}
$h_5,h_7$& $ \epsilon $& &  Positive \\ \cline{1-2} \cline{4-4}
$\tau_2^2$& $ \frac{(c_{\rm eff}/8)-\kappa ( 4C_4+ {\varLambda_3}/{2})-C_5 \theta_A}{(5/4)C_7} $&  &  Positive \\ \cline{1-2} \cline{4-4}
$h_6$&  $\epsilon $ and $\min\{ h_2, h_4,\frac{1}{\sqrt{4C_4 \varLambda_1}} \}$&  &  Positive \\ \cline{1-2} \cline{4-4}
$\tau_3^2$& $\frac{\beta \gamma_0^2}{8(2+C_4C_1h^2)C_8} $ &  &  Positive \\
& $C_8=C_4C_1(C_{\rm rel +1})(1+\frac{1}{\kappa})$&&\\  \cline{1-2} \cline{4-4}

$h_8$&  $\epsilon $ and $\min\{ h_2,h_4, \sqrt{\frac{\alpha}{\beta C_1C_4 }},\sqrt{\frac{\gamma_0^2}{4C_8 }} \}$&  &  Positive \\ \cline{1-2} \cline{4-4}
$h_0$& $\min\{ h_7, h_8 \}$& &  Positive \\ \cline{1-2} \cline{4-4}
$\epsilon_1$& $\tau_1 \xi_\ell $, where $\tau_1:=\min\{ \tau_2, \tau_3 \}$&   &  Positive \\  \hline
  \end{tabular} }
%  \begin{tabular}{| p{1cm} | p{5.5cm} |p{2.0cm}  |p{1.5cm}|}%\footnotesize%\label{tb1}
%  %   \hline
% %Constant & Dependency on & Appears first in &  Value \\
% \hline \hline
%$ \epsilon$ &    & $\min\{1,\frac{C_{\rm eff}}{9C_4}\}$&   & [0, 1] \\ \cline{1-2} \cline{4-4}
%$ h $ &    & $\frac{\gamma_0^2}{24}<\frac{\gamma_0^2}{8(2+C_4C_1h^2)}$& Theorem  4.4 &  [0,1] \\  \cline{1-2} \cline{4-4}
%
%  \end{tabular}
\end{center}
\end{table}




%% end of file `docultex.tex'

1402.3068.extracted_bib
{GoMiSa {\sc M. Goossens, F. Mittelbach, and A. Samarin, %{\em The \LaTeX\ {\em Companion, Addison-Wesley, Reading, MA, 1994. % %
1402.3068.extracted_bib
{Higham {\sc N.~J. Higham, {\em Handbook of Writing for %the Mathematical Sciences, Society for Industrial and Applied %Mathematics, Philadelphia, PA, 1993. % %
1402.3068.extracted_bib
{Lamport {\sc L. Lamport, \LaTeX: {\em A Document %Preparation System, Addison-Wesley, Reading, MA, 1986. % %
1402.3068.extracted_bib
{SerLev {\sc R. Seroul and S. Levy, {\em A %Beginner's Book of \TeX, Springer-Verlag, Berlin, New %York, 1991.
1402.3068.extracted_bib
{oden {\sc M. Ainsworth and J. T. Oden, {\em A posteriori error estimation in finite element analysis, Wiley, New York, 2000. %
1402.3068.extracted_bib
{arnold % {\sc Arnold, N. Douglas, and R. S. Falk, % {\em A uniformly accurate finite element method for the {Reissner-{Mindlin plate, % SIAM J. Numer. Anal., {26 ({1989), pp. {1276-1290.
1402.3068.extracted_bib
{3matlab {\sc C. Bahriawati and C. Carstensen, {\em Three matlab implementation of the lowest order Raviart-Thomas MFEM with a posteriori error control, Comput. Methods Appl. Math., 5 (2005), pp. 333-1361. %
1402.3068.extracted_bib
{Braess % D. Braess, % {Finite elements- theory, fast solvers and applications in solid mechanics, % Cambridge university press, 1992.
1402.3068.extracted_bib
{Brenner-Scott {\sc S. C. Brenner and L. R. Scott, {\em The mathematical theory of finite element methods, third ed., Texts in Applied Mathematics, 15, Springer, New York, 2008. %
1402.3068.extracted_bib
{Brennerenrichment % S.~C. Brenner. % \newblock Two-level additive {Schwarz preconditioners for % nonconforming finite % element methods. % \newblock {\em Math. Comp., 65(215):897--921, 1996.
1402.3068.extracted_bib
{bek_mixed {\sc R. Becker and S. Mao, {\em An optimally convergent adaptive mixed finite element method, {Numer. Math., 111 (2008), pp. 35-54. %FJOURNAL = {Numerische Mathematik, % VOLUME = {111, % YEAR = {2008, % NUMBER = {1, % PAGES = {35--54,
1402.3068.extracted_bib
{bek_non {\sc R. Becker and S. Mao, {\em A convergent nonconforming adaptive finite element method with quasi-optimal complexity, {SIAM J. Numer. Anal., %FJOURNAL = {SIAM Journal on Numerical Analysis, {47 ({2010), pp. {4639-4659.
1402.3068.extracted_bib
{biven {\sc P. Binev, W. Dahmen, and R. DeVore, {\em Adaptive finite element methods with convergence rates, {Numer. Math., %FJOURNAL = {Numerische Mathematik, {97 ({2004), pp. {219-268.
1402.3068.extracted_bib
{Brezzi {\sc F. Brezzi and M. Fortin, {\em Mixed and hybrid finite element methods, Springer Verlag, 1991.
1402.3068.extracted_bib
{Car_apo {\sc C. Carstensen, {\em A~posteriori error estimate for the mixed finite element method, Math. Comp., 66 (1997), pp. 465-476.
1402.3068.extracted_bib
{Car1 {\sc C. Carstensen, {\em A unifying theory of a posteriori finite element error control, Numer. Math., 100 (2005), pp. 617-637.
1402.3068.extracted_bib
{Car_apo1 {\sc C. Carstensen, {\em Convergence of adaptive finite element methods in computational mechanics, Appl. Numer. Math., 59 (2009), pp. 2119-2130.
1402.3068.extracted_bib
{CH {\sc C. Carstensen and J. Hu, {\em A unifying theory of a~posteriori error control for nonconforming finite element methods, Numer. Math., 107 (2007), pp. 473-502. %
1402.3068.extracted_bib
{Carhop % C. Carstensen, M. Eigel, C. L{\"obhard and R.H.W. Hoppe, % {A review of unified a posteriori finite element error control, % Numer. Math. Theory Methods Appl., 5, 2012, 509-558. %
1402.3068.extracted_bib
{CC {\sc C. Carstensen, Asha. K. Dond, N. Nataraj, and Amiya. K. Pani, {\em Error analysis of nonconforming and mixed FEMs for second-order linear non-selfadjoint and indefinite elliptic problems, {\em http://arxiv.org/abs/1401.4810. % %
1402.3068.extracted_bib
{CHZ % P. Ciarlet, Jr., J. Huang and J. Zou, % {Some observations on generalized saddle-point problems, % SIAM J Matrix Anal. Appl. 25, 2003, 224-236. % % %
1402.3068.extracted_bib
{Doug % J. Douglas, Jr. and J. E. Roberts, % {Global estimates for mixed methods for second order elliptic equations, % Math. Comp. 44, 1985, 39-52. % %
1402.3068.extracted_bib
{dauge % M. Dauge, {Elliptic boundary problems on corner domains, % Lecture Notes in Math. Vol. 1341, Springer verlag, Berlin 1988. % %
1402.3068.extracted_bib
{evans % L.C. Evans, % {Partial Differential Equations, % Graduate Studies in Mathematics Vol. 19, AMS, Providence, Rhode Island, 2002. % % %
1402.3068.extracted_bib
{GT % D. Gilberg, N. S. Trudinger, % {Elliptic partial differential equations of second order, % Springer Verlag, 1982. % %
1402.3068.extracted_bib
{G % P. Grisvard, % {Elliptic Problems in Nonsmooth Domains, Pitman, Boston, 1985.
1402.3068.extracted_bib
{car_cgt2 {\sc C. Carstensen and R. Hella, {\em An optimal adaptive mixed finite element method, {Math. Comp., %FJOURNAL = {Mathematics of Computation, {80 ({2011), pp. {649-667.
1402.3068.extracted_bib
{car_non_cgt {\sc C. Carstensen and R.H.W. Hoppe, {\em Convergence analysis of adaptive nonconforming finite element methods, Numer. Math., 103 (2006), pp. 251-266.
1402.3068.extracted_bib
{car_cgt1 {\sc C. Carstensen and R. H. W. Hoppe, {\em Error reduction and convergence for an adaptive mixed finite element method, {Math. Comp., {75 ({2006), pp. {1033-1042.
1402.3068.extracted_bib
{axiom {\sc C. Carstensen, M. Feischl, M.Page, and D. Praetorius, {\em Axioms of Adaptivity, ASC Report 38/2013, Institute for Analysis and Scientific Computing, Vienna University of Technology.
1402.3068.extracted_bib
{quasi {\sc J. Cascon, C. Kreuzer, R. H. Nochetto, and K. G. Siebert, {\em Quasi-optimal convergence rate for an adaptive finite element method, {SIAM J. Numer. Anal., % FJOURNAL = {SIAM Journal on Numerical Analysis, {46 (2008), pp. {2524-2550.
1402.3068.extracted_bib
{cascon {\sc J. Cascon, R. H. Nochetto, and K. G. Siebert, {\em Design and convergence of AFEM in $H({\rm div)$, {Math. Models Methods Appl. Sci., % FJOURNAL = {SIAM Journal on Numerical Analysis, {17 ({2007), pp. {1849-1881.
1402.3068.extracted_bib
{Chen {\sc J. Chen and L. Li, {\em Convergence and domain decomposition algorithm for nonconforming and mixed methods for nonselfadjoint and indefinite problems, Comput. Methods Appl. Engrg., 173 (1999), pp. 1-20. %
1402.3068.extracted_bib
{chen_cgs {\sc L. Chen, M. Holst, and J. Xu, {\em Convergence and optimality of adaptive mixed finite element methods, Math. Comp., 78 (2009), pp. 35-53.
1402.3068.extracted_bib
{chen_hoppe {\sc H. Chen, X. Xu, and R.H.W. Hoppe, {\em Convergence and Quasi-optimality of adaptive nonconforming finite element methods for some nonsymmetric and indefinite problems, Numer. Math., 116 (2010), pp. 383-419.
1402.3068.extracted_bib
{l2 {\sc A. Demlow and R Stevenson, {\em Convergence and quasi-optimality of an adaptive finite element method for controlling {$L_2$ errors, {Numer. Math., %FJOURNAL = {Numerische Mathematik, {117 ({2011), pp. {185-218.
1402.3068.extracted_bib
{daglas {\sc J. Jr. Douglas and J.E. Robert, {\em Global estimates for mixed methods for second order elliptic equations Math. Comput. 44(1985), pp. 39-51. %
1402.3068.extracted_bib
{Du_xie % {\sc S. Du and X. Xie, % {\em Error reduction, convergence and optimality of an adaptive % mixed finite element method, % {J. Syst. Sci. Complex., % %FJOURNAL = {Journal of Systems Science \& Complexity, % {25 (2012), pp. {195-208. %
1402.3068.extracted_bib
{prac {\sc M. Feischl, T. F$\ddot{u$hrer, and D. Praetorius, {\em Adaptive FEM with optimal convergence rates for a certain class of non-symmetric and possibly non-linear problems, ASC Report 43/2012, Institute for Analysis and Scientific Computing, Vienna University of Technology.
1402.3068.extracted_bib
{nochetto {\sc K. Mekchay and R.H. Nochetto, {\em Convergence of adaptive finite element methods for general second order elliptic PDE, SIAM J. Numer. Anal., 43 (2005), pp. 1803-1827.
1402.3068.extracted_bib
{Marini {\sc L. D. Marini, {\em An inexpensive method for the evaluation of the solution of the lowest order {Raviart-{Thomas mixed method, {SIAM J. Numer. Anal., %FJOURNAL = {SIAM Journal on Numerical Analysis, {22 ({1985), pp. {493-496.
1402.3068.extracted_bib
{dof {\sc W. D$\ddot{o$rfler, {\em A convergent adaptive algorithm for Poisson's equation, SIAM J. Numer. Anal. 33 (1996), pp. 1106-1124. %
1402.3068.extracted_bib
{mao % {\sc S. Mao, X. Zhao, and Z. Shi, % {\em Convergence of a standard adaptive nonconforming finite % element method with optimal complexity, % {Appl. Numer. Math., % FJOURNAL = {Applied Numerical Mathematics. An IMACS Journal, %{60 ({2010), pp. {673-688.
1402.3068.extracted_bib
{morin1 {\sc P. Morin, R.H. Nochetto, and K.G. Siebert, {\em Data oscillation and convergence of adaptive fem, SIAM J. Numer. Anal., 38 (2000), pp. 466-488.
1402.3068.extracted_bib
{morin2 {\sc P. Morin, R.H. Nochetto, and K.G. Siebert, {\em Convergence of adaptive finite element methods, SIAM Rev., 44 (2002), pp. 631-658.
1402.3068.extracted_bib
{veeser {\sc P. Morin, K. G. Siebert, and A. Veeser, {\em A basic convergence result for conforming adaptive finite elements, {Math. Models Methods Appl. Sci., % FJOURNAL = {Mathematical Models \& Methods in Applied Sciences, {18 ({2008), pp. {707-737.
1402.3068.extracted_bib
{schatz1 {\sc A. H. Schatz, {\em An observation concerning Ritz-Galerkin methods with indefinite bilinear forms, Math. Comput., 28 (1974), pp. 959-962.
1402.3068.extracted_bib
{schatz {\sc A. H. Schatz and J. Wang, {\em Some new error estimates for Ritz-Galerkin methods with minimal regularity assumptions, Math. Comp., 65 (1996), pp. 19-27. %
1402.3068.extracted_bib
{scottzhang %{\sc L. R. Scott and S.Zhang, {\em Finite element interpolation of nonsmooth function satisfying boundary conditions, Math. Comp., 54 (1990), pp. 483-493.
1402.3068.extracted_bib
{steven {\sc R. Stevenson, {\em Optimality of a standard adaptive finite element method, Found. Comput. Math., 7 (2007), pp.245-269.
1402.3068.extracted_bib
{verf {\sc R. Verf\"{urth, {\em A posteriori error estimation and adaptive mesh-refinement techniques, J. Comput. Appl. Math., 50 (1994), pp. 67-83.