Comments (3)
I'm in favor of all the proposed changes because the CPU/memory costs are trivial and I think it's important to have a consistent, easy to understand convention for users. Once we have that, if there's interest, we can work towards exposing more customizability, but it's important to get the basics right and to make the user-facing API as simple as reasonably possible.
It's very important that we get feedback now, especially from users who oppose these changes. Ultimately, a preference for on-grid vs off-grid boundaries is a somewhat subjective tradeoff and there is no globally optimal choice. We need to make a decision one way or the other for consistency, but there are two ways this can go so if you're a user of this package and prefer a different set of conventions, please speak up now.
from diffeqoperators.jl.
Here's a visualization of the different choices and the effect on the operator:
Dirichlet0: Boundaries off the array (special case). On-grid
A = [2 1 0 0 0
1 2 1 0 0
0 1 2 1 0
0 0 1 2 1
0 0 0 1 2]
Intuition: The "1" off the grid is multiplied by the BC (0) so it's safely ignored.
Dirichlet: Boundaries off the array. On-grid
A = [2 1 0 0 0
1 2 1 0 0
0 1 2 1 0
0 0 1 2 1
0 0 0 1 2]
+ B which is zeros(u) with B[1] and B[end] the BCs
Intuition: The BC is off the grid so we need to add it to the end points for the "1" that was off the grid.
Dirichet: Boundaries on the array. On-grid
A = [1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1]
Intuition: The BC is in the array, so `A*u` with just the stencil works.
Dirichet: Boundaries off the array. Center-grid (constant extrapolation)
Doesn't work.
A = [1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1]
Seems like it would, but it doesn't actually have the BC in the calculation so
it would be unstable.
Neumann0: Boundaries on the array. On-grid
A = [2 2 0 0 0 0 0
1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1
0 0 0 0 0 2 2]
Intuition: Reflect the "1" the falls off over the boundary to get "2 2"
Neumann0: Boundaries off the array. Center-grid
A = [2 2 0 0 0 0 0
1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1
0 0 0 0 0 2 2]
Intuition: Reflect the "1" the falls off over the boundary to get "2 2"
Periodic: Boundaries on the array. Center-grid
A = [2 1 0 0 0 0 1
1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1
1 0 0 0 0 1 2]
Intuition: The "1" wraps around to the other side
Periodic: Boundaries on the array. On-grid
A = [2 1 0 0 0 1 0
1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1
0 1 0 0 0 1 2]
Intuition: The "1" wraps around to the other side, but u[1]=u[end] so you need
to skip it.
Periodic: Boundaries on the left side of the array. On-grid
A = [2 1 0 0 0 0 1
1 2 1 0 0 0 0
0 1 2 1 0 0 0
0 0 1 2 1 0 0
0 0 0 1 2 1 0
0 0 0 0 1 2 1
1 0 0 0 0 1 2]
Intuition: The "1" wraps around to the other side.
Some points to note:
- For Neumann, on-grid and center grid is always the same. Neumann derivatives always need "hidden values" for the derivatives.
- For Dirichlet, center grid doesn't work.
- For Dirichlet with non-zero BC with BCs off the array, the operation is not a matrix multiplication.
- The memory difference having the boundaries on and off is just
2
for the arrays. Since the actual stencils we are using are lazy, there's no cost. So here the difference shown is 5x5 vs 7x7, but that's pseudo because we don't make the matrix! - Dirichlet with values on grid is no longer tridiagonal (but still banded)
- Periodic with bcs on grid is a little weird looking but all of them are not banded, so it doesn't really hurt the sparsity pattern.
As you can see from the visualization, Dirichlet with values off is the "standard" version people use, Neumann under either interpretation is "standard", periodic with a center discretization is "standard".
from diffeqoperators.jl.
Closing this because it's incorrect. Solving for points on the boundary will make the derivative numerically unstable, so we instead want to do SciML/DifferentialEquations.jl#260 (comment)
from diffeqoperators.jl.
Related Issues (20)
- Performance issues with nonlinear_diffusion! HOT 4
- Fix symbolic arrays test to allow MTK v5.21.0 HOT 7
- UndefVarError: dereference not defined when precompiling DiffEqOperators (Win10, julia 1.6.1-2) HOT 7
- Version Compatibility? HOT 1
- Error defining second-order differential HOT 8
- The latest versions of LoopVectorization (0.12.80-81) are not compatible and cause errors HOT 1
- Noob Help: Simple, Linear, Inhomogenous PDE HOT 3
- concretization of High Dimension PDE HOT 2
- LoopVectorization break HOT 5
- Unexpected concretization of Laplace operator HOT 3
- v4.35.0 failed to release HOT 2
- Feature request: Support Summation-By-Parts operators HOT 4
- Support Unitful HOT 6
- Incorrect boundary padded vector with composed PeriodicBC on 2d data HOT 6
- Implementing Drift-Diffusion model in higher dimensions HOT 3
- Uninitialized field `opnorm` in JacVecOperator HOT 4
- SplitODEProblem broken for MatrixFreeOperator HOT 3
- ERROR: LoadError: Some tests did not pass: 30 passed, 0 failed, 4 errored, 4 broken. HOT 18
- Warning related to "LoopVectorization.check_args"
- Lots of precompilation errors and noise HOT 7
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from diffeqoperators.jl.