I am implementing a MatLab program to solve the equation given in this paper, which involves solving integrals coming from the variational formulation of the problem. One of them is $$\int_{\partial \Omega} f\left(\sum_{i=1}^Nu_i\phi_i\right) \cdot \phi_j d\Omega$$ where $\phi_j$ is a function such that, in the node $(x_j, y_j)$ of the triangular partition of $\Omega$, its value is 1 and 0 elsewhere. That is, $\{\phi_i\}$ is the standard nodal basis defined for a triangular partition on FEM.
My confusion is how to evaluate this integral. Here is what I think I need to do:
At each iteration, we want to evaluate $u^{(k)} = (u_1^{(k)}, ..., u_N^{(k)})^T$ an approximation of $u$ at each node from our triangular partition with $N$ nodes. This integral is only evaluated at the nodes which lie on the boundary $\partial \Omega$. Suppose, for example, that $(x_1, y_1)$ is one of those nodes and we want the value of $u$ at this node (ie, we want u_1 at some iteration). Then, we should have to evaluate:
$$\int_{\partial \Omega} f\left(\sum_{i=1}^Nu_i\phi_i(x_1, y_1)\right) \cdot \phi_1(x_1, y_1) d\Omega = \int_{\partial \Omega} f(u_1) d\Omega$$
I am really confused if this is correct, and if it is, how do I calculate this integral?
I am sorry if the question is confusing, but the paper can give a clear idea of what should be done, if needed.
Thanks!