# Fast and Efficient Lunar Finite Element Gravity Model
## 2 Orthogonal Approximation
With this, the above equation is generalized to \(n\)-dimensions by taking the \(n\)th root of the inverted matrix dimension:  
\[a=(\Phi^{T}W\Phi)^{-1}\Phi^{T}Wf\equiv Cf\equiv(\bar{\Phi}^{T}\bar{\Phi})^{-1 }\bar{\Phi}^{T}\bar{f}\equiv\bar{C}\bar{f} \tag{17}\]  
where \(C\) is the Cholesky square root of \(W^{1/2}\) of the weighted matrix \(W=W^{1/2}W^{1/2}\). and \(\bar{\Phi}\equiv W^{1/2}\Phi\), \(\bar{f}=W^{1/2}f\) which makes \(\Phi\) Kronecker factorable. This final result allows the same advantages to be used in the identity-weighted matrix case. With this, Equation (17) can be expressed as:  
\[a=(\bar{\Phi}^{T}\bar{\Phi})^{-1}\bar{\Phi}^{T}\bar{f}=\Big{\{}\Big{[}(\Phi_{ x}^{T}\bar{\Phi}_{x})^{-1}\bar{\Phi}_{x}^{T}\Big{]}\otimes\Big{[}(\bar{\Phi}_{ y}^{T}\bar{\Phi}_{y})^{-1}\bar{\Phi}_{y}^{T}\Big{]}\otimes\Big{[}(\bar{ \[PERSON]}_{z}^{T}\bar{\[PERSON]}_{z})^{-1}\bar{\[PERSON]}_{z}^{T}\Big{]}\Big{\}}\bar{f} \tag{18}\]  
where \(\bar{\[PERSON]}\equiv W^{1/2}\[PERSON]\), \(\bar{\[PERSON]}_{x}\equiv W^{1/2}\[PERSON]_{x}\), \(\bar{\[PERSON]}_{y}\equiv W^{1/2}\[PERSON]_{y}\), \(\bar{\[PERSON]}_{z}\equiv W^{1/2}\[PERSON]_{z}\), and \(\bar{f}\equiv W^{1/2}f\).  
It is useful to note that \(W,W_{x},W_{y}\), and \(W_{z}\) have a Kronecker relationship: \(W=W_{x}\otimes W_{y}\otimes W_{z}\), which allows for us to sample locations corresponding to the Chebyshev polynomials. The Kronecker product trick for constructing the two-dimensional Chebyshev least-squares matrix can be provided by:  
\[\bar{C}\equiv(\bar{\Phi}^{T}\bar{\Phi})^{-1}\bar{\Phi}^{T}=(\Phi^{T}W\Phi)^{- 1}\Phi^{T}W^{1/2}=CW^{-1/2}\equiv\bar{C}_{x}\otimes\bar{C}_{y} \tag{19}\]  
where \(\bar{C}_{x}=(\bar{\Phi}_{x}^{T}\bar{\Phi}_{x})\). The most important takeaway from this step is that \(C\) is also Kronecker factorable, just as weighted matrix \(W\) is. From the above conclusions, this formulation is applied to \(n\)-dimensions. Take the following example for approximating a function in three dimensions:  
\[\bar{C}=(\Phi^{T}W\Phi)^{-1}\Phi^{T}W=\bar{C}_{x}\otimes\bar{C}_{y}\otimes \bar{C}_{z},W\equiv W_{x}\otimes W_{y}\otimes W_{z},\bar{f}=W^{1/2}f \tag{20}\]  
For final generalized function, consider a function of \(n\) independent variables:  
\[g(x_{1},x_{2},...,x_{n}),\{x_{i,min}\leq x_{i}\leq x_{i,max}\},i=1,2,...,n \tag{21}\]  
with the following forward and inverse functions:  
\[ \begin{split} Forward\;\big{\{}\;\bar{\xi}_{i}(x_{i})& =-1+2(x_{i}-x_{i,min})/(x_{i,max}-x_{i,min})\\ Inverse\;\big{\{}\;x_{i}(\bar{\xi}_{i})&=x_{i,min}+ \;(\bar{\xi}_{i}+1)\;(x_{i,max}-x_{i,min})/2\end{split}  \tag{22}\]The function to be approximated is written by substituting Equations (22) into (21) using an independent variable (as performed before):  
\[f(\xi_{1},\xi_{2},...,\xi_{n})=f(\xi_{i})=g(x_{i,min}+(\xi_{i}+1)(x_{i,max}-x_{i, min})/2) \tag{23}\]  
Similarly, for discrete measurements, a set of sample points (nodes) is introduced: \(\{\xi_{i,0},\xi_{i,i},...,\xi_{i,Mx_{i}};M_{xi}\geq N_{xi}\}\). This produces the following residual approximation error at each measurement node:  
\[r_{i_{k}}=f(\xi_{1_{k}},\xi_{2_{k}},...,\xi_{n_{k}})-\sum_{a_{1}=0}^{N_{x_{1}}} \sum_{a_{2}=0}^{N_{x_{2}}}\cdots\sum_{a_{n}=0}^{N_{x_{n}}}a_{a_{1}a_{1}...a_{n }}\phi_{a_{1}a_{1}...a_{n}}(\xi_{1_{k}},\xi_{2_{k}},...,\xi_{n_{k}}) \tag{24}\]  
where the \(k\) indices represent the measurement for the \(i\)th variable. Equation (24) can written in vector-matrix notation, \(r=f-\Phi a\). Finally, as mentioned in the two- and three-dimensional examples, the weighted least-squares solution can be written as the Kronecker product of \(n\) smaller least-squares matrices:  
\[a=(\Phi^{T}W\Phi)^{-1}\Phi^{T}Wf\equiv\mathcal{C}_{x1}\otimes\mathcal{C}_{x2} \otimes...\otimes\mathcal{C}_{xn}^{-}W^{1/2}f \tag{25}\]  
where \(\mathcal{C}_{xi}=C_{xi}W_{xi}^{-1/2}\), \(W=W_{x1}\otimes W_{x2}\otimes...\otimes W_{xn}\), and \(a=vec\{a_{1k,2k,...,nk}\}\), and \(f=vec\{f_{1k,2k,...,nk}\}\).