Compressible DNS

273 views
Skip to first unread message

Tomás Chor

unread,
May 27, 2019, 10:49:24 PM5/27/19
to dedalu...@googlegroups.com
Hello, all, I've been trying to run a compressible DNS case but I keep getting the error "RuntimeError: Factor is exactly singular". 

The equations aren't complete and I'm putting the simplest boundary and initial conditions as possible just to see if I can get this code running first. Has anyone faced a similar problem?

Also, I have based myself off of other example codes I found in the gallery, but I haven't seen any compressible simulation examples, so I had to improvise. 
So if anyone has any working code for simulating a compressible flow I'd appreciate it, since I can start from there. 

My code is attached if anyone care to see it.

Cheers
Tomas.
rb4.py

Daniel Lecoanet

unread,
May 27, 2019, 11:30:20 PM5/27/19
to dedalu...@googlegroups.com
Hi Tomas,

There are several different ways to implement the compressible equations in Dedalus. One way I've found success with is detailed in the appendices of this paper:


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 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/CABUFBtDjejqQ0XXzhbUefQRKnfmqacheQ8%2Bs2yqhgeKkPthW-g%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.

toma...@gmail.com

unread,
May 28, 2019, 12:04:05 PM5/28/19
to Dedalus Users
Hi, Daniel, thanks for the reply.

Yes, I'd already seen that paper. But I'm not sure how to reproduce that in Dedalus. We are trying to model equations (11)-(15), right? I see no good way of implementing them in Dedalus without having two nonlinear terms multiplying in the heat diffusion term (which should be in LHS). Basically here's how I'm implementing those equations:

problem.add_equation("P = (50 + ρ0*g*z)*(1+γ)*(S/Cp + (ρ-ρ0)/ρ0)")
problem.add_equation("T = P / R*ρ")

problem.add_equation("dt(ρ) = - dz(w*ρ) - dx(u*ρ)")

problem.add_equation("dt(u) - ν*dx(Πxx) - ν*dz(Πxz) = - u*dx(u) - w*dz(u) + dx(P) + ν*Πxx*dx(ρ)/ρ + ν*Πxz*dz(ρ)/ρ")

problem.add_equation("dt(w) + ν*dz(Πzx) - ν*dx(Πzz) = - u*dx(w) - w*dz(w) + dz(P) + ν*Πzx*dx(ρ)/ρ + ν*Πzz*dz(ρ)/ρ")

problem.add_equation("dt(S) - κ*(dz(Tz) + dx(dx(T)))/T = - u*dx(T) - w*dz(T) + (κ/T)*(Tz*dz(ρ)/ρ + dx(T)*dx(ρ)/ρ) + (ν/T)*(Πzz*wz + Πxx*dx(u) + Πxz*uz + Πxz*dx(w))")

But I get an error because T can't be multiplying on the LHS. Whatever other workaround I used to get rid of T on the LHS led to other errors.

What would be the best way to model those equations?

Thanks!

Daniel Lecoanet

unread,
May 28, 2019, 1:56:47 PM5/28/19
to dedalu...@googlegroups.com
Hi Tomas,

You can use equations D1-D8 in the 2014 paper I linked to before. This isn't the only way to solve the fully compressible equations, but it seems to work reasonably well.

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

Ben Brown

unread,
May 28, 2019, 2:03:19 PM5/28/19
to dedalus-users
Tomas,
     Evan Anders did some work with these equations (the ones Daniel is referring to) in these two open-access papers:

"Convective heat transport in stratified atmospheres at low and high Mach number"

and

"Predicting the Rossby number in convective experiments"

They worked really well in a variety of settings.
--Ben

toma...@gmail.com

unread,
May 29, 2019, 3:16:40 PM5/29/19
to Dedalus Users
Hi, Daniel, thanks again for the response.

I apologize, I had missed the appendix describing the equations. That was very helpful. However, I'm having trouble with something that I suspect it's silly, which is how to define the mean fields. Here's how I've tried doing it:

problem.parameters['T0'] = T0
problem.parameters['ρ0'] = ρ0
problem.parameters['Lz'] = Lz
problem.parameters['H'] = H
problem.substitutions["T_mean"] = "T0*(Lz+H-z)/H"

I've also tried doing it just like this: `problem.substitutions["T_mean"] = "T0*(Lz+H-z)/H"`

