refined Hodge complement for null-basis

chakravala · Oct 4, 2019 · 4eab351 · 4eab351
1 parent 29bf98b
commit 4eab351
Show file tree

Hide file tree

Showing 3 changed files with 15 additions and 7 deletions.
diff --git a/src/algebra.jl b/src/algebra.jl
@@ -164,7 +164,7 @@ for side ∈ (:left,:right)
     c = Symbol(:complement,side)
     p = Symbol(:parity,side)
     @eval @pure function $c(b::Basis{V,G,B}) where {V,G,B}
-        d = getbasis(V,complement(ndims(V),B,diffvars(V)))
+        d = getbasis(V,complement(ndims(V),B,diffvars(V),hasinf(V)+hasorigin(V)))
         mixedmode(V)<0 && throw(error("Complement for mixed tensors is undefined"))
         typeof(V)<:Signature ? ($p(b) ? SBlade{V}(-value(d),d) : d) : SBlade{V}($p(b)*value(d),d)
     end
@@ -759,15 +759,15 @@ function generate_products(Field=Field,VEC=:mvec,MUL=:*,ADD=:+,SUB=:-,CONJ=:conj
             @eval begin
                 function $c(b::$Chain{T,V,G}) where {T<:$Field,V,G}
                     mixedmode(V)<0 && throw(error("Complement for mixed tensors is undefined"))
-                    $(insert_expr((:N,:ib,:D),VEC)...)
+                    $(insert_expr((:N,:ib,:D,:P),VEC)...)
                     out = zeros($VEC(N,G,T))
                     D = diffvars(V)
                     for k ∈ 1:binomial(N,G)
                         @inbounds val = b.v[k]
                         if val≠0
                             @inbounds p = $p(V,ib[k])
                             v = typeof(V)<:Signature ? (p ? $SUB(val) : val) : p*val
-                            @inbounds setblade!(out,v,complement(N,ib[k],D),Dimension{N}())
+                            @inbounds setblade!(out,v,complement(N,ib[k],D,P),Dimension{N}())
                         end
                     end
                     return $Chain{T,V,N-G}(out)
@@ -777,7 +777,7 @@ function generate_products(Field=Field,VEC=:mvec,MUL=:*,ADD=:+,SUB=:-,CONJ=:conj
         @eval begin
             function $c(m::MultiVector{T,V}) where {T<:$Field,V}
                 mixedmode(V)<0 && throw(error("Complement for mixed tensors is undefined"))
-                $(insert_expr((:N,:bs,:bn),VEC)...)
+                $(insert_expr((:N,:bs,:bn,:P),VEC)...)
                 out = zeros($VEC(N,T))
                 D = diffvars(V)
                 for g ∈ 1:N+1
@@ -786,7 +786,7 @@ function generate_products(Field=Field,VEC=:mvec,MUL=:*,ADD=:+,SUB=:-,CONJ=:conj
                         @inbounds val = m.v[bs[g]+i]
                         if val≠0
                             v = typeof(V)<:Signature ? ($p(V,ib[i]) ? $SUB(val) : val) : $p(V,ib[i])*val
-                            @inbounds setmulti!(out,v,complement(N,ib[i],D),Dimension{N}())
+                            @inbounds setmulti!(out,v,complement(N,ib[i],D,P),Dimension{N}())
                         end
                     end
                 end

diff --git a/src/parity.jl b/src/parity.jl
@@ -22,7 +22,10 @@ for side ∈ (:left,:right)
     p = Symbol(:parity,side)
     @eval begin
         @pure $p(V::Bits,B::Bits,N::Int) = $p(count_ones(V&B),sum(indices(B,N)),count_ones(B),N)
-        @pure $p(V::Signature,B,G=count_ones(B)) = $p(count_ones(value(V)&B),sum(indices(B,ndims(V))),G,ndims(V)-diffvars(V))
+        @pure function $p(V::Signature,B,G=count_ones(B))
+            o = hasorigin(V) && (hasinf(V) ? iszero(B&UInt(1))&(!iszero(B&UInt(2))) : isodd(B))
+            $p(count_ones(value(V)&B),sum(indices(B,ndims(V))),G,ndims(V)-diffvars(V))⊻o
+        end
         @pure function $p(V::DiagonalForm,B,G=count_ones(B))
             ind = indices(B,ndims(V))
             g = prod(V[ind])
@@ -32,7 +35,11 @@ for side ∈ (:left,:right)
     end
 end
 
-@pure complement(N::Int,B::UInt,D::Int=0)::UInt = ((~B)&(one(UInt)<<(N-D)-1))|(B&((one(UInt)<<D-1)<<(N-D)))
+@pure function complement(N::Int,B::UInt,D::Int=0,P::Int=0)::UInt
+    UP,ND = UInt(1)<<P-1, N-D
+    C = ((~B)&(UP⊻(UInt(1)<<ND-1)))|(B&(UP⊻((UInt(1)<<D-1)<<ND)))
+    count_ones(C&UP)>1 ? C⊻UP : C
+end
 
 ## product parities
 

diff --git a/src/utilities.jl b/src/utilities.jl
@@ -171,6 +171,7 @@ Base.@pure promote_type(t...) = Base.promote_type(t...)
     assign_expr!(e,x,:di,:(dualindex(V)))
     assign_expr!(e,x,:D,:(diffvars(V)))
     assign_expr!(e,x,:μ,:(diffvars(V)≠0))
+    assign_expr!(e,x,:P,:(hasinf(V)+hasorigin(V)))
     return x
 end