Skip to content

Commit

Permalink
generalized support for symbolic scalars Reduce.RExpr, SymPy.Sym #6
Browse files Browse the repository at this point in the history
  • Loading branch information
@chakravala
chakravala committed Mar 4, 2019
1 parent 00330ca commit 1edc6b5
Show file tree
Hide file tree
Showing 4 changed files with 82 additions and 104 deletions.
6 changes: 4 additions & 2 deletions src/Grassmann.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,6 @@ include("utilities.jl")
include("multivectors.jl")
include("algebra.jl")
include("forms.jl")
include("symbolic.jl")
include("generators.jl")

## fundamentals
Expand Down Expand Up @@ -227,6 +226,9 @@ end

# ParaAlgebra

__init__() = @require Reduce="93e0c654-6965-5f22-aba9-9c1ae6b3c259" import Reduce
function __init__()
@require Reduce="93e0c654-6965-5f22-aba9-9c1ae6b3c259" include("symbolic.jl")
@require SymPy="24249f21-da20-56a4-8eb1-6a02cf4ae2e6" generate_product_algebra(:(SymPy.Sym),:(SymPy.:*),:(SymPy.:+),:(SymPy.:-),:svec,:(SymPy.conj))
end

end # module
92 changes: 54 additions & 38 deletions src/algebra.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,55 +19,62 @@ function add!(out::MultiVector{T,V},val::T,A::Vector{Int},B::Vector{Int}) where
return out
end

for (op,set) ∈ ((:add,:(+=)),(:set,:(=)))
sm = Symbol(op,:multi!)
sb = Symbol(op,:blade!)
for (s,index) ∈ ((sm,:basisindex),(sb,:bladeindex))
for (i,B) ∈ ((:i,Bits),(:(bits(i)),Basis))
@eval begin
@inline function $s(out::MArray{Tuple{M},T,1,M},val::T,i::$B) where {M,T<:Field}
@inbounds $(Expr(set,:(out[$index(intlog(M),$i)]),:val))
return out
end
@inline function $s(out::Q,val::T,i::$B,::Dimension{N}) where Q<:MArray{Tuple{M},T,1,M} where {M,T<:Field,N}
@inbounds $(Expr(set,:(out[$index(N,$i)]),:val))
return out
end
end
end
end
for s ∈ (sm,sb)
@eval begin
@inline function $(Symbol(:join,s))(V::VectorSpace{N,D},m::MArray{Tuple{M},T,1,M},A::Bits,B::Bits,v::T) where {N,D,T<:Field,M}
if v ≠ 0 && !(hasdual(V) && isodd(A) && isodd(B))
$s(m,parity(A,B,V) ? -(v) : v,A ⊻ B,Dimension{N}())
end
return m
end
@inline function $(Symbol(:join,s))(m::MArray{Tuple{M},T,1,M},v::T,A::Basis{V},B::Basis{V}) where {V,T<:Field,M}
if v ≠ 0 && !(hasdual(V) && hasdual(A) && hasdual(B))
$s(m,parity(A,B) ? -(v) : v,bits(A) ⊻ bits(B),Dimension{ndims(V)}())
const Sym = :(Reduce.Algebra)
const SymField = Any

set_val(set,expr,val) = Expr(:(=),expr,set≠:(=) ? Expr(:call,:($Sym.:+),expr,val) : val)

