The tutorial consists in the state reconstruction of a 1D heat transfer problem where we have some measurement points inside the domain but partially known boundary conditions.
The equations of the true problem are
\begin{eqnarray*} \rho C_p\frac{\partial T}{\partial t} + k \Delta T &=& 0 &\text{ in } (0,1) \times (0, t_f]\\ \nabla T \cdot \mathbf{n} &=& g_{true}(t) &\text{ on } x=0 \times (0, t_f]\\ T &=& 0.5 &\text{ on } x=1 \times (0, t_f]\\ T &=& T_0 &\text{ on } (0,1) \times t=0 \end{eqnarray*}
where \(\rho\), \(C_p\), \(k\) and \(T_0\) are defined by the user in the ITHACAdict and 0/T file while \(g_{true}(t) = 5 t + 2\) and \(T_0 = 1\).
The objective of the tutorial is to reconstruct the true solution having some noisy temperature measurements defined in the measurementsDict and the model
\begin{eqnarray*} \rho C_p\frac{\partial T}{\partial t} + k \Delta T &=& 0 &\text{ in } (0,1) \times (0, t_f]\\ \nabla T \cdot \mathbf{n} &=& g(t) &\text{ on } x=0 \times (0, t_f]\\ T &=& 0.5 &\text{ on } x=1 \times (0, t_f]\\ T &=& T_0 &\text{ on } (0,1) \times t=0 \end{eqnarray*}
where the boundary condition \(g(t)=5t\) is different from the right one \(g_{true}\).
To solve this problem, we use the Ensemble Kalman Filter implemented in the ITHACAmuq::muq2ithaca::EnsembleKalmanFilter method. We start by assuming a Gaussian prior for the state (the temperature in this case)
\[ \omega_{prior} \sim \mathcal{N}(T_0, \sigma_{prior}I). \]
Then, we sample it creating the prior ensemble
\[ E_{prior}^0. \]
Now, at each timestep, we perform a forecasting step. In the forecasting step, we solve the direct problem with the "wrong" heat flux for each of the member of the ensemble adding to the solution the model error
\[ \omega_{model} \sim \mathcal{N}(0, \sigma_{model}I). \]
This way, we produced at each timestep a posterior ensamble
\[ E_{post}^n. \]
If there are no measurements available at this timestep
\[ E_{post}^n = E_{prior}^{n-1}. \]
However, when we have available measurements, we perform a Kalman filter correction step. The Kalman filter updates the posterior ensamble such that
\[ \hat{E}_{post}^n = E_{post}^{n} + Z M, \]
where
\[ Z = \frac{1}{Nseeds - 1} A (HA)^T, \]
\[ A = E_{prior}(i) - \mu_{prior}, \]
\mu_{prior} being a matrix having the mean value on the ensemble repeated on each column. \(HA\) is similar to \(A\) but with the value of the state in the observation point,
\[ M = P Y, \]
\[ P = \frac{1}{Nseeds - 1} (HA) (HA)^T + \sigma_{meas} I, \]
and \(Y\) is the matrix containing the difference between the measured and computed state at the measurement points.
The image below shows the true end reconstructed temperature at \(x=0.5\) together with the 95% quantile
Here there's the plain code
1.10.0
