Restoring the Discontinuous Heat Equation Source Using Sparse Boundary Data and Dynamic Sensors

This study focuses on addressing the inverse source problem associated with the parabolic equation. We rely on sparse boundary flux data as our measurements, which are acquired from a restricted section of the boundary. While it has been established that utilizing sparse boundary flux data can enable source recovery, the presence of a limited number of observation sensors poses a challenge for accurately tracing the inverse quantity of interest. To overcome this limitation, we introduce a sampling algorithm grounded in Langevin dynamics that incorporates dynamic sensors to capture the flux information. Furthermore, we propose and discuss two distinct sensor migration strategies. Remarkably, our findings demonstrate that even with only two observation sensors at our disposal, it remains feasible to successfully reconstruct the high-dimensional unknown parameters.


Introduction 1.Mathematical model
We give the mathematical statement of our interested inverse problem.Firstly the heat equation is given as follows: (∂ t − ∆)u(x, t) = χ D , (x, t) ∈ Ω × (0, T ), u(x, t) = 0, (x, t) ∈ ∂Ω × (0, T ) ∪ Ω × {0}. (1.1) Here Ω ⊂ R 2 is the two-dimensional unit disc.In equation (1.1), the discontinuous source χ D is the characteristic function and the support D ⊂ Ω is unknown.More precisely, if x is included by D, the value of χ D is one; otherwise, it equals zero.We assume that ∂D is sufficiently smooth.In this work, we will use the boundary flux data observed by the dynamic moving sensors to recover the unknown support D.
Here z(t) is the observation point and located in the boundary ∂Ω; Γ(t) is the observed area.We employ the notations z(t) and Γ(t) to emphasize that the sensor locations are time-dependent.Our analysis centers on sparse boundary data, where Γ represents a very small subset of the boundary.The limited availability of observations on this sparse boundary presents substantial challenges in recovering χ D .To address these challenges and verify the proposed Algorithm 1, we opt to rely solely on flux information obtained from two sensors, further enhancing the complexity of the problem.

Background and literature
In practical applications of inverse problems, sparse boundary data holds great importance.Consider Equation (1.1), which can describe the diffusion of pollutants in a system, where the solution u represents pollutant concentration.In such scenarios, the region of the pollutant source, denoted as D, should be heavily contaminated.To safeguard the health of engineers, the utilization of boundary measurements becomes imperative.
Furthermore, in the realm of practical inverse problems, the project's cost cannot be overlooked.This includes expenses related to equipment, labor, computation, and more.This cost-conscious perspective has motivated numerous researchers to seek solutions to inverse problems using sparse data.
However, when opting for sparse boundary data, it raises a pertinent question: how do we determine the sensors' location?Although one can establish the well-posedness of the inverse PDE problem, the sensors' locations significantly influence the quality of numerical simulation.Unfortunately, in the context of inverse problems, we often have to select observation areas randomly.In practice, it is entirely possible to initially choose suboptimal observation points, which can pose challenges.
To address this dilemma and ensure the quality of our approximations, we propose the concept of using moving observation points.Specifically, if the initial observation points are situated in less favorable locations, we can relocate the facilities to more suitable positions.Consequently, the data obtained from these improved locations will compensate for the adverse effects of the suboptimal points and enhance the accuracy of our approximations.
In line with this idea, the fundamental question we must address is how to determine whether a location is advantageous or not; that is, how to decide when and where to move the observation points.
The inverse source problem of the heat equation (∂ t −∆)u = F is a classical field in the literature of inverse problems, and abundant academic achievements are generated.
See [4,14,29,30,17,15,16,12,23] and the references therein.Nowadays, due to the significant application value of the sparse data, more and more researchers have paid attention to this field.Here we list several references on the inverse problems with sparse data.In [12], the authors consider the inverse source problem of the heat equation on the two-dimensional unit disc.They recover the space-dependent source f (x) by the flux data at finite points of the boundary.The conclusion in [12] is promoted by [30], in which the variable separable source f 1 (x)f 2 (t) can be uniquely determined by the sparse boundary data.Note that in [30], the spatial term f 1 (x) and temporal term f 2 (t) are both unknown.In [20,23], the authors further recover the semi-discrete unknown source, which is investigated in this work, and the considered equation is the parabolic equation on a general domain.Also, for the geometric inverse problems, the authors in [18,11] recover the manifold by the image of a single specific function.We should note that for such geometric inverse problems, people usually use a whole operator as the measurements.

