of the solution operator $S$. Then we have the following well-known conclusion about error estimates.

\[lem:error estimates\][([@Ciarlet2002The], Thm. 4.4.6)]{} For a given $u\in L^{2}(\Omega)$, let $y$ be the unique weak solution of the state equation (\[eqn:pde\]) and $y_{h}$ be the unique solution of (\[eqn:weak formulation\]). Then there exists a constant $c>0$ independent of $h$, $u$ and $y_{r}$ such that $$\|y-y_{h}\|_{L^{2}(\Omega)}+h\|\nabla y-\nabla y_{h}\|_{L^{2}(\Omega)}\leq ch^{2}(\|u\|_{L^{2}(\Omega)}+\|y_{r}\|_{L^{2}(\Omega)}).$$ In particular, this implies $\|S-S_{h}\|_{L^{2}(\Omega)\rightarrow L^{2}(\Omega)}\leq ch^{2}$ and $\|S-S_{h}\|_{L^{2}(\Omega)\rightarrow H^{1}(\Omega)}\leq ch$.

For the given triangulation $\mathcal{T}_{h}$ with nodes $\{x_{i}\}_{i=1}^{N_{h}}$, let $\{\phi_{i}(x)\}_{i=1}^{N_{h}}$ be a set of nodal basis functions associated with nodes $\{x_{i}\}_{i=1}^{N_{h}}$ which spans $Y_{h}$, $U_{h}$ and satisfies the properties: $$\phi_{i}(x)\geqslant0,\ \ \|\phi_{i}(x)\|_{\infty}=1\ \ \ \forall{i=1,...,N_{h}}, \ \ \sum\limits_{i=1}^{N_{h}}\phi_{i}(x)=1.$$ Then $y_{h}\in Y_{h}$, $u_{h}\in U_{h}$ can be represented as $y_{h}=\sum\limits_{i=1}^{N_{h}}y_{i}\phi_{i},$ $u_{h}=\sum\limits_{i=1}^{N_{h}}u_{i}\phi_{i}$ and we have $y_{h}(x_{i})=y_{i}$, $u_{h}(x_{i})=u_{i}$. The other variables and operators in the subproblems of the inexact ADMM in function space are all discretized by piecewise linear functions similarly.

Before we give the proposed algorithm to solve the problem (\[eqn:reduced problem with linear constraint\]), we first introduce the definition of node interpolation operator $I_{h}$ for the following convergence analysis.

For a given regular and quasi-uniform triangulation $\mathcal{T}_{h}$ of $\Omega$ with nodes $\{x_{i}\}_{i=1}^{N_{h}}$, let $\{\phi_{i}(x)\}_{i=1}^{N_{h}}$ denotes a set of 