cylindrical geometry

274 views
Skip to first unread message

mg.de...@gmail.com

unread,
Sep 30, 2014, 12:29:55 PM9/30/14
to dedalu...@googlegroups.com
Hello,

after getting dedalus running, I would like to generate a set-up for cylindrical geometry. A finite cylindrical gap or Taylor-Couette system with non-axisymmetric flow. Anybody here, who has done this already?

I found an older script for axisymmetric Taylor Couette flows. What would be a good basis to extend it first to be 3D (Chebyshev/Fourier/Fourier?) and in a second step to have a finite height?

Is there a description of possible basis functions, time stepping methods, boundary conditions and how to apply them? What does left_bc/right_bc mean in a 2D domain? Is there top_bc/bottom_bc as well?

Thank you.

Marcus

Ben Brown

unread,
Sep 30, 2014, 12:59:26 PM9/30/14
to mg.de...@gmail.com, dedalus-users
Marcus,
     Jeff Oishi has a lot of experience with Taylor Couette, so I'll leave him to answer most of your e-mail.

Regarding the boundary conditions: left/right are the same as bottom/top.  The only direction we can impose arbitrary boundary conditions currently is the chebyshev directrion, and left/right refer to the bottom/top of the chebyshevs.

Timesteppers etc. will be documented in the forthcoming instrument paper, but you can see many of the timesteppers we've implemented by following the documented literature trail in dedalus2/pde/timesteppers.py.  You might look in particular at Ascher et al 1995 for the basic implicit/explicit multistep timesteppers (CNAB2, SBDF2, SBDF3, etc.).  The actual implemented variable timestep variants come from Wang & Ruuth 2008, but I find Ascher et al a bit clearer for understanding the underlying methods.
--Ben


--
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 post to this group, send email to dedalu...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/a5a69e01-c32a-4fdd-a102-8d125f11dbe4%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Keaton Burns

unread,
Sep 30, 2014, 1:03:02 PM9/30/14
to mg.de...@gmail.com, dedalu...@googlegroups.com
Hi Marcus,

As with the installation instructions, much of the documentation is forth-coming, but the best place to start is probably with the tutorial notebooks, which can be found here:

http://dedalus-project.readthedocs.org/en/latest/getting_started.html

We currently support N-dimensional domains, where the first (N-1) dimensions are represented by Fourier bases, and the last is represented by a Fourier or Chebyshev basis. Boundary conditions are only applicable to the Chebyshev direction (the Fourier directions assume periodicity), and so left/right BC’s correspond to the start/end of the Chebyshev interval.

You could do vertically periodic, non-axisymmetric Taylor Couette as you describe, with a Fourier-Fourier-Chebyshev domain for (z, θ, r). Here the left/right BC’s are on the inner and outer cylinders (the Chebyshev direction).

You can’t explicitly set top and bottom boundary conditions since z must be a Fourier basis — but there are some trick you can play with explicitly setting the symmetry of the Fourier series in that direction to mimic no-slip or stress-free boundaries on the Fourier interval. However, we’re currently working on a sin/cos basis which will do this sort of thing explicitly (and more efficiently).

-Keaton

Keaton Burns

unread,
Sep 30, 2014, 1:18:10 PM9/30/14
to Ben Brown, mg.de...@gmail.com, dedalus-users
To make sure we’re all on the same page — I think Ben was referring to something like a vertically bounded box here (e.g. for parallel Couette flow) that is horizontal periodic, where the Chebyshev basis is the vertical direction (z), so the “left”/“right” BC’s are setting the conditions on the bottom/top of the box.  

For a cylindrical Taylor-Couette setup, the radial expansion is done in Chebyshev polynomials, so the “left”/“right” BC’s are setting the conditions on the inner and outer cylinders.  Currently, the vertical/axial direction (z) must be represented by a Fourier basis and is considered periodic, with no other explicitly available BC’s.


Ben Brown

unread,
Sep 30, 2014, 1:22:59 PM9/30/14
to Keaton Burns, mg.de...@gmail.com, dedalus-users
Keaton & Marcus,
     Yes, that's what I was referring to.  Thank you for the clearer and more detailed note.

--Ben

Daniel Lecoanet

unread,
Sep 30, 2014, 1:27:40 PM9/30/14
to Keaton Burns, Ben Brown, mg.de...@gmail.com, dedalus-users
Hi Marcus,

To reiterate what has already been mentioned, it would be probably be easiest to first try doing periodic-in-z non-axisymmetric Taylor-Couette (using Fourier-Fourier-Chebyshev), and then afterwards play with the boundary conditions (basically imposing a symmetry across the middle of the z-domain) to do a finite box with stress-free boundary conditions.

Let me know if you have any questions about generalizing from the axisymmetric case to the non-axisymmetric case.  I would be happy to send you a script I have for a 3D simulation in a box if that would help you to go from 2D to 3D.

Daniel

mg.de...@gmail.com

unread,
Sep 30, 2014, 5:13:46 PM9/30/14
to dedalu...@googlegroups.com, keaton...@gmail.com, bpb...@gmail.com
Thank you to all of you for the extensive answers. As soon as one knows that everything except one direction is Fourier, the BCs become much more clear.

Daniel, of course I would appreciate if you would like to share your 3D TC code. In the beginning changing something already existing is usually easier than writing completely new.

