Take advantage of new Quantity-KVector support (#130)

Project.toml

@@ -13,7 +13,7 @@ QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
 Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [compat]
-CliffordNumbers = "0.1.4"
+CliffordNumbers = "0.1.9"
 CoordRefSystems = "0.12, 0.13, 0.14, 0.15"
 FastGaussQuadrature = "1"
 HCubature = "1.5"

src/differentiation.jl

@@ -125,15 +125,6 @@ function differential(
         ts::Union{AbstractVector{T}, Tuple{T, Vararg{T}}},
         diff_method::DifferentiationMethod = _default_method(G)
 ) where {G <: Geometry, T <: AbstractFloat}
-    # Calculate the Jacobian, convert Vec -> KVector
-    J = jacobian(geometry, ts, diff_method)
-    J_kvecs = Iterators.map(_kvector, J)
-
-    # Extract units from Geometry type
-    Dim = Meshes.paramdim(geometry)
-    units = _units(geometry)^Dim
-
-    # Return norm of the exterior products
-    element = foldl(∧, J_kvecs)
-    return LinearAlgebra.norm(element) * units
+    J = Iterators.map(_KVector, jacobian(geometry, ts, diff_method))
+    return LinearAlgebra.norm(foldl(∧, J))
 end

src/utils.jl

@@ -46,11 +46,12 @@ _default_method(g::G) where {G <: Geometry} = _default_method(G)
 #                        CliffordNumbers and Units
 ################################################################################
 
-# Meshes.Vec -> ::CliffordNumber.KVector
-function _kvector(v::Meshes.Vec{Dim, T}) where {Dim, T}
+# Meshes.Vec -> Unitful.Quantity{CliffordNumber.KVector}
+function _KVector(v::Meshes.Vec{Dim, T}) where {Dim, T}
     ucoords = Iterators.map(Unitful.ustrip, v.coords)
-    return CliffordNumbers.KVector{1, VGA(Dim)}(ucoords...)
+    return CliffordNumbers.KVector{1, VGA(Dim)}(ucoords...) * _units(v)
 end
 
 # Extract the length units used by the CRS of a Geometry
 _units(::Geometry{M, CRS}) where {M, CRS} = Unitful.unit(CoordRefSystems.lentype(CRS))
+_units(::Meshes.Vec{Dim, T}) where {Dim, T} = Unitful.unit(T)

test/utils.jl

@@ -1,9 +1,9 @@
 @testitem "Utilities" setup=[Setup] begin
     using LinearAlgebra: norm
 
-    # _kvector
+    # _KVector
     v = Meshes.Vec(3, 4)
-    @test norm(MeshIntegrals._kvector(v)) ≈ 5.0
+    @test norm(MeshIntegrals._KVector(v)) ≈ 5.0u"m"
 
     # _units
     p = Point(1.0u"cm", 2.0u"mm", 3.0u"m")