Both of these get accepted by the code without any errors, but every time I try to define equation D1 of the paper like this:

problem.add_equation("dt(w) + dz(T) + T_mean*dz(Y) + T*dz(lnρ) - ν*(D1p1) \
= -T*dz(Y) - u*dx(w) - w*wz + ν*(D1p2)")

I get this error: "NameError: name 'T_mean' is not defined"

Am I missing something here?

Thank you very much!

Daniel Lecoanet

unread,
May 29, 2019, 3:33:41 PM5/29/19
to dedalu...@googlegroups.com
Hi Tomas,

Hmm... That's very strange. Looks like it should work to me. Could you send your script so I can try running it myself?

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

Tomás Chor

unread,
May 29, 2019, 3:38:50 PM5/29/19
to dedalu...@googlegroups.com
Yeah, sure. Here's the code attached.

Thanks for the help!

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/TrBYJqfhUeg/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dedalus-user...@googlegroups.com.

To post to this group, send email to dedalu...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.


--
Tomás L. Chor
Atmospheric and Oceanic Sciences graduate program
University of California, Los Angeles
rb5.py

Daniel Lecoanet

unread,
May 29, 2019, 3:59:26 PM5/29/19
to dedalu...@googlegroups.com, Keaton Burns
Hi all,

There seems to be some sort of bug in the equation parser. For some reason:

problem.add_equation("ρ = exp(lnρ)")

seems to be causing problems. Keaton -- any thoughts on what's going on?

If you put this past equation 8, (uz=dz(u)), the code understands what T_mean is.

However, there are still error messages because I think lnρ is supposed to be the log of the background density, not the log of the perturbation density.

Daniel

On Wed, May 29, 2019 at 3:39 PM Tomás Chor <toma...@gmail.com> wrote:
Yeah, sure. Here's the code attached.

Thanks for the help!

On Wed, May 29, 2019 at 12:33 PM Daniel Lecoanet <leco...@princeton.edu> wrote:
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/TrBYJqfhUeg/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dedalus-user...@googlegroups.com.

To post to this group, send email to dedalu...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.
--
Tomás L. Chor
Atmospheric and Oceanic Sciences graduate program
University of California, Los Angeles

--
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.

Keaton Burns

unread,
May 29, 2019, 4:05:09 PM5/29/19
to dedalu...@googlegroups.com, Daniel Lecoanet
Try rearranging things so that all parameters and substitutions are added before any equations.  I think this is necessary, but unfortunately we don’t have a good system for indicating this.  Might be good to add to the issue tracker.

Tomás Chor

unread,
May 31, 2019, 4:21:25 PM5/31/19
to dedalu...@googlegroups.com
Thank you both Daniel and Keaton. Changing the order of the commands I was able to run with what I think are the correct BC and IC. However, the simulation keeps blowing up after about 20 or 30 iterations. 
I've tried messing with the CFL condition both I can't figure out how to decrease the CFL limit. What are the recommended approaches? I'm pasting my BCs below and I'm also attaching the full code in case anyone wants to take a look at it.

Thanks!

BC:

problem.add_bc("left(Y) = 1")
problem.add_bc("left(dz(T)) = 0")
problem.add_bc("right(dz(T)) = 0")
problem.add_bc("left(u) = 0")
problem.add_bc("right(u) = 0")
problem.add_bc("left(w) = 0")
problem.add_bc("right(w) = 0")


For more options, visit https://groups.google.com/d/optout.
RB.py

Daniel Lecoanet

unread,
May 31, 2019, 4:37:42 PM5/31/19
to dedalu...@googlegroups.com
Hi Tomas,

Sorry you've uncovered a bug in the imposition of boundary conditions! Your system should only have 6 BC's, but Dedalus thinks you need seven. You can fix this by changing equation D3 to:

problem.add_equation("dt(Y) + w*dz(lnρ) + dx(u) + wz \
                     = - u*dx(Y) - w*dz(Y)",tau=False)

If you remove the Y boundary condition, it seems to run ok at least to 52 timesteps (when the simulation ends).

Daniel

