diff --git a/Project.toml b/Project.toml index 1150ba8..acd4093 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "Grassmann" uuid = "4df31cd9-4c27-5bea-88d0-e6a7146666d8" authors = ["Michael Reed"] -version = "0.8.19" +version = "0.8.20" [deps] AbstractTensors = "a8e43f4a-99b7-5565-8bf1-0165161caaea" diff --git a/src/algebra.jl b/src/algebra.jl index 4603046..d2db780 100644 --- a/src/algebra.jl +++ b/src/algebra.jl @@ -35,7 +35,7 @@ Geometric algebraic product: ω⊖η = (-1)ᵖdet(ω∩η)⊗(Λ(ω⊖η)∪L(ω @pure ⟑(a::Submanifold{V},b::Submanifold{V}) where V = mul(a,b) *(a::X,b::Y,c::Z...) where {X<:TensorAlgebra,Y<:TensorAlgebra,Z<:TensorAlgebra} = *(a*b,c...) -@pure function mul(a::Submanifold{V},b::Submanifold{V},der=derive_mul(V,UInt(a),UInt(b),1,true)) where V +@pure function mul(a::Submanifold{V},b::Submanifold{V},der=derive_mul(V,UInt(a),UInt(b),1,true),::Val{field}=Val(false)) where {V,field} if isdiag(V) ba,bb = UInt(a),UInt(b) (diffcheck(V,ba,bb) || iszero(der)) && (return Zero(V)) @@ -938,18 +938,20 @@ adder(a,b,op=:+) = adder(typeof(a),typeof(b),op) return AntiSpinor{V}(out) end end end end - @noinline function product(a::Type{S},b::Type{<:Chain{V,G,T}},swap=false) where S<:TensorGraded{V,L} where {V,G,L,T} + @noinline function product(a::Type{S},b::Type{<:Chain{V,G,T}},swap=false,field=false) where S<:TensorGraded{V,L} where {V,G,L,T} MUL,VEC = mulvecs(a,b) + vfield = Val(field) anti = isodd(L) ≠ isodd(G) type = anti ? :AntiSpinor : :Spinor + args = field ? (:g,) : () if G == 0 return S<:Chain ? :(Chain{V,L}(broadcast($MUL,a.v,Ref(@inbounds b[1])))) : swap ? :(Single(b)⟑a) : :(a⟑Single(b)) elseif S<:Chain && L == 0 return :(Chain{V,G}(broadcast($MUL,Ref(@inbounds a[1]),b.v))) elseif (swap ? L : G) == mdims(V) && !istangent(V) - return swap ? (S<:Single ? :(⋆(~b)*value(a)) : S<:Chain ? :(@inbounds a[1]*⋆(~b)) : :(⋆(~b))) : :(@inbounds ⋆(~a)*b[1]) + return swap ? (S<:Single ? :(⋆(~b,$(args...))*value(a)) : S<:Chain ? :(@inbounds a[1]*⋆(~b,$(args...))) : :(⋆(~b,$(args...)))) : :(@inbounds ⋆(~a,$(args...))*b[1]) elseif (swap ? G : L) == mdims(V) && !istangent(V) - return swap ? :(b[1]*complementlefthodge(~a)) : S<:Single ? :(value(a)*complementlefthodge(~b)) : S<:Chain ? :(@inbounds a[1]*complementlefthodge(~b)) : :(complementlefthodge(~b)) + return swap ? :(b[1]*complementlefthodge(~a,$(args...))) : S<:Single ? :(value(a)*complementlefthodge(~b,$(args...))) : S<:Chain ? :(@inbounds a[1]*complementlefthodge(~b,$(args...))) : :(complementlefthodge(~b,$(args...))) elseif binomial(mdims(V),G)*(S<:Chain ? binomial(mdims(V),L) : 1)<(1<(G+mingrade(b),mdims(V)))) && (!istangent(V)) return Zero(V) @@ -1255,18 +1259,18 @@ for (op,product) ∈ ((:∧,:exteradd),(:*,:geomadd), if (swap ? maxgrade(b)>> : :⊘)(a,multispin(b),$(args...)))) else - :(isodd(grade(BB))≠isodd(N) && return :(a⊘multispin(b))) + :(isodd(grade(BB))≠isodd(N) && return :($(swap ? :>>> : :⊘)(a,multispin(b),$(args...)))) end) inspin = iseven(grade(BB)) VECS = isodd(G) ? VEC : string(VEC)*"s" @@ -1537,14 +1544,14 @@ for input ∈ (:Couple,:PseudoCouple) if S<:Chain @inbounds for j ∈ 1:bn[G+1] @inbounds B = ib[j] - $preproduct!(V,out,A,B,derive_pre(V,A,B,val,:(@inbounds a[$j]),MUL)) + $preproduct!(V,out,A,B,derive_pre(V,A,B,val,:(@inbounds a[$j]),MUL),vfield) end else B = UInt(basis(a)) if S<:Single - $preproduct!(V,out,A,B,derive_pre(V,A,B,val,:(a.v),MUL)) + $preproduct!(V,out,A,B,derive_pre(V,A,B,val,:(a.v),MUL),vfield) else - $preproduct!(V,out,A,B,derive_pre(V,A,B,val,false)) + $preproduct!(V,out,A,B,derive_pre(V,A,B,val,false),vfield) end end end @@ -1556,7 +1563,7 @@ for input ∈ (:Couple,:PseudoCouple) @inbounds val = out[bs2[g]+i] !isnull(val) && for (B,val2) ∈ ((UInt(BB),(swap ? parityclifford(grade(BB)) : false) ? :(-value(imaginary(b))) : :(value(imaginary(b)))),(indexbasis(N,$pg)[1],(swap ? parityclifford($pg) : false) ? :(-value($$calar(b))) : :(value($$calar(b))))) @inbounds A = ia[i] - $preproduct2!(V,out2,A,B,derive_pre(V,A,B,val,val2,MUL)) + $preproduct2!(V,out2,A,B,derive_pre(V,A,B,val,val2,MUL),vfield) end end end @@ -1608,8 +1615,9 @@ for com ∈ (:spinor,:s_m,:m_s,:anti,:a_m,:m_a,:multivector,:s_a,:a_s) left,leftspin,right,rightspin,br = leftrightsym(com) VEC = outspin ? :svecs : :svec genloop = Symbol(:generate_loop_,com) - @eval @noinline function $genloop(V,a,b,t,MUL,product!,preproduct!,d=nothing) + @eval @noinline function $genloop(V,a,b,t,MUL,product!,preproduct!,field=false,d=nothing) $(insert_expr((:N,br...,:bn),:mvec)...) + vfield = Val(field) if mdims(V)diffmode(V) end -@pure function parityinterior(V::M,a,b,::Val{lim}=Val(false)) where {W,M<:Manifold{W},lim} +@pure function parityinterior(V::M,a,b,::Val{lim}=Val(false),::Val{field}=Val(false)) where {W,M<:Manifold{W},lim,field} A,B,Q,Z = symmetricmask(V,a,b); N = rank(V) diffcheck(V,A,B) && (return lim ? (Values{0,Tuple{UInt,Int}}(),Z) : (1,Zero(UInt),false,Z)) gout,tout,G,diag = 0,false,count_ones(B),isdiag(V) bas = diag ? Values((B,)) : indexbasis(grade(V),G) g = if diag - gg = V[indices(B,N)] - Values((prod(isempty(gg) ? 1 : signbool.(gg)),)) + if field + bi = basisindex(N,B) + Values((:(value(value(g))[$bi]),)) + else + gg = V[indices(B,N)] + Values((prod(isempty(gg) ? 1 : signbool.(gg)),)) + end else value(metrictensor(V,G))[bladeindex(grade(V),B)] end @@ -124,7 +129,7 @@ end p,C,t = parityregressive(Signature(W),A,complement(N,bas[i],diffvars(V)),Val{true}()) CQ,tout = C|Q,tout|t if t - ggg = (p⊻parityright(0,sum(indices(bas[i],N)),G)) ? -(g[i]) : g[i] + ggg = (p⊻parityright(0,sum(indices(bas[i],N)),G)) ? field ? :(-($(g[i]))) : -(g[i]) : g[i] if CQ ∉ bs bs = (bs...,CQ) bg = (bg...,(CQ,ggg)) @@ -132,7 +137,7 @@ end j = findfirst(x->x==CQ,bs) bg = @inbounds (bg[1:j-1]...,(CQ,(bg[j][2]+ggg)),bg[j+1:length(bs)]...) end - gout += ggg + iszero(gout) ? (gout = ggg) : (gout += ggg) end end end @@ -149,16 +154,21 @@ end parity(V,A,B) ? -1 : 1 end -@pure function parityinner(V::M,a::UInt,b::UInt) where {W,M<:Manifold{W}} +@pure function parityinner(V::M,a::UInt,b::UInt,::Val{field}=Val(false)) where {W,M<:Manifold{W},field} A,B = symmetricmask(V,a,b) C = A&B; G = count_ones(C) g = if isdiag(V) || hasconformal(V) - gg = V[indices(C,mdims(V))] - abs(prod(isempty(gg) ? 1 : signbool.(gg))) + if field + bi = basisindex(mdims(V),C) + :(value(value(g))[$bi]) + else + gg = V[indices(C,mdims(V))] + abs(prod(isempty(gg) ? 1 : signbool.(gg))) + end else value(metrictensor(V,G))[bladeindex(mdims(V),C)] end - parity(Signature(W),A,B) ? -(g) : g + parity(Signature(W),A,B) ? fieldneg(g) : g end function parityseq(V,B) @@ -172,85 +182,89 @@ end maxempty(x) = isempty(x) ? 0 : maximum(x) -function paritygeometric(V,A::UInt,B::UInt) +function paritygeometric(V,A::UInt,B::UInt,field=Val(false)) #if iszero(B) # Values(((A,1),)) a,b = splitbasis(V,A),splitbasis(V,B) ga,gb = maxempty(count_ones.(a)),maxempty(count_ones.(b)) if (ga≤1 && gb≤1) ? count_ones(A) ≥ count_ones(B) : ga ≥ gb - paritygeometric(V,A,b) + paritygeometric(V,A,b,field) else - paritygeometric(V,a,B) + paritygeometric(V,a,B,field) end end -function paritygeometric(V,A::UInt,b::Tuple) +function paritygeometric(V,A::UInt,b::Tuple,field) i = length(b) vals = if iszero(i) Values((((UInt(0),1),(A,1)),)) else p = parityseq(V,b) - out = paritygeometricright(V,((UInt(0),p),(A,1)),b[1]) - isone(i) ? out : paritygeometricright(V,out,b,2) + out = paritygeometricright(V,((UInt(0),p),(A,1)),b[1],field) + isone(i) ? out : paritygeometricright(V,out,b,2,field) end combinebasis(V,vals) end -function paritygeometric(V,a::Tuple,B::UInt) +function paritygeometric(V,a::Tuple,B::UInt,field) i = length(a) vals = if iszero(i) Values((((UInt(0),1),(B,1)),)) else p = parityseq(V,a) - out = paritygeometricleft(V,a[i],((UInt(0),p),(B,1))) - isone(i) ? out : paritygeometricleft(V,a,out,i-1) + out = paritygeometricleft(V,a[i],((UInt(0),p),(B,1)),field) + isone(i) ? out : paritygeometricleft(V,a,out,i-1,field) end combinebasis(V,vals) end -function paritygeometricright(V,a,b::Tuple,i) - out = vcat(paritygeometricright.(Ref(V),a,Ref(b[i]))...) - i==length(b) ? out : paritygeometricright(V,out,b,i+1) +function paritygeometricright(V,a,b::Tuple,i,field) + out = vcat(paritygeometricright.(Ref(V),a,Ref(b[i]),field)...) + i==length(b) ? out : paritygeometricright(V,out,b,i+1,field) end -function paritygeometricright(V,aeai,B::UInt) +function paritygeometricright(V,aeai,B::UInt,field) ae,ai = @inbounds (aeai[1],aeai[2]) Ae,Aeg = @inbounds (ae[1],ae[2]) Ai,Aig = @inbounds (ai[1],ai[2]) G = count_ones(B) CCg = if isdiag(V) || hasconformal(V) - g,C,t,Z = interior(V,Ai,B) - Cg = parityclifford(G)⊻isodd(G*count_ones(Ai)) ? -Aig*g : Aig*g + g,C,t,Z = interior(V,Ai,B,Val(false),field) + Aigg = fieldprod(Aig,g) + Cg = parityclifford(G)⊻isodd(G*count_ones(Ai)) ? fieldneg(Aigg) : Aigg t ? ((ae,(C,Cg)),) : () else - Cg,Z = parityinterior(V,Ai,B,Val(true)) - g = parityclifford(G)⊻isodd(G*count_ones(Ai)) ? -Aig*Aeg : Aig*Aeg + Cg,Z = parityinterior(V,Ai,B,Val(true),field) + AigAeg = fieldprod(Aig,Aeg) + g = parityclifford(G)⊻isodd(G*count_ones(Ai)) ? fieldneg(AigAeg) : AigAeg ([((Ae,g),cg) for cg ∈ Cg]...,) end out = if iszero(Ai&B) - p = parity(grade(V),Ai⊻Ae,B) ? -Aeg : Aeg + p = parity(grade(V),Ai⊻Ae,B) ? fieldneg(Aeg) : Aeg (combinegeometric(V,(Ae⊻B,p),ai),CCg...) else CCg end isempty(out) ? Values{0,Tuple{Tuple{UInt,Int},Tuple{UInt,Int}}}() : Values(out) end -function paritygeometricleft(V,a::Tuple,b,i) - out = vcat(paritygeometricleft.(Ref(V),Ref(a[i]),b)...) - isone(i) ? out : paritygeometricleft(V,a,out,i-1) +function paritygeometricleft(V,a::Tuple,b,i,field) + out = vcat(paritygeometricleft.(Ref(V),Ref(a[i]),b,field)...) + isone(i) ? out : paritygeometricleft(V,a,out,i-1,field) end -function paritygeometricleft(V,A::UInt,bebi) +function paritygeometricleft(V,A::UInt,bebi,field) be,bi = @inbounds (bebi[1],bebi[2]) Be,Beg = @inbounds (be[1],be[2]) Bi,Big = @inbounds (bi[1],bi[2]) G = count_ones(A) CCg = if isdiag(V) || hasconformal(V) - g,C,t,Z = interior(V,Bi,A) - Cg = parityreverse(G) ? -Big*g : Big*g + g,C,t,Z = interior(V,Bi,A,Val(false),field) + Bigg = fieldprod(Big,g) + Cg = parityreverse(G) ? fieldneg(Bigg) : Bigg t ? ((be,(C,Cg)),) : () else - Cg,Z = parityinterior(V,Bi,A,Val(true)) - g = parityreverse(G) ? -Big*Beg : Big*Beg + Cg,Z = parityinterior(V,Bi,A,Val(true),field) + BigBeg = fieldprod(Big,Beg) + g = parityreverse(G) ? fieldneg(BigBeg) : BigBeg ([((Be,g),cg) for cg ∈ Cg]...,) end out = if iszero(A&Bi) - p = parity(grade(V),A,Be⊻Bi) ? -Beg : Beg + p = parity(grade(V),A,Be⊻Bi) ? fieldneg(Beg) : Beg (combinegeometric(V,(A⊻Be,p),bi),CCg...) else CCg @@ -309,15 +323,22 @@ function combinegeometric(V,ei) e,i = @inbounds ei[1],ei[2] E,Eg = @inbounds e[1],e[2] I,Ig = @inbounds i[1],i[2] - iszero(E&I) ? (E⊻I,Eg*Ig) : nothing + iszero(E&I) ? (E⊻I,fieldprod(Eg,Ig)) : nothing end function combinegeometric(V,e,i) E,Eg = @inbounds e[1],e[2] I,Ig = @inbounds i[1],i[2] #p = parity(grade(V),E,I) ? -Eg*Ig : Eg*Ig - return ((E,Eg*Ig),(I,1)) + return ((E,fieldprod(Eg,Ig)),(I,1)) end +fieldprod(a::Expr,b::Expr,op=:*) = :($op($a,$b)) +fieldprod(a::Expr,b,op=:*) = isone(b) ? a : isone(abs(b)) ? :(-($a)) : :($op($a,$b)) +fieldprod(a,b::Expr,op=:*) = isone(a) ? b : isone(abs(a)) ? :(-($b)) : :($op($a,$b)) +fieldprod(a,b) = a*b +fieldneg(a::Expr) = :(-($a)) +fieldneg(a) = -a + ### parity cache const parity_cache = Dict{UInt,Vector{Vector{Bool}}}[] @@ -360,6 +381,10 @@ end ### parity product caches +function interior(V,a,b,c::Val{lim},d::Val{field}=Val(false)) where {lim,field} + lim||field ? parityinterior(V,a,b,c,d) : interior(V,a,b) +end + for par ∈ (:conformal,:regressive,:interior) calc = Symbol(:parity,par) T = Tuple{Any,UInt,Bool,UInt} diff --git a/src/products.jl b/src/products.jl index 4d0eef9..041f7fa 100644 --- a/src/products.jl +++ b/src/products.jl @@ -196,7 +196,7 @@ function generate_mutators(M,F,set_val,SUB,MUL) end return false end - @inline function $(Symbol(:join,spre))(V,m::$M,a::UInt,b::UInt,v::S) where {T<:$F,S<:$F,M} + @inline function $(Symbol(:join,spre))(V,m::$M,a::UInt,b::UInt,v::S,field=nothing) where {T<:$F,S<:$F,M} if v ≠ 0 && !diffcheck(V,a,b) A,B,Q,Z = symmetricmask(V,a,b) val = :($$MUL($(parityinner(grade(V),A,B)),$v)) @@ -234,12 +234,12 @@ function generate_mutators(M,F,set_val,SUB,MUL) end return false end - @inline function $(Symbol(:geom,spre))(V,m::$M,a::UInt,b::UInt,v::S) where {T<:$F,S<:$F,M} + @inline function $(Symbol(:geom,spre))(V,m::$M,a::UInt,b::UInt,v::S,vfield::Val{field}=Val(false)) where {T<:$F,S<:$F,M,field} if v ≠ 0 && !diffcheck(V,a,b) A,B,Q,Z = symmetricmask(V,a,b) if isdiag(V) pcc,bas,cc = (hasinf(V) && hasorigin(V)) ? conformal(V,A,B) : (false,A⊻B,false) - g = parityinner(V,A,B) + g = parityinner(V,A,B,vfield) val = :($$MUL($g,$(pcc ? :($$SUB($v)) : v))) if istangent(V) !iszero(Z) && (val = Expr(:call,:*,val,getbasis(loworder(V),Z))) @@ -248,7 +248,7 @@ function generate_mutators(M,F,set_val,SUB,MUL) $spre(m,val,bas|Q,Val(mdims(V))) cc && $spre(m,hasinforigin(V,A,B) ? :($$SUB($val)) : val,(conformalmask(V)⊻bas)|Q,Val(mdims(V))) else - for (bas,g) ∈ paritygeometric(V,A,B) + for (bas,g) ∈ paritygeometric(V,A,B,vfield) val = :($$MUL($g,$v)) if istangent(V) !iszero(Z) && (val = Expr(:call,:*,val,getbasis(loworder(V),Z))) @@ -267,63 +267,93 @@ function generate_mutators(M,F,set_val,SUB,MUL) $(Symbol(prod,S))(V,m,UInt(A),UInt(B),v) end end - @eval begin - @inline function $(Symbol(prod,s))(V,m::$M,A::UInt,B::UInt,val::T) where {T,M} - if val ≠ 0 - if isdiag(V) - g,C,t,Z = $uct(V,A,B) - v = val - if istangent(V) && !iszero(Z) - T≠Any && (return true) - _,_,Q,_ = symmetricmask(V,A,B) - v *= getbasis(loworder(V),Z) - count_ones(Q)+order(v)>diffmode(V) && (return false) - end - t && $s(m,$MUL(g,v),C,Val(mdims(V))) - else - Cg,Z = $(Symbol(:parity,uct))(V,A,B,Val(true)) - v = val - if istangent(V) && !iszero(Z) - T≠Any && (return true) - _,_,Q,_ = symmetricmask(V,A,B) - v *= getbasis(loworder(V),Z) - count_ones(Q)+order(v)>diffmode(V) && (return false) - end - for (C,g) ∈ Cg - $s(m,$MUL(g,v),C,Val(mdims(V))) - end + end + @eval begin + @inline function $(Symbol(:meet,s))(V,m::$M,A::UInt,B::UInt,val::T) where {T,M} + if val ≠ 0 + g,C,t,Z = regressive(V,A,B) + v = val + if istangent(V) && !iszero(Z) + T≠Any && (return true) + _,_,Q,_ = symmetricmask(V,A,B) + v *= getbasis(loworder(V),Z) + count_ones(Q)+order(v)>diffmode(V) && (return false) + end + t && $s(m,$MUL(g,v),C,Val(mdims(V))) + end + return false + end + @inline function $(Symbol(:meet,spre))(V,m::$M,A::UInt,B::UInt,val::T,field::Val=Val(false)) where {T,M} + if val ≠ 0 + g,C,t,Z = regressive(V,A,B) + v = val + if istangent(V) && !iszero(Z) + _,_,Q,_ = symmetricmask(V,A,B) + v = Expr(:call,:*,v,getbasis(loworder(V),Z)) + v = :(h=$v;iszero(h)||$(count_ones(Q))+order(h)>$(diffmode(V)) ? 0 : h) + end + t && $spre(m,Expr(:call,$(QuoteNode(MUL)),g,v),C,Val(mdims(V))) + end + return false + end + @inline function $(Symbol(:skew,s))(V,m::$M,A::UInt,B::UInt,val::T) where {T,M} + if val ≠ 0 + if isdiag(V) + g,C,t,Z = interior(V,A,B) + v = val + if istangent(V) && !iszero(Z) + T≠Any && (return true) + _,_,Q,_ = symmetricmask(V,A,B) + v *= getbasis(loworder(V),Z) + count_ones(Q)+order(v)>diffmode(V) && (return false) + end + t && $s(m,$MUL(g,v),C,Val(mdims(V))) + else + Cg,Z = parityinterior(V,A,B,Val(true),false) + v = val + if istangent(V) && !iszero(Z) + T≠Any && (return true) + _,_,Q,_ = symmetricmask(V,A,B) + v *= getbasis(loworder(V),Z) + count_ones(Q)+order(v)>diffmode(V) && (return false) + end + for (C,g) ∈ Cg + $s(m,$MUL(g,v),C,Val(mdims(V))) end end - return false end - @inline function $(Symbol(prod,spre))(V,m::$M,A::UInt,B::UInt,val::T) where {T,M} - if val ≠ 0 - if isdiag(V) - g,C,t,Z = $uct(V,A,B) - v = val - if istangent(V) && !iszero(Z) - _,_,Q,_ = symmetricmask(V,A,B) - v = Expr(:call,:*,v,getbasis(loworder(V),Z)) - v = :(h=$v;iszero(h)||$(count_ones(Q))+order(h)>$(diffmode(V)) ? 0 : h) - end - t && $spre(m,Expr(:call,$(QuoteNode(MUL)),g,v),C,Val(mdims(V))) - else - Cg,Z = $(Symbol(:parity,uct))(V,A,B,Val(true)) - v = val - if istangent(V) && !iszero(Z) - _,_,Q,_ = symmetricmask(V,A,B) - v = Expr(:call,:*,v,getbasis(loworder(V),Z)) - v = :(h=$v;iszero(h)||$(count_ones(Q))+order(h)>$(diffmode(V)) ? 0 : h) - end - for (C,g) ∈ Cg - $spre(m,Expr(:call,$(QuoteNode(MUL)),g,v),C,Val(mdims(V))) - end + return false + end + @inline function $(Symbol(:skew,spre))(V,m::$M,A::UInt,B::UInt,val::T,field::Val=Val(false)) where {T,M} + if val ≠ 0 + if isdiag(V) + g,C,t,Z = interior(V,A,B,Val(false),field) + v = val + if istangent(V) && !iszero(Z) + _,_,Q,_ = symmetricmask(V,A,B) + v = Expr(:call,:*,v,getbasis(loworder(V),Z)) + v = :(h=$v;iszero(h)||$(count_ones(Q))+order(h)>$(diffmode(V)) ? 0 : h) + end + t && $spre(m,Expr(:call,$(QuoteNode(MUL)),g,v),C,Val(mdims(V))) + else + Cg,Z = parityinterior(V,A,B,Val(true),field) + v = val + if istangent(V) && !iszero(Z) + _,_,Q,_ = symmetricmask(V,A,B) + v = Expr(:call,:*,v,getbasis(loworder(V),Z)) + v = :(h=$v;iszero(h)||$(count_ones(Q))+order(h)>$(diffmode(V)) ? 0 : h) + end + for (C,g) ∈ Cg + $spre(m,Expr(:call,$(QuoteNode(MUL)),g,v),C,Val(mdims(V))) end end - return false end + return false end + + end + end end end @@ -1089,16 +1119,21 @@ end ### Product Algebra +const product_contraction_metric = product_contraction @generated function contraction2(a::TensorGraded{V,L},b::Chain{V,G,T}) where {V,G,L,T} - product_contraction(a,b,false,:product) + product_contraction(a,b,false,false,:product) end -for (op,prop) ∈ ((:⟑,:product),(:contraction,:product_contraction)) +@generated function contraction2_metric(a::TensorGraded{V,L},b::Chain{V,G,T},g) where {V,G,L,T} + product_contraction(a,b,false,true,:product) +end +for (op,prop,field) ∈ ((:⟑,:product,false),(:contraction,:product_contraction,false),(:wedgedot_metric,:product,true),(:contraction_metric,:product_contraction,true)) + args = field ? (:g,) : () @eval begin - @generated function $op(b::Chain{V,G,T},a::TensorTerm{V,L}) where {V,G,L,T} - $prop(a,b,true) + @generated function $op(b::Chain{V,G,T},a::TensorTerm{V,L},$(args...)) where {V,G,L,T} + $prop(a,b,true,$field) end - @generated function $op(a::TensorGraded{V,L},b::Chain{V,G,T}) where {V,G,L,T} - $prop(a,b) + @generated function $op(a::TensorGraded{V,L},b::Chain{V,G,T},$(args...)) where {V,G,L,T} + $prop(a,b,false,$field) end end end @@ -1106,30 +1141,32 @@ for op ∈ (:∧,:∨) prop = Symbol(:product_,op) @eval begin @generated function $op(a::Chain{w,G,T},b::Chain{W,L,S}) where {T,w,S,W,G,L} - $prop(a,b) + $prop(a,b,false) end @generated function $op(b::Chain{Q,G,T},a::TensorTerm{R,L}) where {Q,G,T,R,L} $prop(a,b,true) end @generated function $op(a::TensorTerm{Q,G},b::Chain{R,L,T}) where {Q,R,T,G,L} - $prop(a,b) + $prop(a,b,false) end end end -for (op,product!) ∈ ((:∧,:exteraddmulti!),(:⟑,:geomaddmulti!), - (:∨,:meetaddmulti!),(:contraction,:skewaddmulti!)) +for (op,product!,field) ∈ ((:∧,:exteraddmulti!,false),(:⟑,:geomaddmulti!,false), + (:∨,:meetaddmulti!,false),(:contraction,:skewaddmulti!,false), + (:wedgedot_metric,:geomaddmulti!,true),(:contraction_metric,:skewaddmulti!,true)) preproduct! = Symbol(product!,:_pre) - prop = op≠:⟑ ? Symbol(:product_,op) : :product + prop = op∉(:⟑,:wedgedot_metric) ? Symbol(:product_,op) : :product + args = field ? (:g,) : () @eval begin - @generated function $op(b::Multivector{V,T},a::TensorGraded{V,G}) where {V,T,G} - $prop(a,b,true) + @generated function $op(b::Multivector{V,T},a::TensorGraded{V,G},$(args...)) where {V,T,G} + $prop(a,b,true,$field) end - @generated function $op(a::TensorGraded{V,G},b::Multivector{V,S}) where {V,G,S} - $prop(a,b) + @generated function $op(a::TensorGraded{V,G},b::Multivector{V,S},$(args...)) where {V,G,S} + $prop(a,b,false,$field) end - @generated function $op(a::Multivector{V,T},b::Multivector{V,S}) where {V,T,S} + @generated function $op(a::Multivector{V,T},b::Multivector{V,S},$(args...)) where {V,T,S} MUL,VEC = mulvec(a,b) - loop = generate_loop_multivector(V,:(a.v),:(b.v),promote_type(T,S),MUL,$product!,$preproduct!) + loop = generate_loop_multivector(V,:(a.v),:(b.v),promote_type(T,S),MUL,$product!,$preproduct!,$field) if mdims(V)