chakravala / grassmann.jl Goto Github PK
View Code? Open in Web Editor NEW⟨Grassmann-Clifford-Hodge⟩ multilinear differential geometric algebra
Home Page: https://grassmann.crucialflow.com
License: GNU Affero General Public License v3.0
⟨Grassmann-Clifford-Hodge⟩ multilinear differential geometric algebra
Home Page: https://grassmann.crucialflow.com
License: GNU Affero General Public License v3.0
There seems to be a bug in the constructor when creating a MultiVector from a blade:
julia> using Grassmann
julia> basis"++"
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> a = v1 + v2
1v₁ + 1v₂
julia> typeof(a)
MBlade{Int64,⟨++⟩,1,2}
julia> MultiVector(a)
ERROR: type MArray has no field v
Stacktrace:
[1] getproperty(::Any, ::Symbol) at ./sysimg.jl:18
[2] MultiVector(::MBlade{Int64,⟨++⟩,1,2}) at .julia/packages/Grassmann/H9Zog/src/multivectors.jl:261
[3] top-level scope at none:0
Thanks for making such a great GA lib available.
How do you visualize the results?
I'm trying to analyze geometric structure in GA.
Can you recommend a workflow or any libraries (if you are already aware of by change) ?
Thanks!
Best,
Hongying
P.S. Previously I use "gable" by Leo Dorst, and also noticed GAlgebra (in Python) and its Julia wrapper).
It is sometimes convenient to call a multiplication operator with a single argument, as in \wedge(xs...)
where xs
is a collection of chains. For example, *
returns its argument when called with a single argument. It would be convenient if \wedge
did the same; this would avoid code such as
if length(xs) == 1
vol = abs(xs[1])
else
vol = abs(∧(xs...))
end
Please see comments in the pdf of your JuliaCon submission.
grassmann-juliacon-2019.pdf
Summary of Comments on grassmann-juliacon-2019.pdf
This paper is a description of Grassmann.jl and its supporting packages, which are appear quite impressive and extensive. Unfortunately, it
assumes extensive knowledge of terminology and notation from the reader, without ever defining what that reader should know, whose
notations they should be familiar with, and what that reader can expect to gain from reading.
As someone who has written a lot of geometric and topological code in low dimensions, and been interested in Clifford algebra enough
to read the book of Dorst, Fontijne, & Mann and some of Hestenes, I was eager to review this. I was disappointed because it does not
clearly say what the package does, nor what decisions were made between alternative implementations. (There are many intriguing hints,
e.g., Due to the design of the VectorBundle type
system, these operations make code optimization easier at compile-time by evaluating the bit parameters.)
As Fred Brooks has said, "Since anyone can create something new [in a synthetic field], that alone does not establish a contribution.
Rather, one must show that the creation is better." (Quoted from https://cra.org/resources/best-practice-memos/evaluating-computerscientists-and-engineers-for-promotion-and-tenure/) A paper should do more than quote from the documentation of a package; it
should provide context so the reader (from the author's desired audience) can understand the decisions made by the package author and
why they are better than alternatives.
According to the documentation of \oslash
, this should hold:
help?> ⊘
"⊘" can be typed by \oslash<tab>
search: ⊘
⊘(ω::TensorAlgebra,η::TensorAlgebra)
Sandwich product: ω⊘η = (~η)⊖ω⊖η
julia> v⊘v1 == (~v1)⊖v⊖v1
false
works on a Multivector and Basis, but doesn't seem to give correct result for Blades.
julia> a
1.0v₁ + 2.0v₂ + 0.0v₃
julia> b
3.0v₁
julia> ab = a*b
3.0 - 6.0v₁₂
julia> reverse(ab)
3.0 + 6.0v₁₂
julia> reverse(ab(2))
6.0v₁ + 0.0v₂ + 0.0v₃
julia> typeof(ab(2))
SBlade{Float64,⟨+++⟩,2,3}
The docstring for >>>
seems to document \oslash
. I therefore originally assumed these are just different names for the same function, but that doesn't seem to be so.
I find that the uparrow and downarrow operators for conformal geometries convert rational numbers to floating-point numbers:
julia> using Grassmann
julia> S = Signature(4, 1, 1)
⟨∞∅++⟩
julia> V = SubManifold(S)
⟨∞∅++⟩
julia> B = Λ(S)
DirectSum.Basis{⟨∞∅++⟩,16}(v, v∞, v∅, v₁, v₂, v∞∅, v∞₁, v∞₂, v∅₁, v∅₂, v₁₂, v∞∅₁, v∞∅₂, v∞₁₂, v∅₁₂, v∞∅₁₂)
julia> Chain{V,1}(1//1,2//1,3//1,4//1)
(1//1)v∞ + (2//1)v∅ + (3//1)v₁ + (4//1)v₂
julia> x = Chain{V,1}(1//1,2//1,3//1,4//1)
(1//1)v∞ + (2//1)v∅ + (3//1)v₁ + (4//1)v₂
julia> ↑(x)
-0.0 + 11.5v∞ + 3.0v∅ + 3.0v₁ + 4.0v₂
julia> ↓(x)
0.0 + 1.5v₁ + 2.0v₂
The problem is caused by the terms G.v∞/2
and inv(G.v∞∅)
, which don't have any type information and thus default to Float64
. Multiplying the basis vectors by one(valuetype(\omega))
before the division/inversion solves the problem.
This breaks my code because I like to use rational numbers for testing to avoid round-off errors. I am currently using a corrected copy of these functions that I copied into my code.
julia> using Grassmann
julia> basis"+++"
(⟨+++⟩, v, v₁, v₂, v₃, v₁₂, v₁₃, v₂₃, v₁₂₃)
julia> v + v1 + v2
Internal error: encountered unexpected error in runtime:
MethodError(f=typeof(Base.string)(), args=(Expr(:call, Symbol("!"), Expr(:., :typ, :(:mutable))),), world=0x0000000000000eb9)
rec_backtrace at /Users/mason/julia/src/stackwalk.c:94
record_backtrace at /Users/mason/julia/src/task.c:217
jl_throw at /Users/mason/julia/src/task.c:417
jl_method_error_bare at /Users/mason/julia/src/gf.c:1649
jl_method_error at /Users/mason/julia/src/gf.c:1667
jl_lookup_generic_ at /Users/mason/julia/src/gf.c:2195
jl_apply_generic at /Users/mason/julia/src/gf.c:2215
lift_leaves at ./compiler/ssair/passes.jl:300
getfield_elim_pass! at ./compiler/ssair/passes.jl:658
run_passes at ./compiler/ssair/driver.jl:123
optimize at ./compiler/optimize.jl:164
typeinf at ./compiler/typeinfer.jl:35
typeinf_ext at ./compiler/typeinfer.jl:576
typeinf_ext at ./compiler/typeinfer.jl:613
jfptr_typeinf_ext_1 at /Users/mason/julia/usr/lib/julia/sys.dylib (unknown line)
jl_apply_generic at /Users/mason/julia/src/gf.c:2219 [inlined]
jl_apply at /Users/mason/julia/src/./julia.h:1571 [inlined]
jl_type_infer at /Users/mason/julia/src/gf.c:277
jl_compile_method_internal at /Users/mason/julia/src/gf.c:1819
jl_fptr_trampoline at /Users/mason/julia/src/gf.c:1863
do_call at /Users/mason/julia/src/interpreter.c:323
eval_stmt_value at /Users/mason/julia/src/interpreter.c:362 [inlined]
eval_body at /Users/mason/julia/src/interpreter.c:759
jl_interpret_toplevel_thunk_callback at /Users/mason/julia/src/interpreter.c:885
Interpreter frame (ip: 0)
Core.CodeInfo(code=Array{Any, (2,)}[
Expr(:call, :+, :v, :v1, :v2),
Expr(:return, SSAValue(1))], codelocs=Array{Int32, (2,)}[1, 1], method_for_inference_limit_heuristics=nothing, ssavaluetypes=2, linetable=Array{Any, (1,)}[Core.LineInfoNode(mod=Main, method=Symbol("top-level scope"), file=:none, line=0, inlined_at=0)], ssaflags=Array{UInt8, (0,)}[], slotflags=Array{UInt8, (0,)}[], slotnames=Array{Any, (0,)}[], inferred=false, inlineable=false, propagate_inbounds=false, pure=false)jl_interpret_toplevel_thunk at /Users/mason/julia/src/interpreter.c:894
jl_toplevel_eval_flex at /Users/mason/julia/src/toplevel.c:764
jl_toplevel_eval at /Users/mason/julia/src/toplevel.c:773 [inlined]
jl_toplevel_eval_in at /Users/mason/julia/src/toplevel.c:793
eval at ./boot.jl:328
eval_user_input at /Users/mason/julia/usr/share/julia/stdlib/v1.1/REPL/src/REPL.jl:85
macro expansion at /Users/mason/julia/usr/share/julia/stdlib/v1.1/REPL/src/REPL.jl:117 [inlined]
#26 at ./task.jl:259
jl_apply at /Users/mason/julia/src/./julia.h:1571 [inlined]
start_task at /Users/mason/julia/src/task.c:572
1 + 1v₁ + 1v₂
I'm trying to implement sqrt
for pure scalar multivectors (so for Basis{V, 0, 0}
, MValue{V, 0, v, T}
and SValue{V, 0, v, T}
), but I have a hard time to find the right syntax for the methods type declaration. It is easy for Basis
:
sqrt(t::Basis{B, 0, 0x0000000000000000}) where {B} = 1
But I'm not sure how to define the type of MValue
and SValue
for scalar basis v
. I can get v
into the scope and write
julia> basis"2"
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> sqrt(x::SValue{B, 0, v, T}) where {B, T} = sqrt(x.v)
sqrt (generic function with 23 methods)
julia> sqrt(13*v)
3.605551275463989
but that seems ugly and does not make sure that v
matches the scalar basis for any B
.
I guess I have to use the function indexbasis
but I'm not sure about its syntax.
So this is more like a question than an actual issue (except maybe a documentation issue).
Btw: What is the difference between SValue and MValue?
PS: I hope that its clear, but I'm saying it anyway: I really like this package, thanks for all the work you are putting into it.
When I add two multivectors A
and B
everything works like expected:
julia> using Grassmann
julia> @basis "2"
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> A = 2v1 + v2
2v₁ + 1v₂
julia> B = v1 + v2
1v₁ + 1v₂
julia> A + B
3v₁ + 2v₂
julia> A
2v₁ + 1v₂
julia> B
1v₁ + 1v₂
However when I add one of the basis vectors (v1
or v2
) that are not orthogonal to the multivector, it changes the multivector:
julia> v1 + A
3v₁ + 1v₂
julia> A
3v₁ + 1v₂
julia> v1 + A
4v₁ + 1v₂
julia> v1 + A
5v₁ + 1v₂
Applying abs2
to a Chain
loses type information
julia> S = S"∞∅++"
⟨∞∅++⟩
julia> V = SubManifold(S)
⟨∞∅++⟩
julia> B = Λ(S)
DirectSum.Basis{⟨∞∅++⟩,16}(v, v∞, v∅, v₁, v₂, v∞∅, v∞₁, v∞₂, v∅₁, v∅₂, v₁₂, v∞∅₁, v∞∅₂, v∞₁₂, v∅₁₂, v∞∅₁₂)
julia> x = Chain(B.v1)
0v∞ + 0v∅ + 1v₁ + 0v₂
julia> x |> typeof
Chain{⟨∞∅++⟩,1,Int64,4}
julia> abs2(x) |> typeof
MultiVector{⟨∞∅++⟩,Any,16}
While x
has element type Int
, abs2(x)
has only Any
.
http://lomont.org/math/geometric-algebra/Universal%20Geometric%20Algebra%20-%20Hestenes%20-%201988.pdf is the current url for Hestenes' Universal Geometric Algebra paper
I see that you've made it possible to use REDUCE for symbolic manipulations, but I was hoping to be able to use sympy. I tried the obvious thing:
using SymPy, Grassmann
@vars x y # This defines x and y as SymPy.Sym objects
@basis ℝ^3
x*v + y*v12
But the last line causes an error:
ERROR: setindex!() with non-isbitstype eltype is not supported by StaticArrays. Consider using SizedArray.
Stacktrace:
[1] error(::String) at ./error.jl:33
[2] setindex! at /Users/me/.julia/packages/StaticArrays/mcf7t/src/MArray.jl:134 [inlined]
[3] setmulti! at /Users/me/.julia/dev/Grassmann/src/algebra.jl:31 [inlined]
[4] +(::SValue{⟨+++⟩,0,v,Sym}, ::SValue{⟨+++⟩,2,v₁₂,Sym}) at /Users/boyle/.julia/dev/Grassmann/src/algebra.jl:796
[5] top-level scope at none:0
I tried various other combinations, and it seems that whenever the sympy scalars multiply different Grassmann basis elements, this error occurs. But even if they multiply the same basis element, something weird happens:
julia> (x*v1) + (y*v1)
x + yv₁
This sentence from the README
By default, the coefficients are required to be
<:Number
.
suggests that is should work right out of the box, because
julia> SymPy.Sym <: Number
true
I also tried
Grassmann.generate_product_algebra(SymPy.Sym)
and
Grassmann.generate_product_algebra(SymPy.Sym,:(SymPy.:*),:(SymPy.:+),:(SymPy.:-))
to no avail.
I'm sorry that my julia skills are insufficient to let me figure this out, but this is my first day with this language (being attracted mostly by this exciting package). I would appreciate any help you can give.
Hi @chakravala,
In your JuliaCon paper, I noticed the statement "Mixed-symmetry algebra with Leibniz.jl
and Grassmann.jl
, having the geometric algebraic product chain rule, yields automatic differ- entiation and...", and I found a related simple example in the README.md
file, section "Differential forms and Leibniz tangent algebra"
.
Would you mind providing one or two concrete examples on symbolic differentiation using Grassmann.jl
? For example, with respect to a rotor/versor, or a bivector.
I'm thinking that the key part is about a proper setup of a @mixedbasis tangent ...
. Am I right?
I'm reading Chapter 8 of the book "Geometric Algebra for Computer Science", but there are just a few examples of basic functions. I would like to apply vector/multivector differentiation on some more functions, and verify my results.
I'm trying to differentiate a function like \sum_i weight_i || V Bivector_0 /V \wedge v_i ||^2
with respect to the versor V
.
Thank you!
I'm trying to calculate the inverse of a blade with Grassmann v0.1.3
:
julia> basis"3"
(⟨+++⟩, v, v₁, v₂, v₃, v₁₂, v₁₃, v₂₃, v₁₂₃)
julia> inv(v1+v2)
ERROR: MethodError: no method matching inv(::MBlade{Int64,⟨+++⟩,1,3})
Is this function missing, or am I looking in the wrong place?
julia> using Grassmann
julia> @basis tangent(ℝ^0, 2, 2)
(T²⟨₁₂⟩, v, ∂₁, ∂₂, ∂₁₂)
julia> [1∂1, 2∂2]
Bus error: 10
This is with
julia> versioninfo()
Julia Version 1.4.2-pre.0
Commit ef4fe83698* (2020-04-15 16:24 UTC)
Platform Info:
OS: macOS (x86_64-apple-darwin19.3.0)
CPU: Intel(R) Core(TM) i7-8850H CPU @ 2.60GHz
WORD_SIZE: 64
LIBM: libopenlibm
LLVM: libLLVM-8.0.1 (ORCJIT, skylake)
and
(GrassmannPlayground) pkg> status
Project GrassmannPlayground v0.1.0
Status `~/src/jl/GrassmannPlayground/Project.toml`
[c3fe647b] AbstractAlgebra v0.9.1
[537997a7] AbstractPlotting v0.10.1
[159f3aea] Cairo v1.0.3
[a81c6b42] Compose v0.8.2
[22fd7b30] DirectSum v0.5.6
[8d0d7f98] GaloisFields v1.0.1
[5c1252a2] GeometryBasics v0.2.4
[a2cc645c] GraphPlot v0.4.2
[4df31cd9] Grassmann v0.5.7
[edad4870] Leibniz v0.0.5
[093fc24a] LightGraphs v1.3.1
[ee78f7c6] Makie v0.10.0
[93e0c654] Reduce v1.2.5
Wanted to add an example for barycentric coordinates to your tutorial section ( as suggested on bivector.net). Ran into a problem.
julia> basis"3"
julia> a = 0v1
0v₁
julia> b = 1v1
1v₁
julia> c = 1v2
1v₂
julia> ab = b-a
1v₁
julia> ca = a-c
0v₁ - 1v₂ + 0v₃
julia> bc = c-b
-1v₁ + 1v₂ + 0v₃
julia> A = -ab∧ca
1v₁₂ + 0v₁₃ + 0v₂₃
julia> bc∧ca
1v₁₂ + 0v₁₃ + 0v₂₃
julia> bc∧ca/A
ERROR: DimensionMismatch("No precise constructor for SArray{Tuple{8},Int64,1,8} found. Length of input was 3.")
Stacktrace:
[1] SArray{Tuple{8},Int64,1,8}(::Tuple{Tuple{Tuple{Tuple{Int64,Int64,Int64}}}}) at /Users/mewert/.julia/packages/StaticArrays/1g9bq/src/convert.jl:1
[2] StaticArray at /Users/mewert/.julia/packages/StaticArrays/1g9bq/src/convert.jl:4 [inlined] (repeats 3 times)
[3] macro expansion at /Users/mewert/.julia/packages/StaticArrays/1g9bq/src/convert.jl:4 [inlined]
[4] reverse(::Chain{⟨+++⟩,2,Int64,3}) at /Users/mewert/.julia/packages/Grassmann/OeTVc/src/products.jl:531
[5] ~ at /Users/mewert/.julia/packages/DirectSum/Ofduh/src/generic.jl:242 [inlined]
[6] inv(::Chain{⟨+++⟩,2,Int64,3}) at /Users/mewert/.julia/packages/Grassmann/OeTVc/src/algebra.jl:671
[7] /(::Chain{⟨+++⟩,2,Int64,3}, ::Chain{⟨+++⟩,2,Int64,3}) at /Users/mewert/.julia/packages/AbstractTensors/cHS7X/src/AbstractTensors.jl:156
[8] top-level scope at REPL[31]:1
Seems to be the zero elements
julia> A
1v₁₂ + 0v₁₃ + 0v₂₃
julia> (bc∧ca)/1v12
1.0v⃖
This issue is to track and discuss the new updates related to differential SubManifold algebra.
chakravala/DirectSum.jl@e76bea7
chakravala/DirectSum.jl@1009e80
chakravala/DirectSum.jl@60c0e27
chakravala/AbstractTensors.jl@b7d2ff6
This is related to the finalization of the JuliaCon submission from JuliaCon/proceedings-review#18
(v1.3) pkg> add Grassmann
Updating registry at ~/.julia/registries/General
Updating git-repo https://github.com/JuliaRegistries/General.git
Resolving package versions...
Installed ReplMaker ─────── v0.2.3
Installed Requires ──────── v1.0.1
Installed AbstractTensors ─ v0.4.3
Installed DirectSum ─────── v0.5.3
Installed Reduce ────────── v1.2.4
Installed Grassmann ─────── v0.5.0
Updating ~/.julia/environments/v1.3/Project.toml
[4df31cd9] + Grassmann v0.5.0
Updating ~/.julia/environments/v1.3/Manifest.toml
[398f06c4] + AbstractLattices v0.1.2
[a8e43f4a] + AbstractTensors v0.4.3
[861a8166] + Combinatorics v1.0.0
[459fdd68] + ComputedFieldTypes v0.1.0
[22fd7b30] + DirectSum v0.5.3
[9dda63f9] + ForceImport v0.0.3
[4df31cd9] + Grassmann v0.5.0
[edad4870] + Leibniz v0.0.3
[f27b6e38] + Polynomials v0.6.0
[3cdcf5f2] + RecipesBase v0.7.0
[93e0c654] + Reduce v1.2.4
[b873ce64] + ReplMaker v0.2.3
[ae029012] + Requires v1.0.1
[90137ffa] + StaticArrays v0.12.1
[a4af3ec5] + SyntaxTree v1.0.1
[2a0f44e3] + Base64
[8ba89e20] + Distributed
[b77e0a4c] + InteractiveUtils
[8f399da3] + Libdl
[37e2e46d] + LinearAlgebra
[56ddb016] + Logging
[d6f4376e] + Markdown
[3fa0cd96] + REPL
[9a3f8284] + Random
[ea8e919c] + SHA
[9e88b42a] + Serialization
[6462fe0b] + Sockets
[2f01184e] + SparseArrays
[10745b16] + Statistics
[8dfed614] + Test
[cf7118a7] + UUIDs
[4ec0a83e] + Unicode
Building Reduce → ~/.julia/packages/Reduce/7af5f/deps/build.log
julia> using Grassmann
[ Info: Precompiling Grassmann [4df31cd9-4c27-5bea-88d0-e6a7146666d8]
ERROR: LoadError: UndefVarError: VectorBundle not defined
Stacktrace:
[1] top-level scope at /Users/mewert/.julia/packages/Leibniz/XlDhG/src/Leibniz.jl:53
[2] include at ./boot.jl:328 [inlined]
[3] include_relative(::Module, ::String) at ./loading.jl:1105
[4] include(::Module, ::String) at ./Base.jl:31
[5] top-level scope at none:2
[6] eval at ./boot.jl:330 [inlined]
[7] eval(::Expr) at ./client.jl:425
[8] top-level scope at ./none:3
in expression starting at /Users/mewert/.julia/packages/Leibniz/XlDhG/src/Leibniz.jl:53
ERROR: LoadError: Failed to precompile Leibniz [edad4870-8a01-11e9-2d75-8f02e448fc59] to /Users/mewert/.julia/compiled/v1.3/Leibniz/wNkI4_RK5IG.ji.
Stacktrace:
[1] error(::String) at ./error.jl:33
[2] compilecache(::Base.PkgId, ::String) at ./loading.jl:1283
[3] _require(::Base.PkgId) at ./loading.jl:1024
[4] require(::Base.PkgId) at ./loading.jl:922
[5] require(::Module, ::Symbol) at ./loading.jl:917
[6] include at ./boot.jl:328 [inlined]
[7] include_relative(::Module, ::String) at ./loading.jl:1105
[8] include(::Module, ::String) at ./Base.jl:31
[9] top-level scope at none:2
[10] eval at ./boot.jl:330 [inlined]
[11] eval(::Expr) at ./client.jl:425
[12] top-level scope at ./none:3
in expression starting at /Users/mewert/.julia/packages/Grassmann/5163y/src/Grassmann.jl:46
ERROR: Failed to precompile Grassmann [4df31cd9-4c27-5bea-88d0-e6a7146666d8] to /Users/mewert/.julia/compiled/v1.3/Grassmann/i5fGo_RK5IG.ji.
Stacktrace:
[1] error(::String) at ./error.jl:33
[2] compilecache(::Base.PkgId, ::String) at ./loading.jl:1283
[3] _require(::Base.PkgId) at ./loading.jl:1024
[4] require(::Base.PkgId) at ./loading.jl:922
[5] require(::Module, ::Symbol) at ./loading.jl:917
julia>
julia> @basis 1
(⟨+⟩, v, v₁)
julia> @basis 2
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> A = Chain{V,1}(Chain(v1), Chain(v1))
(1v₁ + 0v₂)v₁ + (1v₁ + 0v₂)v₂
julia> x = Chain(v2)
0v₁ + 1v₂
julia> A ⋅ x
1v₁ + 0v₂
julia> A * x
ERROR: MethodError: no method matching Chain{⟨++⟩,1,Int64,2}(::Float64)
Closest candidates are:
Chain{⟨++⟩,1,Int64,2}(::T) where T<:Number at boot.jl:715
Chain{⟨++⟩,1,Int64,2}(::Simplex{V,0,B,T} where T where B) where {V, G, T, X} at /Users/eschnett/.julia/packages/Grassmann/3DicJ/src/multivectors.jl:44
Chain{⟨++⟩,1,Int64,2}(::Base.TwicePrecision) where T<:Number at twiceprecision.jl:243
...
Stacktrace:
[1] convert(::Type{Chain{⟨++⟩,1,Int64,2}}, ::Float64) at ./number.jl:7
[2] macro expansion at /Users/eschnett/.julia/packages/StaticArrays/1g9bq/src/util.jl:11 [inlined]
[3] convert_ntuple at /Users/eschnett/.julia/packages/StaticArrays/1g9bq/src/util.jl:8 [inlined]
[4] SArray{Tuple{4},Chain{⟨++⟩,1,Int64,2},1,4}(::Tuple{Chain{⟨++⟩,1,Int64,2},Float64,Float64,Chain{⟨++⟩,1,Int64,2}}) at /Users/eschnett/.julia/packages/StaticArrays/1g9bq/src/SArray.jl:28
[5] StaticArray at /Users/eschnett/.julia/packages/StaticArrays/1g9bq/src/convert.jl:4 [inlined]
[6] macro expansion at /Users/eschnett/.julia/packages/Grassmann/3DicJ/src/products.jl:0 [inlined]
[7] *(::Chain{⟨++⟩,1,Chain{⟨++⟩,1,Int64,2},2}, ::Chain{⟨++⟩,1,Int64,2}) at /Users/eschnett/.julia/packages/Grassmann/3DicJ/src/products.jl:856
[8] top-level scope at REPL[67]:1
[9] eval(::Module, ::Any) at ./boot.jl:331
[10] eval_user_input(::Any, ::REPL.REPLBackend) at /Users/eschnett/src/julia-1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
[11] run_backend(::REPL.REPLBackend) at /Users/eschnett/.julia/packages/Revise/0meWR/src/Revise.jl:1075
[12] top-level scope at none:0
Which definition of the inner product between multivectors do you use? I'm asking because I was confused by the fact that v1 ⋅ v12
is giving me 0
as a result, and I would have expected some multiple of v1 (which v12 ⋅ v1
gives me).
Am I wrong or is it a bug ?
using Reduce
using Grassmann
@basis S"∞∅+++"
p = :x*v1 + :y*v2 + :z*v3
P = v∅ + p + 0.5 * p^2 * v∞
ERROR: LoadError: MethodError: no method matching +(::Expr, ::Float64)
Closest candidates are:
+(::Any, ::Any, !Matched::Any, !Matched::Any...) at operators.jl:529
+(!Matched::Bool, ::T<:AbstractFloat) where T<:AbstractFloat at bool.jl:104
+(!Matched::Float64, ::Float64) at float.jl:395
...
Stacktrace:
[1] macro expansion at .julia/packages/StaticArrays/3KEjZ/src/mapreduce.jl:41 [inlined]
[2] _map at .julia/packages/StaticArrays/3KEjZ/src/mapreduce.jl:21 [inlined]
[3] map at .julia/packages/StaticArrays/3KEjZ/src/mapreduce.jl:14 [inlined]
[4] +(::Array{Any,1}, ::StaticArrays.SizedArray{Tuple{5},Any,1,1}) at /home/gacquewi/.julia/packages/StaticArrays/3KEjZ/src/linalg.jl:11
[5] +(::MChain{Any,⟨∞∅+++⟩,1,5}, ::MultiVector{Any,⟨∞∅+++⟩,32}) at .julia/packages/Grassmann/cEoGT/src/algebra.jl:1310
[6] +(::Basis{⟨∞∅+++⟩,1,0x0000000000000002}, ::MChain{Any,⟨∞∅+++⟩,1,5}, ::MultiVector{Any,⟨∞∅+++⟩,32}) at ./operators.jl:529
[7] top-level scope at wgCGAbug.jl:5
[8] include at ./boot.jl:328 [inlined]
[9] include_relative(::Module, ::String) at ./loading.jl:1094
[10] include(::Module, ::String) at ./Base.jl:31
[11] exec_options(::Base.JLOptions) at ./client.jl:295
[12] _start() at ./client.jl:464
in expression starting at wgCGAbug.jl:5
The following is a quote of @enkimute from the discussion in #3
check the way he defined the complement and figure out if this may provide a method that works for degenerate metrics in the same unified fashion..
Opening this issue is meant to track ongoing discussions regarding degenerate metrics and complements
The inverse of a vector a
works like expected for vectors with positive square a² > 0
in the way that baa⁻¹ = b
:
julia> using Grassmann
julia> basis"-+++"
(⟨-+++⟩, v, v₁, v₂, v₃, v₄, v₁₂, v₁₃, v₁₄, v₂₃, v₂₄, v₃₄, v₁₂₃, v₁₂₄, v₁₃₄, v₂₃₄, v₁₂₃₄)
julia> v2^2
v
julia> v1*v2*inv(v2), v1
(v₁, v₁)
It however does not give the correct result, I think, for vectors with negative square a² < 0
:
julia> v1^2
-1v
julia> v2*v1*inv(v1), v2
(-1v₂, v₂)
This is fundamentally tied to the Hodge star dual operation
https://en.wikipedia.org/wiki/Hodge_star_operator
Help or discussion is welcome on this topic
For property testing I like to use random vectors (simplices, chains, multivectors). It would be convenient to have a method for rand()
that does that, called e.g. as rand(Chain{V,1,Float64})
.
This is not urgent.
I find this:
julia> using Grassmann
julia> V=SubManifold(Signature(1))
julia> Chain{V,1}(SVector(BigInt(10))) |> typeof
Chain{⟨+⟩,1,BigInt,1}
julia> +Chain{V,1}(SVector(BigInt(10))) |> typeof
Chain{⟨+⟩,1,Any,1}
Applying a unary +
or -
to a chain converts the type to Any
.
May be helpful to add AbstractPlotting.inline!(true)
to Makie plot examples for newcomers so they're not confused (like I was!). MakieOrg/Makie.jl#392
julia> basis"+++"
(⟨+++⟩, v, v₁, v₂, v₃, v₁₂, v₁₃, v₂₃, v₁₂₃)
julia> v1⨼(v1*v2)
-1v₂
julia> v1*(v1*v2)
v₂
For some reason, I'm getting the following behavior for exponentials:
exp(0 * (v1 + v2)) # 1v⃖ (as expected)
exp(0 * v1) # 1.0 - NaNv₁ (unexpected)
The problem appears to be specific to zero multiples of basis elements.
In geometric algebra at least the 0-vector wedged with another 0-vector should be non-zero if both scalars are non-zero.
julia> 2v∧3v
0v
The result type of abs2
seems to depend on argument values:
julia> @basis 2
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> abs2(MultiVector(v+v1)) |> typeof
MultiVector{⟨++⟩,Int64,4}
julia> abs2(MultiVector(v)) |> typeof
Simplex{⟨++⟩,0,v,Int64}
The return type of a function should usually only depend on the argument types.
I am currently not running benchmarks, so this does currently not affect me.
I am experimenting with mixedbasis
to represent non-trivial metric tensors. (Is there a better way to represent them?) I receive these errors:
julia> using Grassmann
julia> S = S"-+++"
julia> @mixedbasis S
julia> v1w1(v1)
ERROR: MethodError: no method matching evaluate2(::SubManifold{v₁₂₃₄w¹²³⁴,2,0x0000000000000011}, ::SubManifold{⟨-++++---⟩*,1,0x0000000000000001})
Closest candidates are:
evaluate2(::A, ::B) where {V, A<:TensorTerm{V,1}, B<:TensorTerm{V,1}} at /Users/eschnett/.julia/packages/DirectSum/Ng7R9/src/operations.jl:229
Stacktrace:
[1] (::SubManifold{v₁₂₃₄w¹²³⁴,2,0x0000000000000011})(::SubManifold{⟨-++++---⟩*,1,0x0000000000000001}) at /Users/eschnett/.julia/packages/DirectSum/Ng7R9/src/operations.jl:249
[2] interform at /Users/eschnett/.julia/packages/AbstractTensors/zchWS/src/AbstractTensors.jl:124 [inlined]
[3] (::SubManifold{⟨-++++---⟩*,2,0x0000000000000011})(::SubManifold{⟨-+++⟩,1,0x0000000000000001}) at /Users/eschnett/.julia/packages/DirectSum/Ng7R9/src/operations.jl:255
[4] top-level scope at REPL[10]:1
[5] eval(::Module, ::Any) at ./boot.jl:331
[6] eval_user_input(::Any, ::REPL.REPLBackend) at /Users/eschnett/src/julia-1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
[7] run_backend(::REPL.REPLBackend) at /Users/eschnett/.julia/packages/Revise/C272c/src/Revise.jl:1075
[8] top-level scope at none:0
It seems that operations are not defined on pure basis vectors, only on linear combinations of them. However:
julia> (-v1w1+v2w2)(v1)
ERROR: InexactError: Bool(-1)
Stacktrace:
[1] Bool at ./float.jl:73 [inlined]
[2] convert at ./number.jl:7 [inlined]
[3] convert at ./essentials.jl:310 [inlined] (repeats 2 times)
[4] setindex! at ./array.jl:826 [inlined]
[5] dualform(::SubManifold{⟨-++++---⟩*,8,0x00000000000000ff}) at /Users/eschnett/.julia/dev/Grassmann/src/forms.jl:142
[6] (::Chain{⟨-++++---⟩*,2,Int64,28})(::SubManifold{⟨-++++---⟩*,1,0x0000000000000001}) at /Users/eschnett/.julia/dev/Grassmann/src/forms.jl:187
[7] interform at /Users/eschnett/.julia/packages/AbstractTensors/zchWS/src/AbstractTensors.jl:124 [inlined]
[8] (::Chain{⟨-++++---⟩*,2,Int64,28})(::SubManifold{⟨-+++⟩,1,0x0000000000000001}) at /Users/eschnett/.julia/dev/Grassmann/src/forms.jl:176
[9] top-level scope at REPL[11]:1
[10] eval(::Module, ::Any) at ./boot.jl:331
[11] eval_user_input(::Any, ::REPL.REPLBackend) at /Users/eschnett/src/julia-1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
[12] run_backend(::REPL.REPLBackend) at /Users/eschnett/.julia/packages/Revise/C272c/src/Revise.jl:1075
[13] top-level scope at none:0
The latter seems to happen in line 142 of forms.jl
:
@inbounds mV[Q] = (intlog(Y)+1,V[intlog(x)+1])
Should the second tuple element have an != 0
added to handle negative signatures?
Hi,
On notation:
What is the difference between every pair of the three: @basis S"++++-"
, @basis V"++++-"
and @basis "++++-"
?
I'm a little confused in trying the examples in README.md
. It seems they all give me the same results.
On algebra
In Grassmann.js
, is V"++++-"
also the conformal model G^{4,1}
of R^3
?
(I noticed in the section "Null-basis of the conformal split" that we can declare a space with signature ∞ and ∅, but did not figure out how I can construct a signature such that it is in correspondence with "++++-".
Thank you,
Hongying
I get an error for the geometric product of v∅
and -v∞
:
julia> using Grassmann
julia> @basis S"∞∅+"
(⟨∞∅+⟩, v, v∞, v∅, v₁, v∞∅, v∞₁, v∅₁, v∞∅₁)
julia> v∅*v∞
-1 - 1v∞∅
julia> v∅*(-v∞)
ERROR: MethodError: no method matching SValue{⟨∞∅+⟩,G,B,T} where T where B where G(::Int64, ::MultiVector{Int64,⟨∞∅+⟩,8})
Closest candidates are:
SValue{⟨∞∅+⟩,G,B,T} where T where B where G(::T) where {V, T} at .julia/packages/Grassmann/H9Zog/src/multivectors.jl:97
SValue{⟨∞∅+⟩,G,B,T} where T where B where G(::Any, ::SValue{V,G,B,T} where T where B) where {V, G} at .julia/packages/Grassmann/H9Zog/src/multivectors.jl:90
SValue{⟨∞∅+⟩,G,B,T} where T where B where G(::Any, ::MValue{V,G,B,T} where T where B) where {V, G} at .julia/packages/Grassmann/H9Zog/src/multivectors.jl:91
...
Stacktrace:
[1] *(::Basis{⟨∞∅+⟩,1,0x0000000000000002}, ::SValue{⟨∞∅+⟩,1,v∞,Int64}) at .julia/packages/Grassmann/H9Zog/src/algebra.jl:157
[2] top-level scope at none:0
I am testing Grassmann.jl with conformal projective geometry, but the geometric product seems to be broken for those spaces.
Easiest way to reproduce is with your example for the null-basis conformal split in README.md. That gives me:
julia> using Grassmann; @basis S"∞∅++"
[ Info: Recompiling stale cache file /home/karl/.julia/compiled/v1.1/Grassmann/i5fGo.ji for Grassmann [4df31cd9-4c27-5bea-88d0-e6a7146666d8]
(⟨∞∅++⟩, v, v∞, v∅, v₁, v₂, v∞∅, v∞₁, v∞₂, v∅₁, v∅₂, v₁₂, v∞∅₁, v∞∅₂, v∞₁₂, v∅₁₂, v∞∅₁₂)
julia> v∞^2, v∅^2, v1^2, v2^2
ERROR: MethodError: no method matching conformal(::Signature{4,3,0x0000000000000002,0}, ::UInt64, ::UInt64)
Closest candidates are:
conformal(::Any, ::Any, ::W<:(Signature{n,m,s,Diff} where Diff)) where {n, m, s, W<:(Signature{n,m,s,Diff} where Diff)} at .julia/packages/Grassmann/wWv7E/src/parity.jl:164
conformal(::Any, ::Any, ::W<:(DiagonalForm{n,m,s,Diff} where Diff)) where {n, m, s, W<:(DiagonalForm{n,m,s,Diff} where Diff)} at .julia/packages/Grassmann/wWv7E/src/parity.jl:164
Stacktrace:
[1] *(::Basis{⟨∞∅++⟩,1,0x0000000000000004}, ::Basis{⟨∞∅++⟩,1,0x0000000000000004}) at .julia/packages/Grassmann/wWv7E/src/algebra.jl:148
[2] ^(::Basis{⟨∞∅++⟩,1,0x0000000000000004}, ::Int64) at .julia/packages/Grassmann/wWv7E/src/algebra.jl:919
[3] macro expansion at ./none:0 [inlined]
[4] literal_pow(::typeof(^), ::Basis{⟨∞∅++⟩,1,0x0000000000000004}, ::Val{2}) at ./none:0
[5] top-level scope at none:0
It worked in 0.1.4, but doesn't in 0.1.5
Hi.
Thanks for helping to push geometric algebra into the mainstream with your hard work.
I am reading introductory material on GA and using Clifford (Python) to add, multiply, reflect, etc.
For the most part, I bypassed linear algebra and just took up with GA because it looks interesting.
I would like to use Julia (and Grassmann) to learn GA, but I don't understand most of the README, yet.
In particular, I couldn't tell how to simply add, multiply or reflect a vector in R3.
I was able to use Clifford because I could follow their docs, but the Grassmann README seems geared to people who already know what they are doing...
Thanks, again.
I encounter this error:
julia> basis"+++"
julia> (v1+v2) + (v1+v2)*(v1+v2)
ERROR: DimensionMismatch("No precise constructor for StaticArrays.MArray{Tuple{3},Int64,1,3} found. Length of input was 8.")
Stacktrace:
[1] StaticArrays.MArray{Tuple{3},Int64,1,3}(::Tuple{Tuple{Tuple{NTuple{8,Int64}}}}) at /Users/eschnett/.julia/packages/StaticArrays/VyRz3/src/convert.jl:1
[2] Type at /Users/eschnett/.julia/packages/StaticArrays/VyRz3/src/convert.jl:4 [inlined] (repeats 3 times)
[3] convert at /Users/eschnett/.julia/packages/StaticArrays/VyRz3/src/convert.jl:10 [inlined]
[4] value at /Users/eschnett/.julia/packages/Grassmann/LeGo1/src/multivectors.jl:327 [inlined]
[5] +(::MBlade{Int64,⟨+++⟩,1,3}, ::MultiVector{Int64,⟨+++⟩,8}) at /Users/eschnett/.julia/packages/Grassmann/LeGo1/src/algebra.jl:1045
[6] top-level scope at none:0
The documentation for \oslash reads:
Sandwich product: ω⊘η = (~ω)⊖η⊖ω
However, it seems to implement η⊘ω
instead (with arguments commuted):
julia> R = exp(π/4*v12)
0.7071067811865476 + 0.7071067811865475v₁₂
julia> (~R) ⊖ v1 ⊖ R
0.0 + 2.220446049250313e-16v₁ + 1.0v₂
julia> R ⊘ v1
-0.7071067811865476 + 0.7071067811865475v₁₂
julia> v1 ⊘ R
0.0 + 2.220446049250313e-16v₁ + 1.0v₂
This fails:
julia> @basis 0
(⟨⟩, ⟨⟩)
julia> Chain{V,1}(v)
ERROR: DimensionMismatch("No precise constructor for SArray{Tuple{0},SubManifold{⟨⟩,0,0x0000000000000000},1,0} found. Length of input was 1.")
Stacktrace:
[1] SArray{Tuple{0},SubManifold{⟨⟩,0,0x0000000000000000},1,0}(::Tuple{Tuple{Tuple{Tuple{SubManifold{⟨⟩,0,0x0000000000000000}}}}}) at /Users/eschnett/.julia/packages/StaticArrays/MB5hl/src/convert.jl:1
[2] StaticArray at /Users/eschnett/.julia/packages/StaticArrays/MB5hl/src/convert.jl:4 [inlined] (repeats 4 times)
[3] Chain{⟨⟩,1,𝕂,253} where 253 where 𝕂(::SubManifold{⟨⟩,0,0x0000000000000000}) at /Users/eschnett/.julia/packages/Grassmann/Of0BW/src/multivectors.jl:33
[4] top-level scope at REPL[209]:1
[5] eval(::Module, ::Any) at ./boot.jl:331
[6] eval_user_input(::Any, ::REPL.REPLBackend) at /Users/eschnett/src/julia-1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
[7] run_backend(::REPL.REPLBackend) at /Users/eschnett/.julia/packages/Revise/0meWR/src/Revise.jl:1075
[8] top-level scope at none:0
This is of low priority to me.
I encountered this problem:
julia> using Grassmann
julia> D = 2
2
julia> S = Signature(D)
⟨++⟩
julia> V = SubManifold(S)
⟨++⟩
julia> x = Chain{V,1}(1,2)
1v₁ + 2v₂
julia> xs = [x]
1-element Array{Chain{⟨++⟩,1,Int64,2},1}:
1v₁ + 2v₂
julia> @show xs
ERROR: _show_nonempty(::IO, ::AbstractVector, ::String) is not implemented
Stacktrace:
[1] error(::String) at ./error.jl:33
[2] _show_nonempty(::IOContext{Base.GenericIOBuffer{Array{UInt8,1}}}, ::Array{Chain{⟨++⟩,1,Int64,2},1}, ::String) at ./arrayshow.jl:420
[3] show(::Base.GenericIOBuffer{Array{UInt8,1}}, ::Array{Chain{⟨++⟩,1,Int64,2},1}) at ./arrayshow.jl:437
[4] sprint(::Function, ::Array{Chain{⟨++⟩,1,Int64,2},1}; context::Nothing, sizehint::Int64) at ./strings/io.jl:105
[5] #repr#353 at ./strings/io.jl:227 [inlined]
[6] repr(::Array{Chain{⟨++⟩,1,Int64,2},1}) at ./strings/io.jl:227
[7] top-level scope at show.jl:634
The problem seems to be that Grassmann defines (see multivectors.jl
, line 139) the method
@pure Base.ndims(::Vector{Chain{V,G,T,X}} where {G,T,X}) where V = ndims(V)
which confuses the base library when showing vectors.
I think it's wrong to overload ndims
for any Vector{T}
. All Vector{T}
are one-dimensional arrays. The corresponding definition for Base.parent
probably also shouldn't be there.
This breaks my code since using @show
is important for debugging for me.
The source code for this library has a very large number of functions marked @pure
. My understanding is that the use of this macro is emphatically discouraged. In this issue JuliaLang/julia#20704, Keno said
@pure
meansI am @vtjnash and I certify that this function doesn't break any of the implicit assumptions I use in the compiler
;)
In this thread: https://discourse.julialang.org/t/pure-macro/3871/3, Matt says that
Improper
@pure
annotations can introduce bugs. The optimizations it enables rely on an extremely strict definition of pure. It really should be named something like@hyperpure
. Some of the restrictions include:
- It must always return exactly (
===
) the same result for a given input. Watch out for mutable types. I think constant globals are okay, though.- The function it’s used on cannot be further extended by other methods after it gets called.
- It cannot recurse.
- It’s undocumented and not exported (for good reason), but this means the complete list of preconditions is really only in a few people’s heads.
For instance,
Line 42 in f4e4f78
@pure
since an array is returned and pure functions must return outputs that are equal in the ===
sense.
Furthermore, adding methods to @pure
functions is said to not be kosher and this repository is replete with @pure
functions with methods added on.
In this thread https://discourse.julialang.org/t/can-pure-functions-throw-an-error/18459, Jameson advises that @pure
functions throwing errors is "Probably, a very very bad idea", but I see that there are functions that throw
in multivectors.jl
and Grassmann.jl
.
Given this, I think cutting back on most if not all uses of the uses of @pure
in this package is probably a good idea.
Perhaps this deserves a separate issue, but I noticed that in utilities.jl
there are struct definitions
struct Grade{G}
@pure Grade{G}() where G = new{G}()
end
## Dimension{N}
struct Dimension{N}
@pure Dimension{N}() where N = new{N}()
end
I believe the above inner constructors are no-ops that don't need to exist. For instance:
julia> struct MyStruct1{T} end
julia> struct MyStruct2{T}
MyStruct2{T}() where {T} = new{T}()
end
julia> @code_lowered MyStruct1{1}()
CodeInfo(
1 ─ %1 = (Core.apply_type)(Main.MyStruct1, $(Expr(:static_parameter, 1)))
│ %2 = %new(%1)
└── return %2
)
julia> @code_lowered MyStruct2{1}()
CodeInfo(
1 ─ %1 = (Core.apply_type)(Main.MyStruct2, $(Expr(:static_parameter, 1)))
│ %2 = %new(%1)
└── return %2
)
The equality check between Multivectors and Numbers currently does not work correctly for numbers that are not zero:
julia> using Grassmann
julia> basis"2"
(⟨++⟩, v, v₁, v₂, v₁₂)
julia> a = v + v1 - v1
1v⃖
julia> typeof(a)
MultiVector{Int64,⟨++⟩,4}
julia> a == 1
false
julia> b = a - 1
0v⃖
julia> b == 0
true
julia> a - 1 == 0
true
A simple fix would probably be changing lines 309 and 310 in multivectors.jl to:
==(a::Number,b::MultiVector{S,V,G} where {S,V}) where G = prod(0 .== value(b - a))
==(a::MultiVector{S,V,G} where {S,V},b::Number) where G = prod(0 .== value(a - b))
The readme contains this example
using Grassmann, Makie; @basis S"∞+++"
streamplot(vectorfield(exp((π/4)*(v12+v∞3)),V(2,3,4),V(1,2,3)),-1.5..1.5,-1.5..1.5,-1.5..1.5,gridsize=(10,10))
With just the statements above, vectorfield
is undefined. If I also ] add GeometryTypes
and using GeometryTypes
, then vectorfield
is found, but I'm then stuck at the next error:
julia> streamplot(vectorfield(exp((π/4)*(v12+v∞3)),V(2,3,4),V(1,2,3)),-1.5..1.5,-1.5..1.5,-1.5..1.5,gridsize=(10,10))
ERROR: MethodError: no method matching ptype(::GeometryBasics.Point{3,Float64})
Closest candidates are:
ptype(::GeometryTypes.Point{N,T} where N) where T at /Users/eschnett/.julia/packages/Grassmann/9buUZ/src/Grassmann.jl:307
Stacktrace:
[1] (::Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}})(::GeometryBasics.Point{3,Float64}) at /Users/eschnett/.julia/packages/Grassmann/9buUZ/src/Grassmann.jl:310
[2] (::AbstractPlotting.var"#apply_f#505"{Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}}})(::GeometryBasics.Point{3,Float64}, ::Type{T} where T) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/basic_recipes/basic_recipes.jl:1061
[3] streamplot_impl(::Type{T} where T, ::Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}}, ::GeometryBasics.HyperRectangle{3,Float32}, ::Tuple{Int64,Int64}, ::Float64, ::Int64, ::Float64) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/basic_recipes/basic_recipes.jl:1081
[4] (::AbstractPlotting.var"#508#512")(::Function, ::GeometryBasics.HyperRectangle{3,Float32}, ::Tuple{Int64,Int64}, ::Float64, ::Int64, ::Float64) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/basic_recipes/basic_recipes.jl:1138
[5] lift(::Function, ::Observables.Observable{Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}}}, ::Observables.Observable{GeometryBasics.HyperRectangle{3,Float32}}, ::Vararg{Any,N} where N) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/interaction/nodes.jl:8
[6] plot!(::StreamPlot{...}) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/basic_recipes/basic_recipes.jl:1132
[7] plot!(::Scene, ::Type{StreamPlot{...}}, ::Attributes, ::Tuple{Observables.Observable{Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}}},Observables.Observable{IntervalSets.Interval{:closed,:closed,Float64}},Observables.Observable{IntervalSets.Interval{:closed,:closed,Float64}},Observables.Observable{IntervalSets.Interval{:closed,:closed,Float64}}}, ::Observables.Observable{Tuple{Grassmann.var"#707#709"{MultiVector{⟨∞+++⟩,Float64,16},SubManifold{⟨∞+++⟩,3,0x000000000000000e},SubManifold{⟨∞+++⟩,3,0x0000000000000007}},GeometryBasics.HyperRectangle{3,Float32}}}) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/interfaces.jl:633
[8] plot!(::Scene, ::Type{StreamPlot{...}}, ::Attributes, ::Function, ::Vararg{Any,N} where N; kw_attributes::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/interfaces.jl:563
[9] plot!(::Scene, ::Type{StreamPlot{...}}, ::Attributes, ::Function, ::IntervalSets.Interval{:closed,:closed,Float64}, ::IntervalSets.Interval{:closed,:closed,Float64}, ::Vararg{IntervalSets.Interval{:closed,:closed,Float64},N} where N) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/interfaces.jl:532
[10] streamplot(::Function, ::Vararg{Any,N} where N; attributes::Base.Iterators.Pairs{Symbol,Tuple{Int64,Int64},Tuple{Symbol},NamedTuple{(:gridsize,),Tuple{Tuple{Int64,Int64}}}}) at /Users/eschnett/.julia/packages/AbstractPlotting/7mERO/src/recipes.jl:15
[11] top-level scope at REPL[8]:1
[12] eval(::Module, ::Any) at ./boot.jl:331
[13] eval_user_input(::Any, ::REPL.REPLBackend) at /Users/eschnett/src/julia-1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
[14] run_backend(::REPL.REPLBackend) at /Users/eschnett/.julia/packages/Revise/C272c/src/Revise.jl:1075
[15] top-level scope at none:0
I am using projective or conformal bases, such as
julia> S = S"∞++"
julia> B = Λ(S)
I want to write dimension-generic code, and want to access the i
-th basis vector. Unfortunately, Grassmann counts starting from v∞
, so that v₁
is B.v(2)
etc. It would be more convenient if B.v(1)
returned v1
, independent of whether v∞
or v∅
were present. I currently have to write B.v(i+1)
or B.v(i+2)
, depending on the basis.
I realize this would be a breaking change, and obviously this is not an urgent feature request.
I wanted to do a test your package by implementing Quaternions and tried to create a quaternion q for a rotation around a vector v = (a, b, c) by an angle alpha using exponentials:
q = exp(alpha/2*(ai + bj + ck))
However the exponential function doesn't seem to be fully implemented yet for multivectors? At least i get an error in the following MWE:
julia> using Grassmann
julia> i, j, k = hyperplanes(ℝ^3)
3-element Array{SValue{⟨+++⟩,2,B,Int64} where B,1}:
-1v₂₃
1v₁₃
-1v₁₂
julia> alpha = 0.5π
1.5707963267948966
julia> exp(alpha/2*(i))
0.7071067811865476 - 0.7071067811865475v₂₃
julia> a, b, c = 1/sqrt(2) * [1, 1, 0]
3-element Array{Float64,1}:
0.7071067811865475
0.7071067811865475
0.0
julia> exp(alpha/2*(a*i + b*j + c*k))
ERROR: MethodError: no method matching abs(::StaticArrays.MArray{Tuple{8},Float64,1,8})
Closest candidates are:
abs(::Bool) at bool.jl:91
abs(::Float16) at float.jl:520
abs(::Float32) at float.jl:521
...
Stacktrace:
[1] abs(::MultiVector{Float64,⟨+++⟩,8}) at .julia/packages/AbstractTensors/fUnNY/src/AbstractTensors.jl:55
[2] _broadcast_getindex_evalf at ./broadcast.jl:578 [inlined]
[3] _broadcast_getindex at ./broadcast.jl:551 [inlined]
[4] getindex at ./broadcast.jl:511 [inlined]
[5] copyto_nonleaf!(::Array{Float64,1}, ::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(abs),Tuple{Base.Broadcast.Extruded{Array{TensorAlgebra{⟨+++⟩},1},Tuple{Bool},Tuple{Int64}}}}, ::Base.OneTo{Int64}, ::Int64, ::Int64) at ./broadcast.jl:928
[6] copy at ./broadcast.jl:791 [inlined]
[7] materialize at ./broadcast.jl:753 [inlined]
[8] exp(::SBlade{Float64,⟨+++⟩,2,3}) at .julia/packages/Grassmann/wWv7E/src/algebra.jl:961
[9] top-level scope at none:0
This fails:
julia> using Grassmann
julia> @basis tangent(ℝ^0, 1, 1)
(T¹⟨₁⟩, v, ∂₁)
julia> ∂1, ∂1 ⊗ ∂1
(∂₁, 0v)
julia> ∂1 ⊗ ∂1 ⊗ ∂1
ERROR: StackOverflowError:
Stacktrace:
[1] interop(::typeof(⊗), ::Simplex{T¹⟨₁⟩,0,v,Int64}, ::SubManifold{T¹⟨₁⟩,1,0x0000000000000001}) at /Users/eschnett/.julia/packages/AbstractTensors/zchWS/src/AbstractTensors.jl:115
[2] ⊗(::Simplex{T¹⟨₁⟩,0,v,Int64}, ::SubManifold{T¹⟨₁⟩,1,0x0000000000000001}) at /Users/eschnett/.julia/packages/AbstractTensors/zchWS/src/AbstractTensors.jl:148
... (the last 2 lines are repeated 39990 more times)
[79983] interop(::typeof(⊗), ::Simplex{T¹⟨₁⟩,0,v,Int64}, ::SubManifold{T¹⟨₁⟩,1,0x0000000000000001}) at /Users/eschnett/.julia/packages/AbstractTensors/zchWS/src/AbstractTensors.jl:115
It seems that, if ∂1^n
is zero, then ∂1^(n+1)
fails with a stack overflow.
I notice that certain constructors for Chain
work when called generically, but not when called with the specific type:
julia> Chain{V,1}(1,2)
1v₁ + 2v₂
julia> Chain{V,1,Int}(1,2)
ERROR: MethodError: no method matching Chain{⟨++⟩,1,Int64,253} where 253(::Int64, ::Int64)
Similarly:
julia> Chain(B.v1)
1v₁ + 0v₂
julia> Chain{Int}(B.v1)
1v₁ + 0v₂
but:
julia> Chain{V}(B.v1)
ERROR: MethodError: no method matching SubManifold{⟨++⟩,0,0x0000000000000000}(::Simplex{⟨++⟩,0,v,Int64})
julia> Chain{V,1}(B.v1)
ERROR: DimensionMismatch("No precise constructor for StaticArrays.SArray{Tuple{2},SubManifold{⟨++⟩,1,0x0000000000000001},1,2} found. Length of input was 1.")
julia> Chain{V,1,Int}(B.v1)
ERROR: MethodError: Cannot `convert` an object of type SubManifold{⟨++⟩,1,0x0000000000000001} to an object of type StaticArrays.SArray{Tuple{2},Int64,1,2}
It would be convenient if one could always call a specific constructor. That is, if a constructor creates a certain type Chain{V,1,T}
, then I'd prefer if I was always able to call Chain{V,1,T}
, and not have to use Chain{V,1}
or Chain{T}
instead.
This is not urgent.
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.