Guide

Creating UnitfulScalars

Unitful scalars can be defined using the u-prefixed non-standard string literals, e. g.

Ry = 1.21104020061 * u"H*Mg*N*K*Pa*S*Tb*Tl*W/(Ba*C*Es*P*Pb*Pm*Ra*V*Yb)"
ħ = 1054.571818 * u"P/PP*pP*ppm*m*M/MM*Mm*mm"

They are parsed by Unitful.jl and then converted to a UnitfulScalar.

Don't try to use units defined on interval or logarithmic scales, such as °C or dBm. Use K and mW instead.

Creating UnitfulTensors

There are two ways to create UnitfulTensors. You can first construct an array of UnitfulScalars

[1u"s" 2u"s"
 3u"s" 4u"s"]

and then convert it to a UnitfulTensor

UnitfulTensor([1u"s" 2u"s"
               3u"s" 4u"s"])

Or you can construct an array of numbers and AxesDimensions separately and then combine them into a UnitfulTensor:

vals = [1. 2.
	3. 4.]
dims = 𝐓 * nodims(2, 2)
UnitfulTensor(vals, dims)

In the latter case, numerical values should be given in SI units.

Unitless quantities and arrays thereof don't have to be wrapped in a UnitfulScalar or UnitfulTensor, but trying to do mixed unitful/unitless arithmetic with something other than Arrays (e. g., sparse matrices or StaticArrays) can result in method ambiguities. They are easily fixed if necessary.

Creating AxesDimensions

AxesDimensions represent the physical dimensions of a UnitfulTensor and can be constructed as a tensor product of one-dimensional AxesDimensions:

AxesDimensions([𝐓, 𝐋, 𝐋, 𝐋]) ⊗ AxesDimensions(inv.([𝐓, 𝐋, 𝐋, 𝐋]))

or using nodims and scalar multiplication

𝐓 * nodims(2, 2)
# or dimensions(u"s") * nodims(2, 2)

for dimensionally homogeneous UnitfulTensors.

Elementwise operations

Elementwise mathematical functions (at least those that can be expanded in a Taylor series) preserve the tensor product structure of AxesDimensions, so this package modifies broadcasting to return UnitfulTensors when the broadcasted expression involves UnitfulScalars or UnitfulTensors:

julia> 1u"m" .* [1 2; 3 4]
2×2 UnitfulMatrix{Float64, SIDimensions, Matrix{Float64}, AxesDimensions{2, SIDimensions}}:
 1.0 m  2.0 m
 3.0 m  4.0 m

User-defined functions also can be used:

julia> A = UnitfulTensor([1u"m" 2u"m"; 3u"m" 4u"m"])
       B = UnitfulTensor([4u"s" 3u"s"; 2u"s" 1u"s"])
       f(x) = x + 2u"m^2/s"
       @. f(A^2 / B)
2×2 UnitfulMatrix{Float64, SIDimensions, Matrix{Float64}, AxesDimensions{2, SIDimensions}}:
 2.25 m^2 s^-1  3.33333 m^2 s^-1
  6.5 m^2 s^-1     18.0 m^2 s^-1

This works for usual mathematical functions from UnitfulScalars to UnitfulScalars. For something like f(x) = (dimensions(x) == NoDims) use map(f, A) instead.

Defining functions of UnitfulTensors

Functions of UnitfulTensors generally act on the numerical values and dimensions independently (with some exceptions, such as pivoted factorizations). If you use a function from some package (or Julia Base) that is not implemented in UnitfulTensors.jl, you can implement it yourself.

For example, reversing a UnitfulTensor along a specific dimension could be done like this:

using UnitfulTensors
import Base: reverse

function reverse(A::UnitfulTensor; kwargs...)
    dims = reverse(dimensions(A); kwargs...)
    vals = reverse(values(A); kwargs...)
    return UnitfulTensor(vals, dims)
end

function reverse(A::AxesDimensions; dims)
    d = dims
    scale, Adims... = dimsplat(A)
    newscale = scale * last(Adims[d]) # because AxisDimensions are normalized
    newdims = (Adims[1:d-1]..., reverse(Adims[d]), Adims[d+1:end]...)
    return AxesDimensions(newdims, newscale)