Main result and outline
We summarize the contributions as follows.
1. We study inverse problems with sparse measurements and moving observation points.
2. We introduce a sampling-based approach and a dynamic sensor migration algorithm designed for source tracing.Moreover, we present two migration strategies and validate their effectiveness through two challenging numerical examples.
The rest of this article is organized as follows.In Section 2, we will provide an overview of the existing problems.We then introduce the Bayesian methodology and the Langevin-based algorithm in Section 3, in which the uniqueness of this inverse problem is also discussed.In particular, we will discuss the proposed sampling methods and proposed dynamical sensors migration algorithm in Section 3.1 and 3.2.Finally, we will verify the performance and introduce the sensors migration strategies in Section 4.

The setting of the inverse problem
This article is an extension of reference [29,30], in which the authors prove the uniqueness theorem of the inverse source problem under sparse data.More precisely, under suitable conditions, the flux data generated from two points of the boundary can uniquely determine the unknown source.However, as we mention in the introduction, when we solve this inverse problem numerically, there is a natural question that how to determine the locations of the observation points.This question holds paramount importance in computational work because the precision of the numerical approximation is profoundly influenced by the placement of observation points.To substantiate this assertion, we test a numerical experiment by the algorithm proposed in [30] and obtain the comparison (please see Figure 2.1).
Figure 2.1: Numerical approximations with different observation points.The target source location is represented by the blue dashed curve, while the brown curve depicts the source recovered using only two sensors.The black dots situated along the boundary mark the positions of these sensors.It is evident from the images that the accuracy of the recovery process is contingent upon the placement of the sensors.
In Figure 2.1, the dotted blue line represents the support of the exact source, while the dashed red line corresponds to the support of the approximated source.The black points located along the boundary represent the observation points.It's evident that the support of the exact source in both figures is the same; the only difference lies in the placement of the observation points.Consequently, the two approximations exhibit significant disparities.
Evidently, in the left figure, the observation points are positioned unfavorably, earning them the label of 'bad' points, while those in the right figure are considered 'good'.
However, it's important to recognize that in practical applications of inverse problems, obtaining the exact solution for the unknown source is often infeasible.Therefore, we are left with the task of selecting observation points largely based on chance.
If, by chance, we initially select 'bad' points and do not alter their positions, it results in the production of poor numerical outcomes, as demonstrated in the left part of Figure 2.1.This underscores the reason why we endeavor to leverage data from moving observation points to offset the approximations adversely affected by these 'bad' observation points.
Following [12,30], the unknown support D of the source can be uniquely determined when the observed area Γ is fixed and only consists of two chosen points on the boundary.Hence, in this work, we will suppose that Γ(t) consists of two points of ∂Ω for each t ∈ (0, T ).With the idea of moving observation points, the measurements we used in this work will be represented as 2) The notation z j (t) is the orbit of each observation point, and the essential difficulty of this work is how to determine the orbit z j .We will provide some rules for determining z j (t) and give the corresponding numerical experiments.These will be discussed in the next sections.

Bayesian framework
A highly effective approach for tackling the inverse problem involves adopting the Bayesian inference framework, as referenced in prior works [31,21,9,5,33].In this context, we define the dataset as d = {d i } n d i=1 and represent the prior distribution of ξ as p(ξ).Utilizing the Bayesian formula, we can express this relationship as follows: where p(•|ξ) is a likelihood function given the parameter ξ.where p(•|ξ) represents a likelihood function conditioned on the parameter ξ.
We denote the mapping from the target parameters ξ to the observations as the forward model, labeled as g.Given that observation data inherently contain noise and that noise is also introduced during the forward evaluation, we make the assumption that the likelihood follows a normal distribution with a mean of d and a variance of σ 2 .Specifically, the likelihood function takes the following form: One of the key advantages of employing the Bayesian inference method is our ability to obtain the distribution of the unknown parameters.This capability allows us to quantitatively assess the uncertainty associated with quantities of interest (QoI), such as the outputs of the forward model.The trajectory of the observation points is subsequently determined in a sequential manner using posterior samples of the QoI.

