Computing the nonlinear residual

[maximilianmetti]

This is a nontrivial part of the code, but we can reference this paper by one of the professors at UCSD who's been really successful at solving a simplified version of PNP.

[maximilianmetti]

Our current implementation projects the residual into the finite element space, but should ideally compute the norm of the residual directly (rather than the norm of the projected residual). We want to make the norm of the residual agnostic to the mesh.... if possible

[maximilianmetti]

@arthbous want to leave a little note describing your work towards this end? We should also search for a paper (or look in a book about finite element solutions nonlinear problems) to get appropriate stopping conditions for this problems.

[maximilianmetti]

@arthbous we should compute the nonlinear residual the same way we compute the local entropy:

• project the residual of each equation into a pw constant space (via our mass_lumping_solver viz. https://github.com/thepnpsolver/modular-pnp/blob/master/src/domains.cpp#L83)
• compute the norm of the projection

This is equivalent to defining the (mesh-dependent and diagonal) mass-lumping matrix for L2-projection, M, and multiplying the residual vector r

eqn_residual = transpose(inv(M) * r) * (inv(M) * r) / (transpose(inv(M) * 1) * (inv(M) * 1))

The total residual is sum(eqn_residual).

Algorithmically, we just need to compute

• the block diagonal inv(M) for the mesh
• the mesh-weight matrix, W = inv(M) / sqrt( trace(M*M) )
• r for each equation (the vector of inner-products with the pw constant basis functions)
• for each equation eqn_residual = transpose(r) * W * r
• total_residual = sum(eqn_residual)

I think this will give us the best approximation to the nonlinear residual. I feel like we tried to do this before but didn't quite get it right... Did you try to do this at some point? What are your thoughts on this approach?

The pain in the ass about this is computing r, which is dependent on the equation at hand (i.e. it's not the same expression for the Poisson and Nernst-Planck equations).