Comments (9)
Issues like these come out of the fact that these are discretization operators meant for efficient PDE solving, not necessarily tailored to general numerical differentiation. The math is the same, but the API preferences diverge to some extent. That's why I've always postulated that when all is said and done, we probably want a separate package for numerical differentiation that shares a common core with the calculations currently done here.
The one-sided issue is simple - we should ultimately let the users choose whether they prefer sacrificing symmetry or approximation order, as you can't have both. The rest are things we should discuss in terms of how we'd like the final APIs (probably plural) to look and then structure the code accordingly.
from diffeqoperators.jl.
Dx*F
is x derivatives, F*Dy
is y derivatives. By hand do the calculation for the discretized Laplace and this becomes apparent (or you might have it the other way around, depending on x as rows or columns)
from diffeqoperators.jl.
So you mean I can use A = LinearOperator(1,2,...)
, and then A*F
is x-derivatives and F*A
is y-derivatives? (And similarly for higher-order...)
That's neat, thank you.
from diffeqoperators.jl.
Yes. And that's if x goes down columns and y goes across rows. If that's swapped, just swap to F*A
and A*F
.
from diffeqoperators.jl.
Hah, this was actually obvious the second one thinks linear algebra...
from diffeqoperators.jl.
I've tried this now, and have some problems:
- The one-sided derivatives at the boundaries means that
A
is not fully symmetric. SoF*A
won't work without a proper transposition onA
- The multiplication operator for
A*F
andF*A
calls operations from LinearMaps as opposed to what happens withA*x
forx
a vector.
julia> @which(F*A)
*(A1::AbstractArray{T,2} where T, A2::LinearMaps.AbstractLinearMap) in LinearMaps at /home/asbjorn/.julia/v0.6/LinearMaps/src/wrappedmap.jl:41
- Currently, the results from applying
F*A
do not work:
using PDEOperators
xarr = linspace(0,1,51)
yarr = linspace(0,1,101)
dy = yarr[2]-yarr[1]
F = [x^2+y for x = xarr, y = yarr]
B = LinearOperator{Float64}(2,2,dy,length(yarr),:None,:None)
# PROBLEM: F*B is neither zero nor constant
from diffeqoperators.jl.
@dextorious explained it nicely. If you check the operator that comes out from LaPlace 2,2, it won't have all of the properties you want with those BCs by design. The discretization error is O(dx^2) though, so it does converge to no error. But these are not "structure-preserving" operators: these are efficient discretizations, which are very different. If you're looking for structure-preserving techniques, that's on the wishlist but not coming soon.
The multiplication operator for AF and FA calls operations from LinearMaps as opposed to what happens with A*x for x a vector.
That's a mistake that should be corrected. @shivin9 we need the other operator defined. Since the matrix is symmetric this should be pretty easy.
from diffeqoperators.jl.
@ChrisRackauckas Yeah I saw that. If F
is multidimensional though then should we recursively apply the operator for every lower dimensional structure till the time we reach it's rows?
I will correct it for 2D domains anyways.
from diffeqoperators.jl.
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.