# Fast and Efficient Lunar Finite Element Gravity Model
## 2 Orthogonal Approximation
\begin{bmatrix}\frac{T_{0}(\xi_{0})}{2} &T_{0}(\xi_{1})&\cdots&T_{0}(\xi_{M_{x-1}})&\frac{T_{0}(\xi_{M_{x}})}{2}\\ T_{1}(\xi_{0})&2T_{1}(\xi_{1})&\cdots&2T_{1}(\xi_{M_{x-1}})&T_{1}(\xi_{M_{x}}) \\ \vdots&\vdots&\ddots&\vdots&\vdots\\ T_{N_{x}-1}(\xi_{0})&2T_{N_{x}-1}(\xi_{1})&\cdots&2T_{N_{x}-1}(\xi_{M_{x-1}})&T_ {N_{x}-1}(\xi_{M_{x}})\\ T_{N_{x}}(\xi_{0})&2T_{N_{x}}(\xi_{1})&\cdots&2T_{N_{x}}(\xi_{M_{x-1}})&T_{N_{x} }(\xi_{M_{x}})\end{bmatrix} ,a= \begin{bmatrix}a_{0}\\ a_{1}\\ \vdots\\ a_{N}\end{bmatrix} \]  
Similarly, for more than one independent variable, let us consider the function \(f(x,y)\), with the two independent variables \(x\) and \(y\):  
\[g(x,y),\{x_{min}<x<x_{max}\},\{y_{min}<y<y_{max}\} \tag{10}\]  
As in the one-variable approximation example, let us use a set of new independent variables, \(\xi\) and \(\eta\), to reduce the computational burden of dimensionality which leads to the normalized forward and inverse functions:  
\[ \begin{split} Forward\left\{\begin{aligned} \xi(x)&=-1+2(x-x_{min})/(x_{max}-x_{min})\\ \eta(x)&=-1+2(y-y_{min})/(y_{max}-y_{min})\\ \end{aligned}\right.\\ Inverse\left\{\begin{aligned} x(\xi)&=x_{min}+\ (\xi+1)\ (x_{max}-x_{min})/2\\ y(\eta)&=y_{min}+\ (\xi+1)\ (y_{max}-y_{min})/2\end{aligned}\right. \end{split}  \tag{11}\]  
Substituting these equations into Equation (11) provides us with the function to be approximated:  
\[f(\xi,\eta)=g(x(\xi),y(\eta))=g\Big{(}x_{\min}+(\xi+1)\frac{(x_{\max}-x_{\min })}{2},y_{\min}+(\eta+1)\frac{(y_{\max}-y_{\min})}{2}\Big{)} \tag{12}\]  
and in the discrete measurement case tested at sample nodes, \(\{\xi_{0},\xi_{1},...,\xi_{M_{x}},\)\(M_{x}\geq N_{x}\}\), \(\{\eta_{0},\eta_{1},...,\eta_{M_{y}};M_{y}\geq N_{y}\}\), the residual approximation error at each node is:  
\[r_{ij}=f(\xi_{i},\eta_{j})-\sum_{\alpha=0}^{N_{x}}\sum_{\beta=0}^{N_{y}}a_{ \alpha\beta}\phi_{\alpha\beta}(\xi_{i},\eta_{j});\ i=0,1,...,M_{x},\,j=0,1,...,M_{y} \tag{13}\]  
which can be displayed in vector-matrix notation: \(r=f-\Phi a\) where \(\Phi=\sum a_{ij}\). When in matrix form, a common issue that needs to be considered is dimensionality. To avoid this problem, efficient choices of basis functions and associated sample points must be made. This technique is examplified in [46]. Referring to section 1.6.2, \(\Phi\) is formed as the Kronecker product of two basic matrices:  
\[\Phi= \begin{bmatrix}1&\xi_{0}&\xi_{0}^{2}\\ 1&\xi_{1}&\xi_{1}^{2}\\ 1&\xi_{2}&\xi_{2}^{2}\\ 1&\xi_{3}&\xi_{3}^{2}\end{bmatrix} \otimes \begin{bmatrix}1&\eta_{0}\\ 1&\eta_{1}\\ 1&\eta_{2}\end{bmatrix} =\Phi_{x}\otimes\Phi_{y} \tag{14}\]For our purposes, the Kronecker product, \(C=A\otimes B\), is implemented in MATLAB [47] via the \"kron\" command. Subsequently, the weighted least-square solution matrix becomes:  
\[a=(\Phi^{T}W\Phi)^{-1}\Phi^{T}f=\Big{\{}\Big{[}(\Phi_{x}^{T}\Phi_{x})^{-1}\Phi_{ x}^{T}\Big{]}\otimes\Big{[}(\Phi_{y}^{T}\Phi_{y})^{-1}\Phi_{y}^{T}\Big{]}\Big{\}}f \tag{15}\]  
As Equation (15) is \"Kronecker favorable\", it can be further represented as \(\Phi=\Phi_{x}\otimes\Phi_{y}\). Furthermore, this result can be generalized to higher dimensions as the following: \(\Phi=\Phi_{x}\otimes\Phi_{y}\otimes\Phi_{z}\), which allows for us to generalize Equation (15) into a three-dimensional approximation space \((x,y,z)\) or \((\xi,\eta,\zeta)\):  
\[a=(\[PERSON]^{T}W\Phi)^{-1}\Phi^{T}f=\Big{\{}\Big{[}(\Phi_{x}^{T}\Phi_{x})^{-1}\Phi _{x}^{T}\Big{]}\otimes\Big{[}(\Phi_{y}^{T}\Phi_{y})^{-1}\Phi_{y}^{T}\Big{]} \otimes\Big{[}(\Phi_{z}^{T}\Phi_{z})^{-1}\Phi_{z}^{T}\Big{]}\Big{\}}f \tag{16}\]  
Using the Kronecker product of three small least-squares operators to produce a larger operator allows us to approximate in higher dimensions by taking the cube root of the size of the three inverted matrices. Going up dimensions would require us only to introduce a fourth variable, which is the beauty of this Kronecker product trick. for more details about the properties of the Kronecker product and how lower-dimensional analysis can be easily \"upgraded\" to higher dimensions.