Calculating the jacobian of norm and square root terms in the Finite Element Method

Question

In the code that my group is writing (Lethe) we use a stabilized approach to solve the Navier-Stokes equation. The GLS stabilized method we use has a stabilisation term which contains a stabilization parameter such as: $$\tau =\left[\left( \frac{2\rho |\mathbf{u}|}{h} \right)^2+\left( \frac{4 \mu}{\gamma h^2} \right) \right]^{0.5}$$

where $\gamma$ is a constant that depends on the type of element, $\rho$ is the density, $h$ is the element size and $\mu$ is the dynamic viscosity.

This stabilization parameter does not depend on the flow variable except for the first term which contains the absolute value of the velocity. I find myself unable to calculate the jacobian of that term, so right now we assemble an analytical jacobian for everything but this term.

Is there anyway to proceed to establish the jacobian of this term analytically? I understand that the $u^2$ term can be treated straightforwardly, but I am unsure on how to calculate the Frechet derivative of the square root.

Answer 1

You can arrive at the Jacobian analytically it just takes a few steps

So assuming we have our typical FE field values:

$$ u_i = \sum_j \phi^j u_i^j $$

Where $i$ represents coordinate direction, $\phi$ our shape function, and $j$ index to DOF.

Start with our stabilization parameter, assuming you are using Euclidian norm:

$$ \tau = \sqrt{\alpha \lVert u \rVert + \beta} $$

Define some functions for ease of writing derivatives

$$ g(\mathbf{u}) = \alpha f(\mathbf{u}) + \beta $$ $$ f(\mathbf{u}) = \lVert u \rVert = \sqrt{u_x^2 + u_y^2 + u_z^2} $$

So now stabilization parameter can be written as

$$ \tau = \sqrt{g(\mathbf{u})} $$

And our derivatives will be in the form

$$ \frac{\partial\tau}{\partial u_i} = \frac{\partial g / \partial u_i}{2 \tau}$$

And so we find derivative of $g$

$$ \frac{\partial g}{ \partial u_i} = \alpha \frac{\partial f}{\partial u_i}$$

and finally derivative of $f$

$$ \frac{\partial f}{ \partial u_i} = \frac{u_i \phi^j}{f(\mathbf{u})} $$

Putting it all together:

$$ \frac{\partial\tau}{\partial u_i} = \frac{\alpha u_i \phi^j}{2 \tau \lVert u \rVert}$$

You may want to double-check my math but it should give you the idea of how to do it.