Reconstruction of quasi-local numerical effective models from low-resolution measurements

We consider the inverse problem of reconstructing an effective model for a prototypical diffusion process in strongly heterogeneous media based on coarse measurements. The approach is motivated by quasi-local numerical effective forward models that are provably reliable beyond periodicity assumptions and scale separation. The goal of this work is to show that an identification of the matrix representation related to these effective models is possible. On the one hand, this provides a reasonable surrogate in cases where a direct reconstruction is unfeasible due to a mismatch between the coarse data scale and the microscopic quantities to be reconstructed. On the other hand, the approach allows us to investigate the requirement for a certain non-locality in the context of numerical homogenization. Algorithmic aspects of the inversion procedure and its performance are illustrated in a series of numerical experiments.


Introduction
This paper focuses on the computational solution of multiscale inverse problems, i.e., where the quantities to be sought are related to an unknown microscopic mathematical model, while measurement data are available only on a much coarser scale. Due to the scale mismatch between given data and unknowns, the direct recovery of microscopic quantities, e.g., in the form of a coefficient of a partial differential equation (PDE), is not only expensive but may also lead to unsatisfactory results. Therefore, the aim of this paper is to introduce a computational framework -inspired by numerical homogenization methods -to reconstruct an effective model, i.e., an alternative quantity on the coarse scale of the available data.
Effective models are the key to bridge the discrepancy between a microscopic coefficient and a coarse scale of interest in the forward setting. They provide models that compute reliable approximations of the solution of a PDE even in the presence of microscopic quantities. If structural assumptions such as (local) periodicity or scale separation hold, classical homogenization methods (see, e.g., [HW97, MS02, EE + 03, EE05]) based on analytical homogenization theory can be used. These methods are local in the sense that the communication among the degrees of freedom is only between neighbors. In a more general setting where these structural assumptions do not hold a priori, numerical homogenization methods (see, e.g., [MP14,HP13,GGS12,OZB14,EGH13,Owh17,Mai20]) provably provide an alternative. These methods are based on a coarse mesh with a characteristic mesh parameter and compute special problem-adapted basis functions with optimal approximation properties. Compared to the locality of classical homogenization methods, numerical homogenization methods typically involve a slight deviation from local communication between the degrees of freedom which, in turn, leads to somewhat increased sparsity patterns of the corresponding system matrices. Since this non-locality can be controlled, we refer to these methods as quasi-local.
This paper follows the pragmatic approach of reconstructing quasi-local effective models (i.e., their representation in term of quasi-local system matrices) that describe the effective behavior of a medium with microstructures based on coarse (i.e., low-resolution) measurement. On the one hand, this provides a surrogate to a direct reconstruction in the above-mentioned multiscale context. On the other hand, the reconstruction of an effective model allows us to investigate the requirement for a certain quasi-locality in the context of numerical homogenization approaches.
The goal of this work is to promote this idea along with algorithmic aspects and preparatory numerical experiments. To demonstrate the feasibility and the potential of the approach, we investigate a stationary linear elliptic multiscale diffusion problem. Moreover, we consider a worst-case scenario without any structural a priori knowledge on the underlying diffusion coefficient and do not assume that the heterogeneous coefficient can be parameterized by a few unknown parameters that could more easily be identified. The main novelty of our approach is that it aims to recover information about the microscopic scale in the sense of reproducing the effective behavior of corresponding solutions, instead of aiming at identifying the actual microscopic coefficient. In the multiscale setting, in which measurement are given on a coarse scale and without a priori assumption on the structure of the microscopic coefficient, one cannot hope to reasonably identify the actual coefficient and our approach presents an alternative strategy. Our aim is not to compete with classical inverse strategies but to provide a first step towards reasonable surrogates in the case where a direct recovery of the coefficient is either too expensive or not reliable due to multiscale aspects of the problem.
The remaining parts of the paper are organized as follows. We start with introducing the microscopic forward problem and motivating the reconstruction of an effective model, represented by an effective system matrix (Section 2). This strategy is inspired by numerical homogenization strategies which provably provide reliable effective models for the given forward model. The rigorous motivation is presented in Section 4 and is mainly to emphasize that an effective model indeed exists in the setting of very general coefficients. Based on the considerations for the forward problem, we then tackle the reconstruction of an effective quasi-local model from given measurements. To this end, we prescribe a quasi-local sparsity pattern of the system matrices and rephrase the inverse problem as a non-linear least squares problem for which we apply iterative minimization techniques such as the gradient descent or the Gauß-Newton method (Section 3). In a series of numerical experiments (Section 5), we show that quasi-local effective models can indeed be reconstructed. In particular, we consider the cases where we are given measurements for all possible (coarse) boundary conditions, and also the setting where solutions are only known for a few boundary conditions. The aim of the experiments is to show that allowing the model to deviate from locality improves the inversion process and, thus, justifies the previous discussion.

