This model considers an SPDE over a domain Ω which is partitioned into k subdomains Ωd, d ∈ {1, …, k}, where ∪d = 1kΩd = Ω. A common marginal variance is assumed but the range can be particular to each Ωd, rd.
From Bakka et al. (2019), the precision matrix is $$ \mathbf{Q} = \frac{1}{\sigma^2}\mathbf{R}\mathbf{\tilde{C}}^{-1}\mathbf{R} \textrm{ for } \mathbf{R}_r = \mathbf{C} + \frac{1}{8}\sum_{d=1}^kr_d^2\mathbf{G}_d , \;\;\; \mathbf{\tilde{C}}_r = \frac{\pi}{2}\sum_{d=1}^kr_d^2\mathbf{\tilde{C}}_d $$ where σ2 is the marginal variance. The Finite Element Method - FEM matrices: C, defined as Ci, j = ⟨ψi, ψj⟩ = ∫Ωψi(s)ψj(s)∂s, computed over the whole domain, while Gd and C̃d are defined as a pair of matrices for each subdomain (Gd)i, j = ⟨1Ωd∇ψi, ∇ψj⟩ = ∫Ωd∇ψi(s)∇ψj(s)∂s and (C̃d)i, i = ⟨1Ωdψi, 1⟩ = ∫Ωdψi(s)∂s.
In the case when r = r1 = r2 = … = rk we have $\mathbf{R}_r = \mathbf{C}+\frac{r^2}{8}\mathbf{G}$ and $\mathbf{\tilde{C}}_r = \frac{\pi r^2}{2}\mathbf{\tilde{C}}$ giving $$ \mathbf{Q} = \frac{2}{\pi\sigma^2}( \frac{1}{r^2}\mathbf{C}\mathbf{\tilde{C}}^{-1}\mathbf{C} + \frac{1}{8}\mathbf{C}\mathbf{\tilde{C}}^{-1}\mathbf{G} + \frac{1}{8}\mathbf{G}\mathbf{\tilde{C}}^{-1}\mathbf{C} + \frac{r^2}{64}\mathbf{G}\mathbf{\tilde{C}}^{-1}\mathbf{G} ) $$ which coincides with the stationary case in Lindgren and Rue (2015), when using $\tilde{\mathbf{C}}$ in place of C.
In practice we define rd as rd = pdr, for known p1, …, pk constants. This gives $$ \mathbf{\tilde{C}}_r = \frac{\pi r^2}{2}\sum_{d=1}^kp_d^2\mathbf{\tilde{C}}_d = \frac{\pi r^2}{2} \mathbf{\tilde{C}}_{p_1,\ldots,p_k} \textrm{ and } \frac{1}{8}\sum_{d=1}^kr_d^2\mathbf{G}_d = \frac{r^2}{8}\sum_{d=1}^kp_d^2\mathbf{\tilde{G}}_d = \frac{r^2}{8}\mathbf{\tilde{G}}_{p_1,\ldots,p_k} $$ where C̃p1, …, pk and G̃p1, …, pk are pre-computed.