Langevin diffusion and Bayesian sampling
In this study, we will employ the Bayesian approach to ascertain the trajectory of observation points and iteratively address the inverse source problem.To attain the desired posterior distribution p(ξ|d), various MCMC methods can be employed for sampling.Among these methods, the Langevin diffusion-based MCMC [7,8,26,3] stands out as a computationally efficient choice, particularly when dealing with highdimensional and multimodal distributions [26,19].Specifically, let us define the energy function U (ξ) as, We are interested in sampling from the following distribution: where ξ ∈ R p , τ is the temperature parameter, and the energy function U is known.
One method for generating samples from π(ξ) is to simulate a stochastic process whose stationary distribution is π.Under certain regularity conditions (see, e.g., [2,28]), the following stochastic differential equation (SDE), known as the preconditioned Langevin diffusion (LD) [27,25], has a stationary distribution π: where L is an arbitrary symmetric, positive-definite matrix, W t is the p-dimensional standard Brownian motion.Motivated by the replica-exchange methods [22,21], [24] extends the preconditioned Langevin diffusion-based methods to multiple chains to accelerate distribution simulation.
In practical sampling scenarios, computing the gradient of the energy function can be computationally expensive, especially when dealing with PDE-based inverse problems that involve time-consuming forward model simulations.Moreover, in this work, the source term exhibits discontinuities concerning the unknown parameters.Therefore, we opt for the gradient-free Metropolis-Hasting (MH) MCMC method [32,10].
To address the inverse problem, we will employ the adaptive preconditioned Crank-Nicolson MCMC (pCN-MCMC) method within a Bayesian inference framework.

Adaptive preconditioned-Crank-Nicolson (pCN)
As the sampling acceptance rate is independent of the dimension of parameters, the method of pCN-MCMC [6] stands as an efficient strategy for addressing large-scale inference challenges.Specifically, it is based on the SDE which is a special case of equation (3.3),where the prior is the Gaussian distribution N (0, B), the gradient of the log-likelihood vanishes in the drift term, and the temperature is τ = 1.The proposal is generated by applying the Crank-Nicolson (CN) schemes [1,6] to the SDE (3.4), where w ∼ N (0, I p ) represents samples from a multivariate normal distribution, and δ denotes the step size.When L = B is chosen, the pCN proposal is resulted, where Moreover, the equation (3.5) can also be rewritten as [13] ξ where , yielding the adaptive pCN scheme, it is used to accelerate the convergence of the pCN MCMC method.One of the optimal choices of C is the posterior covariance matrix, which is not directly available but can be determined by the historical samples empirically.For both the pCN and adaptive pCN scheme, the acceptance rate has the form where Φ(ξ) = log p(d|ξ), please refer to Algorithm 2 for more details on the adaptive pCN sampling method.

Methodology
As mentioned earlier, the uniqueness of this inverse problem has undergone comprehensive scrutiny in prior works [29,30].In this study, we align ourselves with the principles outlined in these references [29,30] to articulate the uniqueness theorem.
Theorem 1.Given two unknown supports D 1 and D 2 which have sufficiently smooth boundaries, we denote the solutions corresponding to D j by u j , j = 1, 2. Recalling the data (2.2), we assume that there exists ϵ 0 > 0 such that z j (t) = θ j for t ∈ (0, ϵ 0 ), j = 1, 2, and where Q is the set of rational numbers.If Proof.The proof of Theorem 1 follows from [29, Theorem 3.1] and the analyticity of the natural exponential function straightforwardly.
Theorem 1 establishes that, under appropriate conditions, we can ascertain the precise support of the source using the data (2.2).However, extensive experimental findings, as presented in [23,30], suggest that relying solely on fixed-sensor flux observations may not yield accurate source predictions.Furthermore, in both sampling methods [23,31] and the classical Levenberg-Marquardt iteration, the accuracy of source capture is greatly contingent on sensor placement, particularly when only using a limited number of sensors.
Inspired by Theorem 1, we introduce an algorithm that introduces variability in the locations of observation points (sensors).For a comprehensive understanding of the procedure, we direct your attention to Algorithms 1 and 2. It's worth noting that the option to select new sensor locations at random is a possibility.However, in Section 4.1.1 and 4.2.1, we delve into two distinct approaches for altering the sensor locations, which we detail further.
Algorithm 1: Peudo-algorithm for dynamical locations determination.Pick up and set t to be the next observation time T k+1 from t ′ ; 4 Numerical experiments.
In this section, we will conduct numerical experiments to test the performance of the algorithm.We shall delve into scenarios involving two distinct unknown sources, 11 Store the posterior samples of the parameters and PDE solutions at time t.
each characterized by disparate geometries and dimensions pertaining to the unknown parameters.Across both sets of problems, the task involves measuring boundary observation flux at two distinct points along the boundary-an endeavor that has proven to be particularly challenging for conventional implementation methodologies [30,29].Employing a comparative approach alongside the fixed observation points method, our analysis reveals a consistent trend: the methods we propose adeptly capture the source.For the convenience of solving the forward model, we transform the Cartesian plane to the polar coordinates, i.e., the domain is then transformed to be Ω := {(r, θ)|0 < r < 1, 0 ≤ θ ≤ 2π}.The equation becomes The forward model (4.9) is solved by the finite difference method.The measurement locations are sequentially determined during the inference following Algorithm 1.