Microscopic Forward Problem
In this section, we present the forward model and identify an effective discrete model, which is characterized by an appropriate system matrix.
2.1. Problem setting. We consider the prototypical second-order linear elliptic diffusion problem where Ω ⊂ R d , d ∈ {1, 2, 3} is a polyhedral domain and the diffusion coefficient A encodes the microstructure of the medium. We do not make any structural assumptions on the coefficient such as periodicity or scale separation. Admissible coefficients are elements of the following set, which only requires minimal assumptions.
Since solutions to problem (2.1) do not necessarily exist in the classical sense, we are interested in the weak solution of (2.1) in the Sobolev space V := H 1 (Ω) which is characterized by the following variational formulation. Given A ∈ A, u 0 ∈ X := H 1/2 (∂Ω), and f ∈ L 2 (Ω), we seek u ∈ V such that where tr : V → X is the trace operator, (w, v) := (w, v) L 2 (Ω) denotes the L 2 inner product, and a(w, v) := Ω A∇w · ∇v dx. Note that instead of (2.1), we could as well consider a general second-order linear PDE in divergence form with additional lower-order terms. Such a generalization is straight-forward but is omitted for simplicity.
2.2. Coarse discretization. In practice, it is favorable to rewrite problem (2.2) as a problem with homogeneous Dirichlet boundary conditions in V 0 . Let E b : X → V be a linear extension operator, which also defines the restriction operator R : V → V 0 by R := 1 − E b tr. Then, we can decompose u = Ru + (1 − R)u = Ru + E b u 0 and problem (2.2) reduces to finding Ru ∈ V 0 such that Let us now introduce a coarse target scale H (e.g., the resolution of the data available for the inverse problem). We adopt the notation from numerical homogenization where a capital H is used to indicate that the scale is indeed a coarse one. In typical applications, H will be much larger than the microscopic scale, i.e., the scale on which the diffusion coefficient varies.
In order to discretize (2.3), let T H be a mesh of orthotopes with characteristic mesh size H and denote with Q 1 (T H ) the corresponding space of piecewise bilinear functions. Further, we define the discrete spaces V H := Q 1 (T H )∩V , V 0 H := V H ∩V 0 , and X H := tr V H with dimensions m = dim V H , m 0 = dim V 0 H , and n = dim X H , respectively. The choice of these finite element spaces is not unique and other standard finite element spaces could be used.
In the context of inverse problems, it is reasonable to consider that u 0 is defined as the first order finite element approximation of coarse experimental boundary data which approximate the real data up to order H in the H 1/2 norm. That is, in the following we will assume that u 0 ∈ X H . Further, we assume to have a discretized extension oper- H and a corresponding restriction operator H tr. Based on the above spaces, we introduce an injective linear operator and its discretized version The operator L A (and similarly also L G A, ) can be written as For simplicity, we assume in the following that f is a fixed function and only the boundary conditions may change. The generalization to the case where f is variable as well is conceptually straightforward but slightly more involved. The decomposition (2.8) motivates the distance function between operators defined by for all A, B : X H × L 2 (Ω) → V . Note, however, that also other distance functions could be used.
2.3. Characterization in terms of a stiffness matrix. As a next step, we discuss an alternative representation of the operator L G A, using the stiffness matrix corresponding to the discrete formulation (2.5). Given a coefficient A ∈ A and a mapping G : V H → V as above, the stiffness matrix S H = S H (A, G) is defined by where i → z i is a fixed ordering of the m nodes in T H and Λ z denotes the classical finite element hat function associated with the node z ∈ T H . Therefore, we may define the operator (2.11) with the classical finite element mass matrix M H , the restriction S H,0 = R H S H R T H of S H to the inner nodes of T H , and f H := Π H f the L 2 projection of f onto V H . For better readability, we use the notation v H (or B H ) for both the vector v H ∈ R m (or the matrix B H ∈ R m 0 ×m ) and the corresponding function v H ∈ V H (or the mapping B H : V H → V 0 H ). The following theorem is the basis of the inverse strategy discussed in the next section. It is a direct consequence of Lemma 4.3 and Corollary 4.2, which are proven in Section 4. It states the existence of an appropriate coarse model, characterized by an operator G, that is able to capture the effective behavior of the original model (2.6).
Theorem 2.1 (Existence of an effective model). There exists an operator G as in (2.4) and a corresponding stiffness matrix S H (A, G) such that Moreover, there exists a choice of G and ∼ | log H| such that S H (A, G) ∈ M( , T H ) with defined as the set of matrices that may have a non-zero entry at position [i, j] only if the corresponding nodes z i and z j belong to the -neighborhood of each other. Theneighborhood N is defined by We call the operator L S H (A,G) the effective model and S H (A, G) the effective stiffness matrix.
Theorem 2.1 justifies the inverse procedure that is presented in the following section. In particular, we may interpret coarse measurements of a solution operator L A as approximations obtained by an effective model due to the fact that there exists an appropriate effective model which is reasonably close.
2.4. Quasi-locality and connection to numerical homogenization. The operator G defined in (2.4) and the corresponding effective stiffness matrix defined in (2.10) are strongly related to numerical homogenization methods. In contrast to analytical homogenization, these approaches provably work beyond structural assumptions such as (local) periodicity or a clear separation of scales which cannot be guaranteed for general microstructures.
One such method is the Localized Orthogonal Decomposition (LOD) approach which provides effective models that provably cope with arbitrary rough coefficients in a large class of model problems including diffusion problems [MP14,HP13,HMP14], elasticity [HP16, ACM + 20] and wave propagation [GP15, AH17, Pet17, Ver17, GHV18, MP19], without requiring periodicity or scale separation. This method allows us to explicitly characterize an operator G to prove Theorem 2.1 in Section 4. For linear elliptic problems, there are various other numerical homogenization approaches such as the Generalized Finite Element Methods (GFEM) [BL11], AL bases [GGS12], Rough Polyharmonic Splines (RPS) [OZB14], the Generalized Multiscale Finite Element Method (GMsFEM) [EGH13], Gamblets [Owh17], CEM-GMsFEM [CEL18], the higher-order multiscale approach described in [Mai20, Ch. 3], and their variants with similar properties as LOD. All these methods compute special problem-adapted basis functions with optimal approximation properties based on underlying Galerkin methods. To achieve optimal accuracy, a moderate price in terms of the computational complexity has to be paid compared to a standard finite element method (fixed order) on the same mesh in order to account for microscopic information. The computational overhead is either characterized by an increase in the number of degrees of freedom per mesh entity (GFEM, GMsFEM), e.g., elements or nodes, or in enlarging the support of the discrete basis functions (LOD, RPS, Gamblets, AL bases). In both cases, the result is a slightly denser sparsity pattern of the corresponding system matrices which is due an increased communication between the degrees of freedom. This quasi-locality distinguishes the above numerical homogenization methods from classical numerical multiscale methods based on homogenization theory such as the Multiscale Finite Element Method (MsFEM) [HW97], the Two-Scale Finite Element Method [MS02], or the Heterogeneous Multiscale Method (HMM) [EE + 03, EE05] that share the communication pattern of standard finite element methods. Quasi-locality also showed to be viable in connection with the pollution effect in highfrequency time-harmonic wave propagation [BS97] which cannot be avoided unless the mesh size is coupled to, e.g., the polynomial degree [MS10,MS11,MPS13] or the support of the basis functions [Pet17] in a logarithmic way. Promising results using non-local models have also been achieved in the field of peridynamics [Sil00, Lip14,Du17] or in isogeometric analysis [HCB05,CHB09]. Moreover, quasi-locality also includes the increased communication of higher-order finite element approaches.