Cheers
Marcus

j s oishi

unread,
Sep 30, 2014, 5:24:49 PM9/30/14
to mg.de...@gmail.com, dedalus-users, Keaton Burns, Ben Brown
Hi Marcus.

Actually, I think I'm the one with the 3D TC code. I will send it your way later this week. I'm about to teach in a few minutes, and I can't remember which computer it's on. Daniel has a 3D convection in a Cartesian box that might also be helpful.

Jeff

Ben Brown

unread,
Sep 30, 2014, 6:10:49 PM9/30/14
to mg.de...@gmail.com, dedalus-users, Keaton Burns
Marcus,
     As another heads up: when handling cylindrical geometries, an initial approach including the 1/r geometric terms in radial divergence terms and in the theta or phi (azimuthal) derivatives require a very broad band representation in our matrix structures.  Namely, the "order" parameter in the parsed problems needs to be large (say order=20 or possibly higher).  This would vastly slow down your computations. 

We have found that it is very, very effective to pre-multiply through both the LHS and RHS of the equation by r or r^2 so that there are no 1/r terms anywhere in your equations.  With at most r^2 terms, you can now fully represent the system with order=2 or 3 or so.  I'm a little fuzzy on the exact details, as it's been about 9 months since I've touched similar problems.

Jeff's TC codes have this implemented.
--Ben

mg.de...@gmail.com

unread,
Oct 1, 2014, 6:59:55 AM10/1/14
to dedalu...@googlegroups.com, keaton...@gmail.com
I have seen this in one of the scripts. What does order=20 mean? Has nothing to do with the order of the Chebyshev basis, right?

On Wednesday, October 1, 2014 12:10:49 AM UTC+2, Ben Brown wrote:
> Marcus,
>      As another heads up: when handling cylindrical geometries, an initial approach including the 1/r geometric terms in radial divergence terms and in the theta or phi (azimuthal) derivatives require a very broad band representation in our matrix structures.  Namely, the "order" parameter in the parsed problems needs to be large (say order=20 or possibly higher).  This would vastly slow down your computations. 
> ...

geoff

unread,
Oct 1, 2014, 7:25:24 AM10/1/14
to mg.de...@gmail.com, dedalu...@googlegroups.com, keaton...@gmail.com
I can answer that! It means that the non-constant coefficients get expanded in Chebyshev polynomials up to degree of 19. Yes, 20=19. The reason is that it uses T_{0}, T_{1}, … , T_{19}, which is a set of 20. The default is order=1. There’s a good chance that we’ll change this convention in future versions. Cheers, Geoff
> --
> 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 post to this group, send email to dedalu...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/dedalus-users/5153fb73-02e1-40c2-9dbe-421200350ff5%40googlegroups.com.

j s oishi

unread,
Oct 1, 2014, 8:48:54 AM10/1/14
to mg.de...@gmail.com, Keaton Burns, dedalus-users

Hi Marcus,

That is the number of chebyshev polynomials the code expands any non-constant coefficients in.

For example, if you write the curvature terms for cylindrical geometry as 1/r, and r is your chebyshev direction, you'd need to expand that coefficient. This is technically independent of the number of modes you expand the domain in. However, 1/r takes a lot of coefficients to accurately represent it. This is what Ben was talking about: if instead you multiply the equations through to make the coefficients r on the non curvature terms, you can set the order to 2.

I will send you a notebook explaining this later today.

Jeff

--
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 post to this group, send email to dedalu...@googlegroups.com.

j s oishi

unread,
Oct 17, 2014, 11:53:57 AM10/17/14
to Marcus, dedalus-users
Hi Marcus,

My apologies for the delay in my response. This is my first semester teaching, and I am finding myself chronically behind on everything! I will have a look at your script, and send you my own over the weekend.

Regarding the divB = 0 constraint, as long as you have resistivity, you can enforce it without an additional equation due to Dedalus's first order formulation. Simply set Bz_z = -(dx(Bx) + dy(By)), where Bz_z is the z-derivative of the z-field. Then, in the induction equation for Bz, you have simply the first derivative of Bz_z:

 dt(Bz) - eta*dz(Bz) = [non-linear terms]

This works well in my (admittedly very limited) testing thus far. Does this make sense?

Jeff

On Thu, Oct 16, 2014 at 4:06 AM, Marcus <mg.de...@gmail.com> wrote:
Hello Jeff,

you mentioned your 3D TC script two weeks ago. Meanwhile I extended the 2D version on my own. Not really successful, unfortunately. Solutions might be correct, but the code runs very unstable and crashes after always after few tenth of rotations, independent from resolution or order of expansion of parameters or Reynolds number, even the axisymmetric solutions crash after some time. I can not find a reason. The script is attached. A hint from you (or your script?) would be great. Thank you in advance.

BTW I implemented the induction equation as well and checked the TC dynamo (Willis A&A 2006). Onset looks good, but unfortunately suffers from the same 'instability'. If I get it stably running, it might be a nice extension of the dedalus examples. Beside the crashing code, I have also not yet found a clean way to implement the divB=0 condition, because it will be an additional equation without adding an independent variable. Not sure how to solve it. Do you now if somebody already tried it using dedalus?

Cheers
Marcus

Reply all
Reply to author
Forward
0 new messages