function declare_mutating_operations(M,neg,F,set_val)
for (op,set) ∈ ((:add,:(+=)),(:set,:(=)))
sm = Symbol(op,:multi!)
sb = Symbol(op,:blade!)
for (s,index) ∈ ((sm,:basisindex),(sb,:bladeindex))
for (i,B) ∈ ((:i,Bits),(:(bits(i)),Basis))
@eval begin
@inline function $s(out::$M,val::S,i::$B) where {M,T<:$F,S<:$F}
@inbounds $(set_val(set,:(out[$index(intlog(M),$i)]),:val))
return out
end
@inline function $s(out::Q,val::S,i::$B,::Dimension{N}) where Q<:$M where {M,T<:$F,S<:$F,N}
@inbounds $(set_val(set,:(out[$index(N,$i)]),:val))
return out
end
end
return m
end
end
for (prod,uct) ∈ ((:meet,:regressive),(:skew,:interior),(:cross,:crossprod))
for s ∈ (sm,sb)
@eval begin
@inline function $(Symbol(prod,s))(V::VectorSpace{N,D},m::MArray{Tuple{M},T,1,M},A::Bits,B::Bits,v::T) where {N,D,T<:Field,M}
if v ≠ 0
p,C,t = $uct(A,B,V)
t && $s(m,p ? -(v) : v,C,Dimension{N}())
@inline function $(Symbol(:join,s))(V::VectorSpace{N,D},m::$M,A::Bits,B::Bits,v::S) where {N,D,T<:$F,S<:$F,M}
if v ≠ 0 && !(hasdual(V) && isodd(A) && isodd(B))
$s(m,parity(A,B,V) ? $neg(v) : v,A ⊻ B,Dimension{N}())
end
return m
end
@inline function $(Symbol(prod,s))(m::MArray{Tuple{M},T,1,M},A::Basis{V},B::Basis{V},v::T) where {V,T<:Field,M}
if v ≠ 0
p,C,t = $uct(bits(A),bits(B),V)
t && $s(m,p ? -(v) : v,C,Dimension{N}())
@inline function $(Symbol(:join,s))(m::$M,v::S,A::Basis{V},B::Basis{V}) where {V,T<:$F,S<:$F,M}
if v ≠ 0 && !(hasdual(V) && hasdual(A) && hasdual(B))
$s(m,parity(A,B) ? $neg(v) : v,bits(A) ⊻ bits(B),Dimension{ndims(V)}())
end
return m
end
end
for (prod,uct) ∈ ((:meet,:regressive),(:skew,:interior),(:cross,:crossprod))
@eval begin
@inline function $(Symbol(prod,s))(V::VectorSpace{N,D},m::$M,A::Bits,B::Bits,v::T) where {N,D,T,M}
if v ≠ 0
p,C,t = $uct(A,B,V)
t && $s(m,p ? $neg(v) : v,C,Dimension{N}())
end
return m
end
@inline function $(Symbol(prod,s))(m::$M,A::Basis{V},B::Basis{V},v::T) where {V,T,M}
if v ≠ 0
p,C,t = $uct(bits(A),bits(B),V)
t && $s(m,p ? $neg(v) : v,C,Dimension{N}())
end
return m
end
end
end
end
end
end
Expand Down Expand Up @@ -377,6 +384,14 @@ end
### Product Algebra Constructor