Inverse Problem: Reconstruction of an Effective Model
Based on the considerations in the forward setting, we are now able to formulate the inverse problem which reconstructs an effective model characterized by a stiffness matrix with a certain sparsity behavior.
3.1. Problem setting. Let us assume that the diffusion coefficient A is unknown and that structural assumptions such as periodicity, quasi-periodicity, and given parameterization by few degrees of freedom are not satisfied a priori. In an ideal setting, information about solutions to problem (2.2) in the form of a solution operator would be given. In practical applications, however, boundary data and information about the corresponding solutions are only available on some (coarse) scale, possibly much larger than the (micro) scale on which the diffusion coefficient and the corresponding solutions vary. A classical formulation of the inverse problem, for a fixed right-hand side f , consists in recovering A in (2.2) given the mapping However, since the unknown coefficient A includes microscopic features, a direct approach of recovering A by simulations on the microscopic scale is computationally unfeasible and does not provide reasonable reconstructions due to the mismatch between coarse data and a possibly microscopic coefficient. In this section, we present an alternative approach to recover information about the (macroscopic) effective model taking into account the presence of a micro-scale diffusion coefficient. Rather than reconstructing the diffusion coefficient itself, we tackle the reconstruction of an effective stiffness matrix that is able to reproduce the given data related to solutions to (2.2). Therefore, the alternative formulation of the inverse problem reads: givenL eff : X H → V H , find the corresponding stiffness matrixS H .
This alternative inverse problem is built upon the assumption that coarse given measurements on the scale H may be represented by an effective mapping (up to order H) based on Theorem 2.1. The effective model corresponding to the reconstructed matrix provides a reliable surrogate that contains effective properties and that may be used for simulations on the coarse data scale. Moreover, knowledge on numerical homogenization methods can, in a second step, provide further information about, e.g., geometric features.
3.2. Medical background. The considered inverse problems are of internal type, i.e., we assume that measurement are available on the bulk domain. This setting is motivated by the use of medical imaging protocols based on Magnetic Resonance Imaging (MRI), which play a key role in modern diagnostics. Using strong magnetic fields, magnetic field gradients, and radio waves, MRI allows clinicians to obtain, in vivo and non-invasively (and without exposing the body to radiation), both geometrical features of the body, e.g., shapes of the organs, and functional properties, e.g., motion or diffusion processes. As an example, Magnetic Resonance Elastography (MRE) is an imaging technique which is sensitive to mechanical parameters of the tissue. During an MRE examination, the tissue undergoes an external mechanical excitation, imposed at given frequencies by so-called actuators, attached to the surface of human tissues. In parallel, using phase-contrast MRI, i.e., postprocessing the phase of the MRI signal, it is then possible to measure the internal displacement of the tissue and hence to recover how shear and compression waves propagate into the body [ME96, SBH + 08, HSB17]. In practice, MRE allows to obtain average displacement fields on each element of a three-dimensional Cartesian mesh (a voxel ), whose resolution is typically of the order of millimeters, and it is practically limited by the examination time and by the properties of the MRI scanner.
Combining these data with a mechanical tissue model, it is then possible to recover information about the elastic behavior of the tissue (a so-called elastogram). Clinical application of MRE are based on reconstructing tissue properties only on an effective scale, i.e., describing the living tissue as a (visco-)elastic material with mechanical parameters varying only on the coarse data scale. The procedure is currently used to diagnose and monitor tissue diseases such as cancer and fibrosis, that are characterized by different tissue stiffnesses. However, parameters describing the microscopic scales -such as tissue porosity and vascular structures -might play an important role in several applications, especially in the context of poroelastic tissues [HBG + 13, HGR + 14]. In those cases, in order to characterize the microstructural properties from coarse data, mathematical models defined on the effective scale -such as the one proposed in Section 3.1 -might provide a valuable alternative to efficiently take smaller scales into account.
. . , f H ] ∈ R m×n and the identity matrix 1 ∈ R m×m . The matrix L eff S H comprises full information about the forward problem in the sense that it includes the solutions of (2.11) for a complete set of basis functions of X H . The operatorL eff may also be interpreted as a matrix, so that the distance between the operators can be measured in general matrix norms. This is especially useful since a splitting of the form (2.8) is generally not known forL eff . Let µ := dim M( , T H ). Instead of (3.4), based on the matrix representation introduced above we consider a minimization problem for the functional J H : R µ → R defined by At this stage, the choice of the norm in R m×n in (3.6) is arbitrary. The results that we will show in Section 5 have been obtained using the Frobenius norm, which seems to be a natural candidate.
3.4. Iterative minimization. In order to find a minimizer of (3.6), we can now apply standard minimization techniques such as the Newton method or the gradient descent method. Here, we adopt a Gauß-Newton method which, in our numerical computations, showed faster convergence in terms of number of iterations.
In order to compute the descent direction, the most important step concerns the computation of the gradient of J H with respect to the relevant entries {s i } µ i=1 of S H (i.e., the diagonal and the non-zero entries above the diagonal, due to symmetry). Using the chain rule, we obtain For the Gauß-Newton method, only the derivatives of L eff S H are needed, i.e., Here ∂s i are relatively easy to compute, as they are defined as global matrices that only contain at most two entries equal to 1.
For ease of notation, let us interpret L eff S H and S H as vectors in R mn and R m 2 , respectively. The Gauß-Newton method to minimize the functional J H is then defined by the following steps.
• Let an initial matrix S 0 H ∈ M( , T H ) be given. • For k = 0, 1, . . . (until a certain stopping criterion is satisfied), solve where D denotes the derivative with respect to the relevant entries of S H and • Set P k ∈ M( , T H ) as the matrix whose relevant entries are given by p k and define (3.9) S k+1 H = S k H + δ k P k with appropriately chosen step size δ k , for example using backtracking line search based on the Armijo-Goldstein condition. Due to the ill-posedness of the inverse problem, the matrix H k might be singular. A possible approach to overcome this issue consists in replacing (3.8) by with a given parameter η > 0. Another possible strategy consists in adding a regularization term to the functional to be minimized, i.e., in replacing (3.6) by where γ > 0 is a given regularization parameter and S reg is a regularization (or stabilization) matrix. Additionally, the computations of the gradient in (3.7) need to be adapted accordingly. In the presence of multiple minimizers, this regularization enforces the solution to be close (depending on the parameter γ) to the matrix S reg . For example, if the aim of the inverse problem is to find defects in an otherwise homogeneous medium, a suitable choice for S reg could be a standard finite element stiffness matrix for a constant diffusion coefficient. In our practical computations, the regularization approach described in (3.10) is better suited since an appropriate regularization matrix S reg is generally not known. We emphasize that the presented inversion process does not need to resolve any microscopic scales in order to obtain an effective numerical model. The information extracted by this procedure (i.e., the stiffness matrixS H ) may be used to simulate other problems subject to the same (unknown) diffusion coefficient. Furthermore, the information gathered can be seen as an intermediate step towards recovering information concerning the original coefficient itself. This additional recovery step will be studied in a future work.

