## 2 The Mathematical Model
: _We assume that for some known constant \(d>1\), the function \(\varepsilon\in C^{2}\left(\mathbb{R}^{3}\right)\) satisfies  
Figure 1: Domain decomposition and mesh discretization in \(\Omega\). The mesh of \(\Omega\) in (a) is a combination of the quadrilateral finite difference mesh \(\Omega_{FDM}\) presented in (**c**) and the finite element mesh \(\Omega_{FEM}\) presented in (**b**). Domains \(\Omega_{FEM}\) and \(\Omega_{FDM}\) overlap by two layers of structured nodes such that they have common boundaries shown by green and blue. \[\varepsilon(x)\in[1,d],\quad\text{ in }\Omega\setminus\Omega_{FDM} \tag{3}\] \[\varepsilon(x)=1,\quad\quad\text{ in }\Omega_{FDM}.\]  
As an immediate consequence of this assumption we have that  
\[\varepsilon\equiv 1\quad\text{in }\Omega_{FDM}\cap\Omega_{FEM}\quad\text{ also }\quad \varepsilon\equiv 1\quad\text{in }\partial\Omega_{FDM}. \tag{4}\]  
Assumption 1 has a physical meaning. The model problem (2) has a lot of applications in the field of coefficient inverse problems what we have pointed out in the introduction. One of the application of the model problem (2) is determination of dielectric permittivity function from measurements of the electric field at the boundary of the investigated domain, and this is a typical CIP. We refer for applications of such CIPs in [9; 16; 17; 18; 19] and references therein. This CIP is an ill-posed problem and needs several assumptions in order to be able to solve it, see details in [20; 21; 22]. One of the assumptions is the introduction of the set of admissible parameters for the dielectric permittivity function, which we define in the Condition (3). The lower bound in this set corresponds to the homogeneous media and the upper bound \(d\) is chosen as maximal values of the dielectric permittivity function and depends on the concrete physical problem. We refer to [16; 18; 23; 24] where one can find values of dielectric permittivity function for different real-life applications with applications in medical imaging. Using these values can determine the upper and lower bounds in the set of admissible parameters for this function. See works [9; 17; 19] for examples of such sets for some real physical problems. Note further that \(\partial\Omega_{FEM}\subset\Omega_{FDM}\). Condition (3) on \(\varepsilon\) and the relation  
\[\
abla\times\
abla\times E=\
abla(\
abla\cdot E)-\
abla\cdot(\
abla E), \tag{5}\]  
together with divergence-free field \(E\), make the equations in (2) independent of each other in \(\Omega_{FEM}\) and \(\Omega_{FDM}\) so that, in \(\Omega_{FDM}\) (with \(\varepsilon=1\)), we just need to solve the system of wave equations:  
\[\frac{\partial^{2}E}{\partial t^{2}}-\Delta E=0,\quad(x,t)\in\Omega_{FDM}\times (0,T]. \tag{6}\]  
**Remark 1**. : _There are different discretizations of [PERSON]'s equations in the time domain such as the FDTD Yee scheme [8; 25], FE Lee-Madsen scheme [26], node-based least squares FE schemes [27; 28], edge elements [29; 30], and a number of discontinuous Galerkin FE schemes [31; 32; 33; 34; 35]. It is well known that, for stable implementation of the finite element solution of [PERSON]'s equations, divergence-free edge elements are the most satisfactory ones from a theoretical point of view [29; 30; 36]. However, the edge elements are less attractive for the solution of time-dependent problems, since a linear system of equations should be solved at each time iteration step. On the contrary, \(P1\) elements can be efficiently used in a fully explicit finite element scheme with lumped mass matrix, see, e.g., [37; 38; 39]. It is also well known that the numerical solution of [PERSON]'s equations using nodal finite elements is often unstable and results in spurious oscillatory solutions [40; 41]. There are a number of techniques to overcome such instabilities, see, e.g., [27; 28; 41; 42; 43]._  
In [44], a finite element analysis shows stability and consistency of the stabilized finite element method for the solution of (1) with \(\sigma(x)=0\).