Circle-shape source
In this experiment, we consider the following unknown domain D, where η = [η 1 , η 2 ] ⊺ are the unknowns.The corresponding expression under the polar coordinates is where (ρ, ω) are the unknown parameters.As the range of the samples is (−∞, +∞) for this algorithm, while the range of our parameters is finite, we use the bijection map ⊺ to the intervals (0, 1) and (0, 2π), respectively.The ξ become the unknowns endowed with the prior N (0, I p ) to be inferred.The relationship between ξ and the location of the center is The statistical properties of η can be estimated through the posterior samples of ξ.
The forward model is computed on a spatial grid of dimensions 33 × 36, with a predefined strength of b = 50.The observed data is generated using ground truth parameters η = [0, 0.5] ⊺ , implying that the true parameter vector ξ is given by (ξ 1 , ξ 2 ) ⊺ = (0, −1) ⊺ .The measurement error variance is σ 2 = 0.05 2 .We perform a total of N 1 = 0 initial iterations, followed by N = 10 4 iterations, with updates occurring every k 0 = 2.5 × 10 3 iterations.In both the pCN and adaptive pCN schemes, the tunable parameters β 1 and β 2 are set to be equal to ensure an acceptance rate of 30% to 40%.

Sensors migration rule
The parameters ξ 1 and ξ 2 represent the radius and center coordinates of the target circle.As the observation locations approach the target area, the gathered data becomes more informative.Since the target area is circular, the proximity of the observation locations to the domain can be quantified by measuring the distance between the center and the observation locations.Hence, the moving observation locations are determined based on the samples of ξ 2 , specifically, New sensor location = (⌊ω/h θ ⌋, ⌈ω/h θ ⌉).where ω is the posterior mean of ω(ξ 2 ), h θ is the spatial step of the θ direction, yielding ⌊ω/h θ ⌋ and ⌈ω/h θ ⌉ the index of locations along the boundary.
Firstly, we take two measured locations along the edge randomly at time T 1 = 0.5, and draw samples from the posterior density.We then calculate the posterior mean of ω(ξ 2 ) as ω, and determine the new sensors location index as ⌊ω/h θ ⌋ and ⌈ω/h θ ⌉ at the next time layer T 2 .Obviously, the selected observation points p 1 and p 2 are adjacent.The process is repeated until the observation time is t = T 3 , and the forward model simulation and Bayesian inference are then stopped.The results of the dynamic procedure are displayed in Figure 4.2.
For comparison purposes, we maintain the sensor location fixed at ( 2218 π, 30  18 π), and flux measurements are taken at time instances t ′ = [0.5, 1, 1.5] ⊺ , respectively.The measurement error remains consistent for both the fixed location strategy and the strategy involving sequentially determined locations.In other words, for each observed time T i , the noise introduced to the flux at the fixed locations ( 2218 π, 30  18 π) is identical to the noise introduced to the flux at the sequentially determined locations.In Figure 4.3, the posterior mean of the samples is represented by the red circle, offering an accurate estimate of the target region's center.However, Figure 4.4 reveals a distinct pattern: the variance of the samples acquired through sequentially determined sensor placements is significantly reduced compared to those obtained using the fixed measurement strategy.This outcome is a result of the fact that sequentially determined measurements provide more informative data, leading to a substantial reduction in uncertainty.We additionally performed a comparative analysis between our proposed method and a random selection strategy for measurement locations.To gain a more comprehensive understanding of the stochastic nature, we conducted two experiments with distinct measurement locations.The trace plots of posterior samples for ξ are presented in Figure 4.5.It is evident that the dynamical sensors placement method, utilizing both the proposed varying method as detailed in Section 4.1.1,and the random varying strategy successfully capture the inverse Quantities of Interest (QoIs).Nevertheless, it is important to highlight that the random varying strategy results in a broader confidence interval.