Justification of the Proposed Inverse Problem
This section is devoted to proving Theorem 2.1 in Section 2. In particular, we explicitly construct an operator G as introduced in (2.4). We emphasize, however, that the explicit construction is only to justify the reconstruction of an effective model. Apart from its communication behavior, none of the specific properties used below are required for the inverse procedure.
4.1. Effective forward approximation beyond structural assumptions. In this subsection, we use the multiscale technique known as Localized Orthogonal Decomposition (LOD) [MP14,HP13] to obtain a coarse forward model on the scale H which produces a forward operator that can be used as a replacement for L A as given in (2.6).
To this end, we discretize (2.3) in a suitable coarse multiscale space. Since the standard space V H is not suitable for the approximation of u if H is larger than the spatial scale of the microstructure, we enrich the coarse model with microscopic information about the problem via corrections of classical finite element functions. The construction of these corrections is based on a projective quasi-interpolation operator I H : V 0 → V 0 H with standard approximation and stability properties, i.e., for an element T ∈ T H with diameter H T , it holds that for all v ∈ V 0 , where the constant C is independent of H, and N(T ) := N 1 (T ) is the neighborhood of T as defined in (2.14). Note that for shape-regular meshes the above estimate also holds globally. For a particular choice of I H , see [EG17,GP17,KPY18].
Based on I H , we define, for any element T ∈ T H and any function v H ∈ V 0 H , the element correction C T v H ∈ W := KerI H by for all w ∈ W , and the full correction C : By construction, it holds that The corrections C T v H have, in general, global support. However, as shown in [HP13,Pet16] (based on [MP14]) they decay exponentially fast (see also the one-dimensional sketch in Figure 4.1). Therefore, we use localized element corrections C T, v H which are obtained by solving (4.2) on local patches N (T ) as defined in (2.14), i.e., for all w ∈ W with w| Ω\N (T ) = 0. As above, we define the full correction C : V 0 H → W by Hat function Λ Correction of hat function CΛ  As shown in [HP13], we get, for any v H ∈ V 0 H , (4.5) The constant c only depends on the contrast β/α, although this dependence seems pessimistic in many cases of practical relevance [PS16,HM17]. For v H ∈ V H , we set Cv H := CRv H and C v H := C Rv H . Given a discretized extension operator E b H : X H → V H and the corresponding restriction operator R H : V H → V 0 H as defined in Section 2.2, the discretized version of (2.3) reads: H . In (4.6), f H := Π H f is the L 2 projection of f onto V H , and u 0 ∈ X H .