On Fri, May 31, 2019 at 4:22 PM Tomás Chor <toma...@gmail.com> wrote:
Thank you both Daniel and Keaton. Changing the order of the commands I was able to run with what I think are the correct BC and IC. However, the simulation keeps blowing up after about 20 or 30 iterations. 
I've tried messing with the CFL condition both I can't figure out how to decrease the CFL limit. What are the recommended approaches? I'm pasting my BCs below and I'm also attaching the full code in case anyone wants to take a look at it.

Thanks!

BC:

problem.add_bc("left(Y) = 1")
problem.add_bc("left(dz(T)) = 0")
problem.add_bc("right(dz(T)) = 0")
problem.add_bc("left(u) = 0")
problem.add_bc("right(u) = 0")
problem.add_bc("left(w) = 0")
problem.add_bc("right(w) = 0")

On Wed, May 29, 2019 at 1:05 PM Keaton Burns <keaton...@gmail.com> wrote:

For more options, visit https://groups.google.com/d/optout.


--
Tomás L. Chor
Atmospheric and Oceanic Sciences graduate program
University of California, Los Angeles

--
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.

Tomás Chor

unread,
May 31, 2019, 4:59:25 PM5/31/19
to dedalu...@googlegroups.com
Hi, Daniel, thanks again. I can't use the keyword here though. According to the docs the only keyword that problem.add_bc accepts is condition. But if I set condition=False I get this error

TypeError: eval() arg 1 must be a string, bytes or code object

in this line

solver = problem.build_solver(de.timesteppers.RK443)

Conda lists my version of Dedalus as 2.1810, but I can't figure out if this is the latest version. I re-ran the script install_conda.sh just in case and I ended up with the same version.


For more options, visit https://groups.google.com/d/optout.

Keaton Burns

unread,
May 31, 2019, 6:16:47 PM5/31/19
to dedalu...@googlegroups.com
Hi Tomas,

To grab the latest version of Dedalus, activate the conda environment, then do:

conda uninstall dedalus
conda install --pre dedalus

to get the latest prerelease version, which includes the new tau option.

Best,
-Keaton

Tomás Chor

unread,
Jun 3, 2019, 4:38:57 PM6/3/19
to dedalu...@googlegroups.com

Hi, Keaton, thanks for the tip.
However, for some reason (I’m guessing because of the type of installation contained in the conda_install.sh script) uninstalling and installing Dedalus by that process doesn’t work.

After issuing conda uninstall dedalus, it fails with this error:

Collecting package metadata: done
Solving environment: failed

PackagesNotFoundError: The following packages are missing from the target environment:
  - dedalus

And I can confirm that conda is set-up correctly, since if I try to delete any other package in the environment, it works without any problems (I’ve tried with several packages).

In any case, if I manually remove the Dedalus package from the dedalus environment and try to issue your command conda install --pre dedalus it fails with EnvironmentLocationNotFound: Not a conda environment: /home/tc/dedalus. So it seems like this command excepts a location with the environment already set-up.

I also tried to find instructions online about the --pre flag but came up empty. All my other attempts failed!

I’d appreciate if you could double-check to see if those are indeed the correct commands. I’d like to try to avoid compiling from source if possible (I’ve tried a couple of times already but unsuccessfully).

Thanks!


For more options, visit https://groups.google.com/d/optout.

Keaton Burns

unread,
Jun 3, 2019, 5:37:03 PM6/3/19
to dedalu...@googlegroups.com
Sorry, those are supposed to be pip commands instead of conda commands. 

Tomás Chor

unread,
Jun 6, 2019, 1:18:57 AM6/6/19
to dedalu...@googlegroups.com
Hi, Keaton, thanks again for the fast reply. However, I'm still having making this work for some reason. So, if  `run pip install --pre dedalus` I get this error:

```
  Installing build dependencies ... done
  Getting requirements to build wheel ... error
  ERROR: Complete output from command /usr/bin/python3 /home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py get_requires_for_build_wheel /tmp/tmpf80az2ya:
  ERROR: Traceback (most recent call last):
    File "/home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py", line 207, in <module>
      main()
    File "/home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py", line 197, in main
      json_out['return_val'] = hook(**hook_input['kwargs'])
    File "/home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py", line 48, in get_requires_for_build_wheel
      backend = _build_backend()
    File "/home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py", line 39, in _build_backend
      obj = getattr(obj, path_part)
  AttributeError: module 'setuptools.build_meta' has no attribute '__legacy__'
  ----------------------------------------
ERROR: Command "/usr/bin/python3 /home/tomaschor/.local/lib/python3.6/site-packages/pip/_vendor/pep517/_in_process.py get_requires_for_build_wheel /tmp/tmpf80az2ya" failed with error code 1 in /tmp/pip-install-3b6pu3fr/dedalus
```

