Potential problem with curl calculations

[Andrei Igoshev]

Dear Dedalus developers,

I have encountered a weird problem while trying to calculate three consecutive curls in a ball. 

If I define:

A  = dist.VectorField (coords, name='A', bases=ball)

A['g'][0] = r**5*np.sin(theta)

and calculate three consecutive curls numerically, I get a weird jump at r = 1.

J2_old = (d3.curl(d3.curl(d3.curl(A)))).evaluate()

Please find a figure attached. The jump gets much worse when I increase the resolution.

It seems that this problem becomes noticeable when I compute the first curl already.

In this case Br and Btheta seems to coincide with their analytical expressions, but the spectral coefficients describing their decomposition into radial basis diverge following ~N.

Please find the script attached.

[Andrei Igoshev]

According to Calum, the difference in spectra is explained because I have different basis for B and B_a. 

However, the problem with triple curl is still real. Please find an example script attached.

Best regards,

Andrei Igoshev

[Keaton Burns]

Hi Andrei,

The problem is that when you set a field’s values in grid space, then the spectral transform will put ~1e-16 noise in all coefficients.  Then if you start differentiating the field, the noise in mode N will grow by something like a factor of N (for Fourier/spherical harmonics) or even up to N^2 (for Chebyshev). So after enough derivatives, and for high enough N, the high coefficients will become O(1) and contaminate the solution.

The fix, if you need to take many derivatives, is to carefully filter out very small amplitude modes by explicitly setting them to zero at the beginning.  For instance, see the attached script which modifies your triple-curl example in this manner.

Best,

-Keaton

--
You received this message because you are subscribed to the Google Groups "Dedalus Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dedalus-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/aa67dd58-53ae-4dd7-992b-2f9412d756f1n%40googlegroups.com.

[Andrei Igoshev]

Hi Keaton,

Thank you very much! It works indeed.

Cheers,

Andrei

You received this message because you are subscribed to a topic in the Google Groups "Dedalus Users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dedalus-users/_awkZHjIV78/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dedalus-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/CADZXxBjU%2Bg5FmK4QZ4ezb4SWmTNRbL9EdrG4QxUOcVUoj2XGZw%40mail.gmail.com.