Error estimates.
The following theorem shows that the approximation error of the presented approach scales optimally with H and that it is independent of the variations of the diffusion coefficient.
Theorem 4.1 (Error of the forward effective model). Let u ∈ V be the solution of (2.2) and u H ∈ V H the solution of (4.6), for given boundary data u 0 ∈ X H , a right-hand side f ∈ L 2 (Ω), as well as an oversampling parameter . For g ∈ L 2 (Ω), letû(g) ∈ V denote the solution of (2.2) with right-hand side g and boundary condition u 0 = 0, and let us introduce the worst-case best-approximation error It holds by the orthogonality property (4.3). For the next steps, we require the identity where we use (4.8) and (4.1). As a next step, we boundē H :=ũ H −ū H . We note that With v H =ē H and similar arguments as above, we obtain (4.10) ē H L 2 (Ω) The error u −ū H can be estimated using [GP17, Prop. 1] which also holds for inhomogeneous Dirichlet boundary conditions, i.e., (4.11) u −ū H L 2 (Ω) The triangle inequality, (4.9), (4.10), and (4.11) yield the desired estimate.  To illustrate the advantage of the LOD, Figure 4.2 shows the error between the numerical solution on a microscopic scale and the numerical solutions using the LOD and a classical finite element approximation on a coarse scale, respectively. The finite element method suffers from pre-asymptotic effects when the micro scale is not resolved, while the LOD produces a finite element function with much better approximation properties.
We emphasize that, choosing large enough (i.e., | log H|), it holds e −c H or even e −c H 2 . As discussed in [GP17], the worst-case best-approximation error is at least O(H), and it scales possibly even better with H for certain pre-asymptotic regimes. In this work, we are mainly interested in solving the inverse problem and do not focus on optimizing the error estimates derived above.
We can now define the discretized operator corresponding to (4.6) by where u H solves (4.6).
Using Theorem 4.1 and the distance function defined in (2.9), we obtain the following result.
Corollary 4.2 (Error of the effective forward operator). Let | log H|. Then it holds Reformulation using the effective stiffness matrix. As described in Section 2.3, we can alternatively represent the operator L eff A, using the stiffness matrix corresponding to the discrete formulation (4.6). Given a coefficient A ∈ A, the corresponding LOD stiffness matrix S H (A) is defined by with the ordering i → z i as in (2.10). Further, the set of LOD stiffness matrices with oversampling parameter based on admissible coefficients is given by We can now prove the following lemma using the operator defined in (2.11).
Lemma 4.3 (Alternative representation of the effective forward operator). Let S H (A) ∈ S( , T H ) be the LOD stiffness matrix corresponding to (4.6). Assume that Proof of Lemma 4.3. Write u H = m j=1 u j Λ z j and observe that (4.6) is equivalent to H , and the definition (4.13), we can write equation (4.16) as H u 0 , which shows (4.15).