High dimensional peanut-shape source
In this example, we examine a problem of higher dimensions with asymmetrical characteristics [29].Specifically, the unknown domain has the boundary ∂D that can be parameterized as with a smooth, periodic function q(θ; ξ) ∈ (0, 1), which has the form To ensure the smoothness of the approximation, we follow the paper [29] and set the penalty term to be the H 2 norm of q(θ; ξ), which implies the Gaussian prior density N (0, B) for ξ, where the covariance matrix B ∈ R (2M +1)×(2M +1) is diagonal, with entries The true parameter is ξ = [1, 0, 0, 0, 0.3] ⊺ , the strength is preset as b = 10, the variance of the measurement error is σ 2 = 0.01 2 .The iteration numbers are N 1 = 1 × 10 3 , N = 1.5 × 10 4 , the update frequent number is k 0 = 2.5 × 10 3 .The parameters in the pCN and adaptive pCN scheme are tuned separately, both lead to the acceptance rate 30% ∼ 40%.
We first assess the random sensor relocation method and showcase the corresponding outcomes in Figure 4.6 (b).We compare the proposed method with the fixed-sensor approach, and as shown in  the data collected by randomly relocating sensors proves to be sufficiently informative for capturing the intricate shape of the high-dimensional source.As depicted in Figure 4.7, the samples derived from the randomly moving sensors strategy tend to cluster around or in proximity to the inverse Quantities of Interest (QoIs) denoted as ξ, whereas the fixed-sensor approach struggles to identify most of the QoIs, with the exception of ξ 1 .

Sensors migration rule
The random sensor relocation strategy demonstrates superior performance compared to the fixed-sensor approach.However, it remains challenging for us to determine the optimal time to relocate the sensors.To address this issue, we introduce a second effective approach that can help us identify when it is appropriate to cease relocating the sensors.
Our approach revolves around relocating the two sensors to locations characterized by the highest flux variance observed across all sampling iterations.To be precise,  (4.12) The second sensor is subsequently moved to the position exhibiting the second largest variance.It's important to emphasize that this strategy incurs no additional costs.
To compute the flux at the sensor locations, we must evaluate the equation's solution using the sampled source, guided by the acceptance rate (3.8).This enables us to compute the flux variance across all sampling iterations at any given point along the domain boundary.Notably, higher variance corresponds to a larger confidence interval and a less confident prediction.Consequently, our strategy involves relocating the sensors to positions characterized by the greatest variance.
Following the aforementioned rule, we sequentially ascertain the measurement locations and present the posterior mean of the prediction samples in Figure 4.8.Notably, the locations at time T 5 coincide with those at time T 3 .Specifically, the suggested locations at time T 5 are ( 21 18 π, 20  18 π).This recurrence of the location 20  18 π also appears at T 4 , indicating a pattern of repeated cycling between the locations ( 4 18 π, 5 18 π) and ( 1918 π, 20 18 π).Consequently, we choose to conclude the sequential process at time T 5 .This method offers a means to determine when to cease taking measurements, providing insights into the richness of the collected data to some extent.The subfigure includes the 1.5σ confidence interval of the posterior samples for q(θ; ξ), within which the actual q(θ; ξ) value is predicted.

Concluding remarks
In this work, we consider the inverse source problem of the heat equation.We use the boundary measurements where the observation sensors can be moved.The Bayesian method is used and we attempt some numerical examples.The numerical results illustrate that the Bayesian approach is feasible to solve such an inverse problem which uses the data from the moving observation points.
We outline several directions for future research.First, it is worth exploring the application of the PDE (1.1) in a more general domain.In this study, we confined our attention to the unit disc in R 2 , which possesses advantageous geometric properties.Expanding our analysis to encompass general domains in R 2 or R 3 would undoubtedly introduce additional complexities to the uniqueness argument and reconstruction algorithms.
Second, a significant open question pertains to the optimal strategy for relocating observation sensors.While we proposed a migration rule (4.12) in this work, we acknowledge that it lacks rigorous mathematical proofs, rendering it somewhat empirical in nature.Consequently, an essential avenue for future research involves a thorough investigation of this migration rule to establish its validity and optimize its performance.

1 5 6 7 Set initial time t 0 = t; 8
Initialization: two observation points (sensors) p 1 and p 2 , constant n, observation time stamps t ′ = [T 1 , ..., T n ] ⊺ , current observation time t = T 1 , initial PDE simulation time t 0 = 0, and the total sampling iteration limits K;2 for k = 1 to K − 1 do 3Observe the flux at p 1 and p 2 ; 4 Simulate the PDE solution until t from t 0 for the computation of Φ(•); Obtain posterior samples {ξ i } N i=1 and {U i } N i=1 according to the sampling schemes presented in Algorithm 2 ;Change two observation points p 1 and p 2 basing on posterior samples ; Set the mean of {U i } N i=1 to be the initial value of the PDE ;

