dedalus
dedalus copied to clipboard
Laplace equation in 2D Dirichlet BCs
Hi, I am trying to implement the Laplace equation in 2D as follows: d2f/dx2 + d2f/dy2 = 0, with Dirichlet boundary conditions: f(x,1.0) = -1.0, f(0.0,y) = f(1.0,y) = f(x,0.0) = 0.0.
Using the Poisson script already developed for BVPs in the github repo threw the following error:
Traceback (most recent call last):
File "poisson.py", line 26, in <module>
problem.add_bc("u(y='left') = sin(8*x)")
File "/home/pavan/miniconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/problems.py", line 139, in add_bc
return self.add_equation(*args, **kw)
File "/home/pavan/miniconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/problems.py", line 131, in add_equation
self._build_object_forms(temp)
File "/home/pavan/miniconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/problems.py", line 152, in _build_object_forms
temp['LHS'] = field.Operand.parse(temp['raw_LHS'], self.namespace, self.domain)
File "/home/pavan/miniconda3/envs/dedalus/lib/python3.8/site-packages/dedalus/core/field.py", line 102, in parse
expression = eval(string, namespace)
File "<string>", line 1, in <module>
TypeError: 'Field' object is not callable
However, when I adapt the script and specify the boundary conditions as add_equations instead:
import os
import numpy as np
import matplotlib.pyplot as plt
from dedalus import public as de
from dedalus.extras import plot_tools
# Create bases and domain
x_basis = de.Fourier('x', 100, interval=(0, 1))
y_basis = de.Chebyshev('y', 100, interval=(0, 1))
domain = de.Domain([x_basis, y_basis], grid_dtype=np.float64)
# Poisson equation
problem = de.LBVP(domain, variables=['u','uy'])
problem.add_equation("dx(dx(u)) + dy(uy) = 0")
problem.add_equation("uy - dy(u) = 0")
problem.add_equation("left(u) = 0")
problem.add_equation("right(u) = 1")
# Build solver
solver = problem.build_solver()
solver.solve()
# Plot solution
u = solver.state['u']
u.require_grid_space()
plot_tools.plot_bot_2d(u)
plt.savefig('laplace.png')
Which now works and gives the following:
I am now only able to specify the desired boundary conditions on the top and bottom of the domain and it looks like no-flux BCs are enforced on the left and right sides of the domain. I understand that the boundary condition is enforced along the y coordinate since I read somewhere that the Chebyshev basis function is the axis along which the boundary condition is imposed on in Dedalus. In this case, I am not entirely sure how to go about imposing the right Dirichlet BCs to satisfy the problem.
For additional info. I followed the installation procedure on the documentation exactly.
Hi, If I understand correctly, it sounds like you are trying to impose BCs along both axes. This feature is not included in the current Dedalus release, but I do believe it is an upcoming feature.
EDIT: To clarify, on the left and right, periodic boundary conditions are being used (a consequence of it being a Fourier basis), not no-flux boundary conditions.
Ah, so in that case:
- Solving the problem with the specification of Dirichlet boundary conditions on all four domains (i.e. the model I want to solve) is not possible at the moment. Is there any indication when this might be out?
- Is it possible to employ a different basis function to atleast specify Neumann BCs along the left and right boundaries. I am guessing that specifying a Fourier basis results in periodic BCs along that axis?
- What is the nature of the issue running the orignal poisson problem on the github repo. It seems that the
add_bc
formulation is much cleaner than adding additional equations in terms of formulating the problem. It seems that running other example scripts such as the KDV script or 1-D Lane Emden equation had no issues.
Thanks!
I wanted to follow up on this: is it now possible to specify Dirichlet boundary conditions along both axes? If not, is it possible to specify Dirichlet boundary conditions along one axis and Neumann boundary conditions along another?
Is it still the case that Dirichlet boundary conditions can only be imposed on one axis?
Doing this correctly for general PDEs is still an open research problem, but we're working on it (see this recent paper), and there will be a separate release and additional examples when it's supported in the code, so stay tuned!
Just a quick follow-up: Is the feature with such boundary conditions currently not supported at all (in the sense that such codes will not run properly or return undefined behavior) or is it generally possible to set such problems up and only be able to execute them in serial? The poisson script can indeed be modified in this sense but my question is basically: Is the fact that it's working a coincidence?
This issue is that, like all Dedalus models, enforcing boundary conditions requires tau terms, and these get complicated when you have different boundaries meeting at edges/corners. We understand how to do it in simple cases (e.g. Poisson) as described in the paper linked above, but not in more general cases. And we simply can't offer general support when doing it correctly and reliably it's an open/active mathematical problem.
Also yes, even with the correct tau terms for a given equation set, double-Chebyshev problems will in general be fully coupled and non easily parallelized. Parallelizing them would require using preconditioned iterative solvers or distributed direct solvers. These are both interesting routes to explore, but again they are research grade questions in their own right.
Thank you for the clarification. Of course, I understand that this issue, being an open research problem, is beyond the scope of a general purpose library can be expected to handle.