The Ritz method for solving partial differential equations using number-theoretic grids

Consider the problem
\begin{gather*} L u(\vec x) = f(\vec x), \\ u(\vec x)\big|_{\partial {G_s}}\big.=g(\vec x), \end{gather*}
where $f(\vec x), g(\vec x) \in E_s^{\alpha}$, $L$ is a linear differential operator with constant coefficients, $G_s$ is the unit cube $[0; 1]^s$.
Its solution is reduced to finding the minimum of the functional
\begin{equation*} v(u(\vec x)) =\underset{G_s}{\int\ldots\int} F\left(\vec x, u, u_{x_1}, \ldots, u_{x_s}\right) dx_1\ ldots dx_s \end{equation*}
under given boundary conditions.
The values of the functional $v(u(\vec x))$ in the Ritz method are considered not on the set of all admissible functions $u(\vec x)$, but on linear combinations
$$ u(\vec x) = W_0(\vec x) + \sum_{k=1}^{n}w_kW_k(\vec x), $$
where $W_k(\vec x)$ are some basic functions that we will find using number-theoretic interpolation, and $W_0(\vec x)$ is a function that satisfies the given boundary conditions, and the rest $W_k( \vec x)$ satisfy homogeneous boundary conditions.
On these polynomials, this functional turns into a function $\varphi (\vec w)$ of the coefficients $w_1, \ldots, w_n$. These coefficients are chosen so that the function $\varphi (\vec w)$ reaches an extremum. Under some restrictions on the functional $v(u(\vec x))$ and the basis functions $W_k(\vec x)$, we obtain an approximate solution of the boundary value problem.