Issue with Helmholtz (Hodge) decomposition

TurbJet · February 26, 2019, 17:28

Greetings,

I am doing decaying compressible homogeneous isotropic turbulence. I need Helmholtz decomposition to separate the solenoidal and dilatational velocity fields. This is the decomposition I am following (according to Wikipedia)
$\textbf{u} = \nabla\times\textbf{A} - \nabla\Phi$
where A is a vector field, and Phi is a scalar field.

I wrote my own code in MATLAB. Here is a pseudo-code for solving the Poisson equation resulting from the Helmholtz decomposition

Code:

Here, kx, ky, kz are wave numbers in each direction. Also, I use 4th-order difference for evaluating gradient and divergence: central difference for interior grid points, one-sided difference for boundary points.

I apply the code to the solution of the turbulence. The turbulence is solved with 4th-order DRP in space and 4th-order RK in time.

Apparently, the solenoidal velocity field is equal to $\nabla\times\textbf{A}$ , and the dilatational velocity field is equal to $\nabla\Phi$ . However, when I compute the divergence of the solenoidal velocity field in order to verify my code, it doesn't return zero divergence; instead, I have maximum divergence around 1.4, average divergence in the domain is around 1.13e-4, seems like the resulting solenoidal field is not divergence-free.

I went over my code over and over again. As far as I could, I can't find any bug. The finite difference for computing the gradient and divergence is verified, it should not be a problem. The spectral method does not look wrong to me. So I am lost. Could anyone provide some hints or suggestions? Or is there any benchmark case that I can use for verifying my code? I tried Taylor-Green vortex, but this is not challenging enough.

PS:
1. Is it possible that this error is due to the one-sided difference near the boundaries?
2. Is it possible this error is due the the mixture usage of finite difference and spectral method?

Appreciate it!

FMDenaro · February 26, 2019, 18:09

Quote:

Originally Posted by TurbJet

Greetings,

I am doing decaying compressible homogeneous isotropic turbulence. I need Helmholtz decomposition to separate the solenoidal and dilatational velocity fields. This is the decomposition I am following (according to Wikipedia)
$\textbf{u} = \nabla\times\textbf{A} - \nabla\Phi$
where A is a vector field, and Phi is a scalar field.

I wrote my own code in MATLAB. Here is a pseudo-code for solving the Poisson equation resulting from the Helmholtz decomposition

Code:

Here, kx, ky, kz are wave numbers in each direction. Also, I use 4th-order difference for evaluating gradient and divergence: central difference for interior grid points, one-sided difference for boundary points.

I apply the code to the solution of the turbulence. The turbulence is solved with 4th-order DRP in space and 4th-order RK in time.

Apparently, the solenoidal velocity field is equal to $\nabla\times\textbf{A}$ , and the dilatational velocity field is equal to $\nabla\Phi$ . However, when I compute the divergence of the solenoidal velocity field in order to verify my code, it doesn't return zero divergence; instead, I have maximum divergence around 1.4, average divergence in the domain is around 1.13e-4, seems like the resulting solenoidal field is not divergence-free.

I went over my code over and over again. As far as I could, I can't find any bug. The finite difference for computing the gradient and divergence is verified, it should not be a problem. The spectral method does not look wrong to me. So I am lost. Could anyone provide some hints or suggestions? Or is there any benchmark case that I can use for verifying my code? I tried Taylor-Green vortex, but this is not challenging enough.

PS:
1. Is it possible that this error is due to the one-sided difference near the boundaries?
2. Is it possible this error is due the the mixture usage of finite difference and spectral method?

Appreciate it!

There are lot of hidden issues when you start from the HH decomposition in continuous form and then use it after a discretization.

The discrete operators must be congruent and a specific setting of the BCs is required.
First of all, I suggest to use standard second order derivative and write the Div Grad operator in discrete form.

More details can be found in my paper here https://www.researchgate.net/publica...ary_conditions

February 26, 2019, 17:28	Issue with Helmholtz (Hodge) decomposition	#1
TurbJet Senior Member Join Date: Oct 2017 Location: United States Posts: 233 Blog Entries: 1 Rep Power: 10	Greetings, I am doing decaying compressible homogeneous isotropic turbulence. I need Helmholtz decomposition to separate the solenoidal and dilatational velocity fields. This is the decomposition I am following (according to Wikipedia) $\textbf{u} = \nabla\times\textbf{A} - \nabla\Phi$ where A is a vector field, and Phi is a scalar field. I wrote my own code in MATLAB. Here is a pseudo-code for solving the Poisson equation resulting from the Helmholtz decomposition Code: $b = \nabla\cdot\textbf{u}$ $\hat{b} = fft(b)$ $\hat{\Phi} = \frac{\hat{b}}{kx^2 + ky^2 + kz^2}$ $\Phi = ifft(\hat{\Phi})$ $\nabla\times\textbf{A} = \textbf{u} - (-\nabla\Phi)$ Here, kx, ky, kz are wave numbers in each direction. Also, I use 4th-order difference for evaluating gradient and divergence: central difference for interior grid points, one-sided difference for boundary points. I apply the code to the solution of the turbulence. The turbulence is solved with 4th-order DRP in space and 4th-order RK in time. Apparently, the solenoidal velocity field is equal to $\nabla\times\textbf{A}$ , and the dilatational velocity field is equal to $\nabla\Phi$ . However, when I compute the divergence of the solenoidal velocity field in order to verify my code, it doesn't return zero divergence; instead, I have maximum divergence around 1.4, average divergence in the domain is around 1.13e-4, seems like the resulting solenoidal field is not divergence-free. I went over my code over and over again. As far as I could, I can't find any bug. The finite difference for computing the gradient and divergence is verified, it should not be a problem. The spectral method does not look wrong to me. So I am lost. Could anyone provide some hints or suggestions? Or is there any benchmark case that I can use for verifying my code? I tried Taylor-Green vortex, but this is not challenging enough. PS: 1. Is it possible that this error is due to the one-sided difference near the boundaries? 2. Is it possible this error is due the the mixture usage of finite difference and spectral method? Appreciate it!

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Issue on zeroGradient BC for system rotation	magnushaese	OpenFOAM Running, Solving & CFD	2	January 28, 2022 14:43
rhoCentralFoam inviscid airfoil issue (Foam::error::printStack(Foam::Ostream&) at ??)	kmkb21	OpenFOAM Running, Solving & CFD	1	March 1, 2018 02:07
Prism Layer - Core Mesh Boundary Issue	InfernoxCJC	STAR-CCM+	7	November 21, 2017 16:46
potential flows, helmholtz decomposition and other stuffs	pigna	Main CFD Forum	1	October 26, 2017 09:34
foamToTecplot360 issue with multiregion solvers	manuc	OpenFOAM Post-Processing	2	November 21, 2016 14:51