function generate_product_algebra(Field=Field,MUL=:*,ADD=:+,SUB=:-,VEC=:mvec,CONJ=:conj)
if Field == Grassmann.Field
declare_mutating_operations(:(MArray{Tuple{M},T,1,M}),:-,Field,Expr)
elseif Field ∈ (SymField,:(SymPy.Sym))
declare_mutating_operations(:(SizedArray{Tuple{M},T,1,1}),SUB,Field,set_val)
end
Field == :(SymPy.Sym) && for par ∈ (:parany,:parval,:parsym)
@eval $par = ($par...,$Field)
end
TF = Field ≠ Number ? :Any : :T
EF = Field ≠ Any ? Field : ExprField
for Value ∈ MSV
Expand Down Expand Up @@ -997,6 +1012,7 @@ function generate_product_algebra(Field=Field,MUL=:*,ADD=:+,SUB=:-,VEC=:mvec,CON
end

generate_product_algebra()
generate_product_algebra(SymField,:($Sym.:*),:($Sym.:+),:($Sym.:-),:svec,:($Sym.conj))

const NSE = Union{Symbol,Expr,<:Number}

Expand Down
20 changes: 12 additions & 8 deletions src/multivectors.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,10 @@ abstract type TensorMixed{T,V} <: TensorAlgebra{V} end

import DirectSum: shift_indices, printindex, printindices, VTI

parany = (Expr,Any)
parsym = (Expr,Symbol)
parval = (Expr,)

## pseudoscalar

using LinearAlgebra
Expand Down Expand Up @@ -94,7 +98,7 @@ for Value ∈ MSV
$Value{V,G}(v,b::MValue{V,G}) where {V,G} = $Value{V,G,basis(b)}(v*b.v)
$Value{V,G,B}(v::T) where {V,G,B,T} = $Value{V,G,B,T}(v)
$Value{V}(v::T) where {V,T} = $Value{V,0,Basis{V}(),T}(v)
show(io::IO,m::$Value) = print(io,(valuetype(m)∉(Expr,Any) ? [m.v] : ['(',m.v,')'])...,basis(m))
show(io::IO,m::$Value) = print(io,(valuetype(m)∉parany ? [m.v] : ['(',m.v,')'])...,basis(m))
end
end

Expand Down Expand Up @@ -136,15 +140,15 @@ for (Blade,vector,Value) ∈ ((MSB[1],:MVector,MSV[1]),(MSB[2],:SVector,MSV[2]))

function show(io::IO, m::$Blade{T,V,G}) where {T,V,G}
ib = indexbasis(ndims(V),G)
@inbounds if T == Any && typeof(m.v[1]) ∈ (Expr,Symbol)
@inbounds typeof(m.v[1])≠Expr ? print(io,m.v[1]) : print(io,"(",m.v[1],")")
@inbounds if T == Any && typeof(m.v[1]) ∈ parsym
@inbounds typeof(m.v[1])∉parval ? print(io,m.v[1]) : print(io,"(",m.v[1],")")
else
@inbounds print(io,m.v[1])
end
@inbounds printindices(io,V,ib[1])
for k ∈ 2:length(ib)
if T == Any && typeof(m.v[k]) ∈ (Expr,Symbol)
@inbounds typeof(m.v[k])≠Expr ? print(io," + ",m.v[k]) : print(io," + (",m.v[k],")")
if T == Any && typeof(m.v[k]) ∈ parsym
@inbounds typeof(m.v[k])∉parval ? print(io," + ",m.v[k]) : print(io," + (",m.v[k],")")
else
@inbounds print(io,signbit(m.v[k]) ? " - " : " + ",abs(m.v[k]))
end
Expand Down Expand Up @@ -271,8 +275,8 @@ function show(io::IO, m::MultiVector{T,V}) where {T,V}
for k ∈ 1:length(ib)
@inbounds s = k+bs[i]
@inbounds if m.v[s] ≠ 0
@inbounds if T == Any && typeof(m.v[s]) ∈ (Expr,Symbol)
@inbounds typeof(m.v[s])≠Expr ? print(io," + ",m.v[s]) : print(io," + (",m.v[s],")")
@inbounds if T == Any && typeof(m.v[s]) ∈ parsym
@inbounds typeof(m.v[s])∉parval ? print(io," + ",m.v[s]) : print(io," + (",m.v[s],")")
else
@inbounds print(io,signbit(m.v[s]) ? " - " : " + ",abs(m.v[s]))
end
Expand All @@ -295,7 +299,7 @@ const VBV = Union{MValue,SValue,MBlade,SBlade,MultiVector}
@pure valuetype(::Union{MValue{V,G,B,T},SValue{V,G,B,T}} where {V,G,B}) where T = T
@pure valuetype(::TensorMixed{T}) where T = T
@inline value(::Basis,T=Int) = one(T)
@inline value(m::VBV,T::DataType=valuetype(m)) = T≠valuetype(m) ? convert(T,m.v) : m.v

This comment has been minimized.

Sign in to view

Copy link
@chakravala

chakravala Mar 4, 2019

Author Owner

strangely, I belive this was the root cause of the mysterious error #7, because conversions to Any seemed to cause trouble here...
@inline value(m::VBV,T::DataType=valuetype(m)) = T∉(valuetype(m),Any) ? convert(T,m.v) : m.v
@inline sig(::TensorAlgebra{V}) where V = V
@pure basis(m::Basis) = m
@pure basis(m::Union{MValue{V,G,B},SValue{V,G,B}}) where {V,G,B} = B
Expand Down
68 changes: 12 additions & 56 deletions src/symbolic.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,62 +2,18 @@
# This file is part of Grassmann.jl. It is licensed under the GPL license
# Grassmann Copyright (C) 2019 Michael Reed

const Sym = :(Reduce.Algebra)
const SymField = Any
import Reduce: RExpr

set_val(set,expr) = Expr(:(=),expr,set≠:(=) ? Expr(:call,:($Sym.:+),expr,:val) : :val)
*(a::RExpr,b::Basis{V}) where V = SValue{V}(a,b)
*(a::Basis{V},b::RExpr) where V = SValue{V}(b,a)
*(a::RExpr,b::MultiVector{T,V}) where {T,V} = MultiVector{promote_type(T,F),V}(broadcast(Reduce.Algebra.:*,a,b.v))
*(a::MultiVector{T,V},b::RExpr) where {T,V} = MultiVector{promote_type(T,F),V}(broadcast(Reduce.Algebra.:*,a.v,b))
*(a::RExpr,b::MultiGrade{V}) where V = MultiGrade{V}(broadcast(Reduce.Algebra.:*,a,b.v))
*(a::MultiGrade{V},b::RExpr) where V = MultiGrade{V}(broadcast(Reduce.Algebra.:*,a.v,b))
∧(a::RExpr,b::RExpr) = Reduce.Algebra.:*(a,b)
∧(a::RExpr,b::B) where B<:TensorTerm{V,G} where {V,G} = SValue{V,G}(a,b)
∧(a::A,b::RExpr) where A<:TensorTerm{V,G} where {V,G} = SValue{V,G}(b,a)

for (op,set) ∈ ((:add,:(+=)),(:set,:(=)))
sm = Symbol(op,:multi!)
sb = Symbol(op,:blade!)
for (s,index) ∈ ((sm,:basisindex),(sb,:bladeindex))
for (i,B) ∈ ((:i,Bits),(:(bits(i)),Basis))
@eval begin
@inline function $s(out::SizedArray{Tuple{M},T,1,1},val::T,i::$B) where {M,T<:SymField}
@inbounds $(set_val(set,:(out[$index(intlog(M),$i)])))
return out
end
@inline function $s(out::Q,val::T,i::$B,::Dimension{N}) where Q<:SizedArray{Tuple{M},T,1,1} where {M,T<:SymField,N}
@inbounds $(set_val(set,:(out[$index(N,$i)])))
return out
end
end
end
end
for s ∈ (sm,sb)
@eval begin
@inline function $(Symbol(:join,s))(V::VectorSpace{N,D},m::SizedArray{Tuple{M},T,1,1},A::Bits,B::Bits,v::T) where {N,D,T<:SymField,M}
if v ≠ 0 && !(Bool(D) && isodd(A) && isodd(B))
$s(m,parity(A,B,V) ? $Sym.:-(v) : v,A ⊻ B,Dimension{N}())
end
return m
end
@inline function $(Symbol(:join,s))(m::SizedArray{Tuple{M},T,1,1},A::Basis{V},B::Basis{V},v::T) where {V,T<:SymField,M}
if v ≠ 0 && !(hasdual(V) && hasdual(A) && hasdual(B))
$s(m,parity(A,B) ? $Sym.:-(v) : v,bits(A) ⊻ bits(B),Dimension{ndims(V)}())
end
return m
end
end
for (prod,uct) ∈ ((:meet,:regressive),(:skew,:interior),(:cross,:crossprod))
@eval begin
@inline function $(Symbol(prod,s))(V::VectorSpace{N,D},m::SizedArray{Tuple{M},T,1,1},A::Bits,B::Bits,v::T) where {N,D,T<:SymField,M}
if v ≠ 0 && !(hasdual(V) && hasdual(A) && hasdual(B))
p,C,t = $uct(A,B,V)
t && $s(m,p ? $Sym.:-(v) : v,C,Dimension{N}())
end
return m
end
@inline function $(Symbol(prod,s))(m::SizedArray{Tuple{M},T,1,1},A::Basis{V},B::Basis{V},v::T) where {V,T<:SymField,M}
if v ≠ 0 && !(hasdual(V) && hasdual(A) && hasdual(B))
p,C,t = $uct(bits(A),bits(B),V)
t && $s(m,p ? $Sym.:-(v) : v,C,Dimension{N}())
end
return m
end
end
end
end
for par ∈ (:parany,:parval,:parsym)
@eval $par = ($par...,RExpr)
end

generate_product_algebra(SymField,:($Sym.:*),:($Sym.:+),:($Sym.:-),:svec,:($Sym.conj))

0 comments on commit 1edc6b5

Please sign in to comment.