end

function reverse(A::AxisDimensions)
    dims = dimsvec(A)
    newdims = reverse(dims) ./ last(dims) # normalization
    return AxisDimensions(newdims)
end

Note that due to the tensor product structure of AxesDimensions dealing with dimensions is usually much easier than with numerical values.

Index notation

You can use @tensor and @tensoropt macros from TensorOperations.jl to write general tensor manipulations in index notation:

using UnitfulTensors, TensorOperations

A = UnitfulTensor([1 2u"s"; 3u"s^-1" 4])

@tensoropt B[i, j, k, l] := A[i, k] * A[l, j]
@tensoropt C[i, j, k, l] := A[i, i'] * A[j', j] * A[k', k] * A[l, l'] * B[i', j', k', l']
@tensoropt D[i, k] := C[i, j, k, j]

# output

2×2 UnitfulMatrix{Float64, SIDimensions, Matrix{Float64}, AxesDimensions{2, SIDimensions}}:
       5735.0  8370.0 s
 12555.0 s^-1   18290.0

Pitfalls

The types defined in this package store references to arrays of AbstractDimensions. Don't mutate them, or there will be surprises:

julia> Z = UnitfulTensor([1u"S" 2u"S"
                          3u"S" 4u"S"])
       V = UnitfulTensor([5u"V", 6u"V"])
       I = Z * V
2-element UnitfulVector{Float64, SIDimensions, Vector{Float64}, AxesDimensions{1, SIDimensions}}:
 17.0 A
 39.0 A
 
julia> dimsvec(normdims(Z)[1])[2] = dimensions(u"V")
       I
2-element UnitfulVector{Float64, SIDimensions, Vector{Float64}, AxesDimensions{1, SIDimensions}}:
           17.0 A
 39.0 kg m^2 s^-3

Make a copy instead if necessary.

Currently dimensions are compared with == for performance reasons, even though they are represented as floating-point numbers. As a consequence, units with integer, half-integer, quarter-integer exponents work, while units like s^(1/3) don't:

julia> 1u"s^(2/3)" + 1u"s" / 1u"s^(1/3)"
ERROR: DimensionMismatch: dimensions (𝐓^11184811/16777216, 𝐓^5592405/8388608) do not match

This will be fixed when I figure out how to switch to ≈ without a significant performance penalty.

Beyond SI

By default, UnitfulTensors.jl converts all quantities into SI units upon construction. These units are provided by the upreferred function from Unitful.jl. It is possible to override this behaviour and define custom units by mechanisms described in the Unitful.jl documentation.

Example:

import Unitful

const var"@uu_str" = Unitful.var"@u_str"
Unitful.@unit ħunit "ħ" ReducedPlanck Unitful.ħ false
Unitful.@unit eunit "e" ElementaryCharge Unitful.q false
Unitful.register(@__MODULE__)
_uT = ħunit/uu"meV"; _uL = uu"nm"; _uM = ħunit*_uT/_uL^2; _uI = eunit/_uT 
Unitful.preferunits(_uT, _uL, _uM, _uI)

using UnitfulTensors

σ = 1u"Ω^-1"; ω = 1u"THz2π"; E = 1u"eV"; vF = 1e6u"m/s"
σ, ω, E / vF

# output

(4108.235902227661 e^2 ħ^-1, 4.135667696923859 meV ħ^-1, 1.5192674478786259 ħ nm^-1)

However, this is not recommended unless really necessary. This feature has not been tested extensively and it is better to stick to SI for consistency.

Simplified UnitfulTensor construction

If you feel lazy to wrap all array literals in UnitfulTensor(...), you can call tensorize_literals(). After that, array literals containing UnitfulScalars will return UnitfulTensors instead of Arrays.

using UnitfulTensors

tensorize_literals()

[1u"s" 2
 3     4u"s^-1"]
 
# output

2×2 UnitfulMatrix{Float64, SIDimensions, Matrix{Float64}, AxesDimensions{2, SIDimensions}}:
 1.0 s       2.0
   3.0  4.0 s^-1