I looked it up and one way to bypass this error is specifying the flag --no-use-pep517 which appears to work, since Dedalus is installed without any errors. However, running the script I then get this error:

ModuleNotFoundError: No module named 'mpi4py'

Which weird since mpi4py appears in the list of packages in conda. I also confirmed that the python I'm using is the miniconda one. I have no idea why this is happening. I understand some of these appear to be bugs/limitations of setuptools or pip, but do you have any idea here? Maybe another way to run the simulation without having to install the pre-release version?

Thanks!



For more options, visit https://groups.google.com/d/optout.

Keaton Burns

unread,
Jun 6, 2019, 1:11:57 PM6/6/19
to dedalu...@googlegroups.com
Hi Tomas,

Based on the path in the error messages, it looks like the “pip” your using is connected to the system python and not the conda environment’s python.  Make sure you’ve activated the dedalus conda environment first.  You can double-check you’re hitting the right python by doing “which python3” and “which pip”, which should both return paths in the dedalus conda env.

Best,
-Keaton

toma...@gmail.com

unread,
Jun 6, 2019, 1:26:43 PM6/6/19
to Dedalus Users
Hi, Keaton, I understand your concern, I had already made sure of that and the problem was another one. In any case, I was able to uninstall Dedalus successfully and install it again using the `--pre` flag. However, it still installed the same version as before I think. Conda lists it as

dedalus 2.1810

And when I try to run it I get the same error as before: not recognizing the `tau` argument.

I'm going to work on this a bit more (maybe even try to compile it from source), but do you guys have any ETA for pushing the latest version as stable? If it's not too far off, maybe it's worth for me to simply wait for it.

Cheers

Tomás Chor

unread,
Jun 6, 2019, 4:20:49 PM6/6/19
to Dedalus Users

Okay, I finally figured out what was going on. Although using pip referred to the correct pip in the Dedalus environment. Using pip install for some reason was referring back to the system-wide installation.
Now I have successfully installed version 2.1905a1 of Dedalus. However, it appears this version also has a bug since I can’t import it. I get the following error:

~/repos/cturbpy/RB.py in <module>                                                                                                                                                                                                             
      9 import time                                                                                                                                                                                                                           
     10                                                                                                                                                                                                                                       
---> 11 from dedalus import public as de                                                                                                                                                                                                      
     12 from dedalus.extras import flow_tools                                                                                                                                                                                                 
     13                                                                                                                                                                                                                                       

~/miniconda3/envs/dedalus/lib/python3.7/site-packages/dedalus/public.py in <module>                                                                                                                                                           
      7 from .tools import logging as logging_setup                                                                                                                                                                                           
      8                                                                                                                                                                                                                                       
----> 9 from .core import future as _future                                                                                                                                                                                                   
     10 from .core.domain import Domain                                                                                                                                                                                                       
     11 from .core import operators                                                                                                                                                                                                           

~/miniconda3/envs/dedalus/lib/python3.7/site-packages/dedalus/core/future.py in <module>                                                                                                                                                      
      6 from functools import partial                                                                                                                                                                                                         
      7                                                                                                                                                                                                                                       
----> 8 from .field import Operand, Data, Scalar, Array, Field                                                                                                                                                                                
      9 from .metadata import Metadata                                                                                                                                                                                                        
     10 from ..tools.general import OrderedSet                                                                                                                                                                                                

~/miniconda3/envs/dedalus/lib/python3.7/site-packages/dedalus/core/field.py in <module>                                                                                                                                                       
     24                                                                                                                                                                                                                                       
     25 # Load config options                                                                                                                                                                                                                 
---> 26 permc_spec = config['linear algebra']['permc_spec']                                                                                                                                                                                   
     27 use_umfpack = config['linear algebra'].getboolean('use_umfpack')                                                                                                                                                                      
     28                                                                                                                                                                                                                                       