9
Consider the problem ∂u ∂t − ∆u = bχ D (x), t ∈ [0, T ] the physical domain is Ω := {(x, y)|x 2 + y 2 < 1}.b represents the predefined intensity of the source term.χ D denotes the characteristic function related to the domain D, where D is the unspecified region in need of reconstruction.We will specify D in the experiments later.The observations are the flux ∂u ∂ν collected at some time and locations along ∂Ω.

Figure 4 . 2 :
Figure 4.2: Estimates of η using posterior sample mean.The center of the blue dashed circle is the true η, while the center of the red dashed circle is the average of all prediction samples.Flux measurements are obtained through two sensors (black dots) situated along the simulation boundary, with their positions dynamically adjusted over time according to Algorithm 1. Specifically, the simulation is first performed to compute the flux at T 1 = 0.5, with sensors located at ( 22 18 π, 30 18 π).Subsequently, flux calculations are performed at T 2 = 1, with sensors positioned at ( 10 18 π, 11 18 π), and at T 3 = 1.5 with sensor placement at ( 10 18 π, 9 18 π).For further details on the sensor migration strategy, refer to Section 4.2.1

Figure 4 . 4 :
Figure 4.4: Trace plots illustrating posterior samples across distinct measurement strategies.The blue curve represents samples for ξ 1 , while the green curve corresponds to ξ 2 .(a) Utilizing sequentially determined locations as recommended by Algorithm 1.(b) Employing fixed sensors.The visual analysis of the plots indicates that the strategy involving fixed sensor positions results in a broader confidence interval.Conversely, the proposed strategy yields predictions with significantly reduced variance.

Figure 4 . 5 :
Figure 4.5: Trace plots illustrating posterior samples across distinct measurement strategies.The blue curve represents samples for ξ 1 , while the green curve corresponds to ξ 2 .Figure (a): utilizing sequentially determined locations as recommended by Algorithm 1 with dynamical sensor placement strategy detailed in Section 4.1.1.Figure(b) and (c): using Algorithm 1 with random sensor placement strategy; please note that we repeat the experiments twice to comprehensively study the performance.It is evident that adopting a dynamic sensor placement strategy (left figure) leads to a narrower confidence band, signifying decreased variance among all the samples.Nevertheless, even randomly varying sensors, as illustrated in (b) and (c), still exhibit considerably lower variance compared to the fixed sensor approach, as demonstrated in Figure 4.4 (b).
Figure 4.5: Trace plots illustrating posterior samples across distinct measurement strategies.The blue curve represents samples for ξ 1 , while the green curve corresponds to ξ 2 .Figure (a): utilizing sequentially determined locations as recommended by Algorithm 1 with dynamical sensor placement strategy detailed in Section 4.1.1.Figure(b) and (c): using Algorithm 1 with random sensor placement strategy; please note that we repeat the experiments twice to comprehensively study the performance.It is evident that adopting a dynamic sensor placement strategy (left figure) leads to a narrower confidence band, signifying decreased variance among all the samples.Nevertheless, even randomly varying sensors, as illustrated in (b) and (c), still exhibit considerably lower variance compared to the fixed sensor approach, as demonstrated in Figure 4.4 (b).

Figure 4 .
6 (a), the limited data available to the fixedlocation strategy hinders its ability to accurately capture the shape of the target source D. Despite some overlap in measurement locations across different time points,
Algorithm 2: The adaptive pCN MCMC.Input: The covariance matrix B of the prior, parameters β 1 , β 2 , noise of the observed data σ, the initial guess ξ 0 , current observation time t.The number of iterations N 1 and N , the update frequent number k 0 .Output: Posterior samples {ξ i } N i=1 and the PDE solution {U i } N i=1 at time t. 1 for k = 1 : N 1 do Generate the candidate ξ * by the pCN scheme (3.6) Accept ξ * as ξ k with the rate α(ξ * , ξ k−1 ) 4 Collect the posterior samples of the parameters {ξ i } N 1 i=1 and the PDE solution {U i } N 1 i=1 at time t. 5 Compute the empirical covariance matrix of the {ξ i } N 1 i=1 as C. 6 for k = N 1 + 1 : N do 3 9 Propose the candidate ξ * by the adaptive pCN scheme (3.7) 10 Accept ξ * as ξ k with the probability α(ξ * , ξ k−1 )