Numerical Experiments
In this section, we present some numerical experiments that illustrate the capability of the proposed method. The inverse problem is based on synthetic data, i.e., the coarse measurements used to feed the inversion algorithm are obtained from finite element functions in V h , defined on a mesh with mesh size h = √ 2 · 2 −9 , that resolve the micro-scale features of the diffusion coefficient. Furthermore, the data are perturbed by multiplicative random noise with intensity up to 5%. 5.1. Example 1: full boundary data. In the first experiment, we assume to have full information on the operator (matrix)L eff , i.e., we assume that measurements in Ω on the coarse scale H = √ 2 · 2 −5 for a complete basis of X H are available. The scalar coefficient A for which the effective behavior should be recovered is constant on a mesh T ε with ε = √ 2 · 2 −7 and the value on each element is independently obtained as a uniformly distributed random number between 1 and 50, i.e., for any T ∈ T ε we have A| T ∼ U (1, 50) (see Figure 5.1, left, for the explicit sample used here). We set f = 1 and start the inverse iteration with the finite element stiffness matrix S 0 H based on the constant coefficient with value 1. The values of the functional J H in the first 20 iterations of the inversion algorithm are given in Figure 5.1 (right). In particular, we compare the performance of a local approach based on matrices in M(0, T H ) with the sparsity pattern of, e.g., a standard first-order finite element method, the HMM, or the Two-Scale Finite Element Method, with a quasi-local approach based on matrices in M( , T H ) for ∈ {1, 2, 3}. One clearly sees that slightly deviating from locality leads to better results in terms of decrease and value of the error functional J H . In particular, for = 0 the functional seems to reach a stagnation relatively quickly, while the results significantly improve when increasing the value of .
A necessary validation step, in order to further investigate the behavior dependent on , consists in solving a diffusion problem using the stiffness matrices reconstructed with using different sparsity patterns, and comparing the resulting numerical solutions with the finite element functions from which the measurements were taken to feed the inversion algorithm. The outcome of this assessment is shown in Figure 5.2, focusing on the cross sections at x 2 = 0.5 (left) and at x 1 = 0.5 (right) of the numerical approximations corresponding to the boundary condition u 0 (x 1 , x 2 ) = x 1 . Figure 5.3 depicts the same cross sections when a random boundary condition u 0 ∈ X H is considered.
Besides the accuracy of the numerical approximations computed based on the recovered stiffness matrices, it is also important to assess the robustness of the reconstructed effective model, i.e., to investigate to which extent the coarsened information about the diffusion coefficient encoded in the stiffness matrix can be used to simulate other scenarios.
For this purpose, we employ the reconstructed stiffness matrices to simulate a diffusion problem with two different right-hand sides, i.e., and g 2 (x 1 , x 2 ) = 10 1 {x 1 ≥0.5} , and compare the numerical results with the corresponding microscopic solution using the diffusion coefficient depicted in Figure 5.1 (left). In both cases, homogeneous Dirichlet boundary conditions are imposed on the outer boundaries. Representative cross sections of the numerical approximations obtained based on the reconstructed stiffness matrices, compared to the corresponding microscopic solutions, are shown in Figure 5.4. The numerical results indicate that robustness can be assured only  with some moderate quasi-locality. Moreover, as in the previous experiments, the quality of the results improves if is increased. 5.2. Example 2: incomplete boundary data. Next, we consider a more realistic case where the operatorL eff is only partially known. In practice, this means that coarse measurements in Ω are available only for q distinct boundary conditions in X H (q < dim X H ). In this setting, the aim is to find an effective model that not only fits the given data, but that is also able to reproduce the coarse behavior for other boundary conditions not considered as input data. The scalar coefficient A whose corresponding stiffness matrix should be recovered is shown in Figure 5.5 (left). We set H = √ 2 · 2 −5 , f = 1, q = 40, and the initial matrix S 0 H is defined as the finite element stiffness matrix based on an independent and uniformly distributed random coefficient on the coarse scale H with values between 0.1 and 10. We adapt the randomized approach used in [OY19] in the context of deep learning. Namely, in each iteration step, we randomly choose half of the available data to compute the new search direction, whereas we use all available data for the line search and for the evaluation of the functional J H . The values of the error functional J H in the first 20 iterations of the inversion algorithm are shown in Figure 5.5 (right). One can observe that classical local stiffness matrices and even the quasi-local approach with = 1 cannot significantly improve the results obtained with the initial guess, while quasi-local matrices with ≥ 2 are able to reduce the values of the functional up to a certain degree.
As in the previous subsection, we validate the outcome of the inversion algorithm by solving a diffusion problem using the reconstructed stiffness matrices and comparing the numerical results with the corresponding microscopic finite element solutions. The cross sections at x 2 = 0.5 and x 1 = 0.5 of the numerical approximations using the different stiffness matrices are shown in Figure 5.6, for the case with boundary condition u 0 (x 1 , x 2 ) = x 1 . We emphasize that, in this setting, neither the reference finite element function (black dotted line in Figure 5.6) nor a coarse measurement from it were part of the input data.
For a further comparison, we also present in Figure 5.7 the same cross sections of the numerical solutions obtained from the stiffness matrices using a full-data approach, i.e., when all available data (40 measurements) are used in every step to compute the new search direction. The reconstructed matrices behave similarly to the ones obtained with the randomized approach. However, it is worth mentioning that the randomized strategy is generally more robust in the case of incomplete boundary data, and additionally requires less computational effort. For further numerical experiments, see also [Mai20,Ch. 4]. 5.3. Discussion. The presented inversion results demonstrate that the reconstruction of the stiffness matrix assuming a local sparsity pattern with communication only between Figure 5.6. Cross sections of reconstructed functions with boundary condition u 0 (x 1 , x 2 ) = x 1 based on local stiffness matrices ( , dotted) and quasi-local ones with = 1 ( ), = 2 ( ), = 3 ( ) for Example 2 obtained from incomplete boundary data and the randomized approach. The corresponding microscopic FE function ( , dashed) is depicted as a reference but was not part of the input data. Left: Cross section at x 2 = 0.5. Right: Cross section at x 1 = 0.5. neighboring degrees of freedom does not allow us to capture effective features of the microscopic problem, while the reconstruction based on a quasi-local approach, especially with ≥ 2, is able to mimic the effective behavior.
Furthermore, some quasi-locality appears to allow for robustness with respect to different right-hand sides, a property which allows us to employ the reconstructed effective model for the simulation of other scenarios, assuming that the microscopic properties remain unchanged.
Our experiments also indicate that a lower bound on seems to be necessary similar to the forward setting, where needs to be increased for smaller values of H ( | log H|) to obtain improvements in the first place. In that sense, our findings also deviate from the Figure 5.7. Cross sections of reconstructed functions with boundary condition u 0 (x 1 , x 2 ) = x 1 based on local stiffness matrices ( , dotted) and quasi-local ones with = 1 ( ), = 2 ( ), = 3 ( ) for Example 2 obtained from incomplete boundary data and the full-data approach. The corresponding microscopic FE function ( , dashed) is depicted as a reference but was not part of the input data. Left: Cross section at x 2 = 0.5. Right: Cross section at x 1 = 0.5. numerical results in [GGS12] which indicate that truly local numerical homogenization might always be possible.

Conclusion
We proposed a strategy to reconstruct the effective behavior of solutions of a multiscale PDE model which involves a coefficient varying on a microscopic scale. The approach is motivated by the effective models (represented by effective stiffness matrices) obtained by numerical homogenization. The aim is to provide a first step towards a reasonable surrogate in the inverse multiscale setting, which is characterized by a mismatch between coarse data scale and microscopic quantities. The method relies on a quasi-local behavior in the sense that the reconstructed system matrices have a slightly denser sparsity pattern than standard finite element matrices, and this allows us to recover the behavior related to characteristic microscopic features of the solutions without requiring numerical computations on the microscopic scale. The method has been numerically validated on a prototypical model problem, considering a stationary linear elliptic diffusion problem with inhomogeneous boundary conditions. Further, even the case of incomplete boundary data can be handled and ideas from learning-type methods may be adopted.
A possible future extension of the approach includes the identification of further information about the underlying coefficient, e.g., geometric features, based on the reconstructed effective models. Further, more involved combinations with learning-type techniques to deal with incomplete data could be studied, as well as adaptivity with respect to the parameter .