~/miniconda3/envs/dedalus/lib/python3.7/configparser.py in __getitem__(self, key)                                                                                                                                                             
   1249     def __getitem__(self, key):                                                                                                                                                                                                       
   1250         if not self._parser.has_option(self._name, key):                                                                                                                                                                              
-> 1251             raise KeyError(key)                                                                                                                                                                                                       
   1252         return self._parser.get(self._name, key)                                                                                                                                                                                      
   1253                                                                                                                                                                                                                                       

KeyError: 'permc_spec'

Any ideas?

thanks

Evan Henry Anders

unread,
Jun 6, 2019, 4:30:24 PM6/6/19
to dedalus-users
Hi everyone,

I just want to chime in to say that I also ran into the 'permc_spec' key error myself when I just updated from an old version of dedalus to a new version. Whether or not it's wise, a workaround is to copy the dedalus.cfg file from the dedalus source (dedalus/dedalus.cfg) to your working directory, and then add these lines:

    # Column permutation option for scipy.sparse.linalg.spsolve
    permc_spec = NATURAL

to the [linear algebra] section (around lines 69-75ish). These used to be in dedalus.cfg by default. If Keaton or anyone can offer some insight into why these lines were removed and what effects adding this back into the .cfg file would have on performance, etc., that'd be great.

Either way, this fixed this bug for me. You could also fix this at the package level (edit dedalus.cfg in the dedalus/ folder), but I like to leave stuff there as untouched as possible.

Best,
Evan

--
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.

Jeffrey S. Oishi

unread,
Jun 6, 2019, 4:34:38 PM6/6/19
to dedalus-users
Hi Tomas,

Yup, that was a total bug. I've fixed it in the latest version of Dedalus. Evan's workaround will also work.

Jeff

Tomás Chor

unread,
Jun 6, 2019, 5:08:05 PM6/6/19
to Dedalus Users
Thanks, guys, indeed this workaround appears to work. I still have to check the results carefully to see if I'm actually reproducing the results correctly, but at least it's running.

Thanks!
Cheers

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/TrBYJqfhUeg/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dedalus-user...@googlegroups.com.

To post to this group, send email to dedalu...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.


--

Jeffrey S. Oishi

unread,
Jun 7, 2019, 11:21:47 AM6/7/19
to dedalus-users
Hi Evan et al,

On Thu, Jun 6, 2019 at 4:30 PM Evan Henry Anders <Evan....@colorado.edu> wrote
to the [linear algebra] section (around lines 69-75ish). These used to be in dedalus.cfg by default. If Keaton or anyone can offer some insight into why these lines were removed and what effects adding this back into the .cfg file would have on performance, etc., that'd be great.

We have recently introduced a new abstraction inside Dedalus that will allow us to easily add different linear algebra solvers and have a unified interface within Dedalus. Keaton and Geoff have been working on new QR solvers for some of our most important linear algebra operations, and having a streamlined way to ensure that these (and other external libraries) can be implemented quickly is very important. Previously, we had the options you saw in the dedalus.cfg file for a two specific linear algebra libraries (UMFPACK and SuperLU). Adding specific options in the cfg file for arbitrarily many libraries becomes a maintenance issue, so we have been moving to a new system.

Unfortunately, in a few places, we forgot to update the code to remove references to the old config file options. There is a pull request on bitbucket (https://bitbucket.org/dedalus-project/dedalus/pull-requests/59/bug-fix-for-issue-60/diff) that should fix the last of these issues. 

To answer your question about performance, adding these lines back in as a workaround has no effect on the performance; in the majority of the code they do nothing (including the file Tomas originally asked about, where the config file was being read, but the options were never being used). In the sparse eigenvalue solver, they will revert to the old default. 

Jeff

Keaton Burns

unread,
Jun 7, 2019, 12:05:16 PM6/7/19
to dedalu...@googlegroups.com
Thanks for tracking this down, Jeff.  I merged your PR and pushed a new prerelease to pypi.

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

Evan Henry Anders

unread,
Jun 7, 2019, 12:13:24 PM6/7/19
to dedalus-users
Thanks for the detailed response, and for patching the bug, Jeff! Appreciate it. That all makes sense, and I'm all in favor of fewer options that I don't really understand in the config file :). Looking forward to the new solvers, too.

Best,
Evan

Reply all
Reply to author
Forward
0 new messages