Merge pull request #14 from tpolakovic/sym

New pointgroup calc.

Project.toml

@@ -1,7 +1,7 @@
 name = "SimpleQuantum"
 uuid = "7c07e4e5-e402-48af-840a-a11b8ea3e28c"
 authors = ["Tomas Polakovic <[email protected]> and contributors"]
-version = "0.2.4"
+version = "0.2.5"
 
 [deps]
 Accessors = "7d9f7c33-5ae7-4f3b-8dc6-eff91059b697"
@@ -11,10 +11,10 @@ GLMakie = "e9467ef8-e4e7-5192-8a1a-b1aee30e663a"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
 Match = "7eb4fadd-790c-5f42-8a69-bfa0b872bfbf"
+MinkowskiReduction = "d617d09d-b3c0-43f1-8886-30d89b900f25"
 RangeHelpers = "3a07dd3d-1c52-4395-8858-40c6328157db"
 SplitApplyCombine = "03a91e81-4c3e-53e1-a0a4-9c0c8f19dd66"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
-SymmetryReduceBZ = "49a35663-c880-4242-bebb-1ec8c0fa8046"
 Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [compat]

src/SimpleQuantum.jl

@@ -10,7 +10,7 @@ using Colors
 using GLMakie
 using Accessors
 using Combinatorics
-using SymmetryReduceBZ
+using MinkowskiReduction
 
 include("crystal.jl")
 include("interface.jl")

src/crystal.jl

@@ -198,13 +198,46 @@ function mpspace(q)
     [(2r - q - 1)/2q for r ∈ 1:q]
 end
 
+function reduced_basis(c::Crystal{2})
+    vs = GaussReduce((eachcol(c.lattice.R))...)
+    hcat(vs...)
+end
+
+function reduced_basis(c::Crystal{3})
+    vs = minkReduce((eachcol(c.lattice.R))...)
+    length(vs) > 3 && return hcat(vs[1:3]...)
+    hcat(vs...)
+end
+
 """
     getPG(c::Crystal)
 
 Calculates the point symmetry group operations of `c`.
 """
 function getPG(c::Crystal)
-    SymmetryReduceBZ.Symmetry.calc_pointgroup(c.lattice.R)
+    R = reduced_basis(c) |> Matrix
+    D = ndims(c)
+    iR = inv(R)
+    norms = mapslices(norm, R; dims=1)
+    vol = c.lattice.V
+    ls = round.(norms)
+    verts = (Iterators.product(ntuple(i -> -ls[i]:ls[i], D)...)
+        .|> collect
+         |> x -> reshape(x, :, 1)
+         |> combinedims)[:,:,1]
+    out = SMatrix{D,D}[]
+    for perm ∈ permutations(1:size(verts, 2), D)
+            vs = R * verts[:, perm]
+            _norms = mapslices(norm, vs; dims=1)
+            _vol = abs(det(vs))
+            if all(norms ≈ _norms) & all(vol ≈ _vol)
+                op = SMatrix{D,D}(vs * iR)
+                if all(op' * op ≈ I)
+                    append!(out, [op])
+                end
+            end
+        end
+    out
 end
 
 """