Boundary conditions and SphericalEllProduct in MHD shell simulations

63 views
Skip to first unread message

Till Harke

unread,
Jul 3, 2025, 11:57:18 AMJul 3
to Dedalus Users

Hello,

I’m currently working on an MHD simulation in a spherical shell geometry and am having trouble understanding the correct implementation of boundary conditions.

The formulation using a potential, as presented in Equation 25 of the part 2 methods paper, seems to work well as implemented in Ben Brown’s code found here:
https://github.com/bpbrown/shell_benchmark/blob/main/c2001_case1.py

In particular, the boundary conditions are set like this:

ell_func_o = lambda ell: ell + 1
A_potential_bc_o = de.radial(de.grad(A)(r=Ro)) + de.SphericalEllProduct(A, coords, ell_func_o)(r=Ro) / Ro

ell_func_i = lambda ell: -ell
A_potential_bc_i = de.radial(de.grad(A)(r=Ri)) + de.SphericalEllProduct(A, coords, ell_func_i)(r=Ri) / Ri
...
problem.add_equation("A_potential_bc_i = 0")
problem.add_equation("A_potential_bc_o = 0")

What I’m struggling to understand is the theoretical or mathematical basis of the SphericalEllProduct operator. Specifically, how it operates on each spherical harmonic degree individually, rather than acting on the full gradient of A.

Any clarification would be greatly appreciated.

Thanks in advance,
Till



Daniel Lecoanet

unread,
Jul 3, 2025, 2:04:29 PMJul 3
to dedalu...@googlegroups.com
Hi Till,

SphericalEllProduct multiplies by a specified function of ell+a, where a is the total regularity of each component of a tensor. The external potential boundary condition is

d A_{ell,i}/dr + (ell + a_i + 1) A_{ell,i}/r_o = 0,

where a_i is the total regularity of the ith component of A_{ell, i}, and r_o is the outer radius. The operators act on each ell separately.

Hope that helps,
Daniel

--
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 visit https://groups.google.com/d/msgid/dedalus-users/a6a99d73-d113-4b02-8bb4-92013f0681e4n%40googlegroups.com.

Till Harke

unread,
Jul 4, 2025, 8:14:45 AMJul 4
to Dedalus Users

Thank you for the quick response, that helps a lot.

I'd like to dive a bit deeper into the background, but unfortunately, I haven’t been able to find any documentation on SphericalEllProduct or the term total regularity in the methods papers.

More concretely, I’m trying to implement a constant external field along with a variable internal field that satisfies the appropriate boundary conditions.

Thanks in advance,
Till

Daniel Lecoanet

unread,
Jul 5, 2025, 8:46:24 AMJul 5
to dedalu...@googlegroups.com
Hi Till,

By total regularity I mean the sum over multi-index that appears in various places in section 3 of Vasil et al 2019 (denoted with \bar in the paper). We write tensors as the sum components which are elements of different regularity spaces, e.g., eqn 91 for rank-2 tensors. You can then work out what the potential BC is for an element of a given regularity space.

Daniel

Till Harke

unread,
Jul 9, 2025, 9:34:36 AMJul 9
to Dedalus Users
Hello Daniel,

I have read further into the documentation and got the impression that the regularity components are more of a 'backend' feature, since the connection to the physical space is not obvious for me. For my application I simply need to set an external axial dipole as a boundary condition. As physical BC I would set

dr(A_10) + (l+1)A_10/Ro = 3/2Ro A_ext_10     for l=1 and m=0
dr(A_lm) +(l+1)A_lm/Ro= 0                          otherwise


In dedalus I realized it like this:

ell_func_o = lambda ell: ell+1
A_potential_bc_o = de.radial(de.grad(A)(r=Ro)) + de.SphericalEllProduct(A, coords, ell_func_o)(r=Ro)/Ro

ell_func_e = lambda ell: 1.5 * Ro if ell == 1 else 0
A_boundary = de.SphericalEllProduct(A0, coords, ell_func_e)(r=Ro)

...

problem.add_equation("A_potential_bc_o = A_boundary ")

while I implemented A0 as a Vectorpotential of an axial dipole. Does this makes sense to you? Or is there a simpler way to implement an external axial dipole field at the boundaries?

Thanks for your help!

Till

Daniel Lecoanet

unread,
Jul 9, 2025, 10:19:38 AMJul 9
to dedalu...@googlegroups.com
Hi Till,

Yes, that looks right to me.

Daniel

Reply all
Reply to author
Forward
0 new messages