preliminary support for Unitful. Solves #68 (#72)

* preliminary support for Unitful. Solves #68

* test printing of unitful particles

* add Unitful note in the documentation

* include tests

.travis.yml

@@ -12,7 +12,7 @@ matrix:
   fast_finish: true
   include:
     - stage: "Documentation"
-      julia: 1.0
+      julia: 1.4
       os: linux
       script:
         - julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd()));

Project.toml

@@ -1,7 +1,7 @@
 name = "MonteCarloMeasurements"
 uuid = "0987c9cc-fe09-11e8-30f0-b96dd679fdca"
 authors = ["baggepinnen <[email protected]>"]
-version = "0.8.9"
+version = "0.8.10"
 
 [deps]
 Distributed = "8ba89e20-285c-5b6f-9357-94700520ee1b"
@@ -34,6 +34,7 @@ Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 SLEEFPirates = "476501e8-09a2-5ece-8869-fb82de89a1fa"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [targets]
-test = ["ControlSystems", "ForwardDiff", "Measurements", "SLEEFPirates", "Test", "Plots"]
+test = ["ControlSystems", "ForwardDiff", "Measurements", "SLEEFPirates", "Test", "Plots", "Unitful"]

docs/src/examples.md

@@ -19,6 +19,14 @@ The application we consider is optimization of a PID controller. Normally, we ar
 ## [Autodiff and Robust optimization](https://github.com/baggepinnen/MonteCarloMeasurements.jl/blob/master/examples/autodiff_robust_opt.jl)
 Another example using MonteCarloMeasurements to perform [robust optimization](https://en.wikipedia.org/wiki/Robust_optimization), this time with automatic differentiation. We use Optim.jl to solve a linear program with probabilistic constraints using 4 different methods, two gradient free, one first-order and one second-order method. We demonstrate calculation of gradients of uncertain functions with uncertain inputs using both Zygote.jl and ForwardDiff.jl.
 
+## Unitful interaction
+Particles with units can be created using the package [Unitful.jl](https://github.com/PainterQubits/Unitful.jl). The interaction is only supported through the construct `Particles{Quantity}`, whereas the reverse construct `Quantity{Particles}` is likely to result in problems. Unitful particles are thus created like this
+```@repl
+using MonteCarloMeasurements, Unitful # hide
+(1±0.1)u"V"
+(1..2)u"m"
+```
+
 
 ## Monte-Carlo sampling properties
 The variance introduced by Monte-Carlo sampling has some fortunate and some unfortunate properties. It decreases as 1/N, where N is the number of particles/samples. This unfortunately means that to get half the standard deviation in your estimate, you need to quadruple the number of particles. On the other hand, this variance does not depend on the dimension of the space, which is very fortunate.

src/MonteCarloMeasurements.jl

@@ -108,6 +108,7 @@ function __init__()
     @require ForwardDiff="f6369f11-7733-5829-9624-2563aa707210" include("forwarddiff.jl")
     @require SLEEFPirates="476501e8-09a2-5ece-8869-fb82de89a1fa" include("sleefpirates.jl")
     @require Measurements="eff96d63-e80a-5855-80a2-b1b0885c5ab7" include("measurements.jl")
+    @require Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" include("unitful.jl")
 end
 
 end

src/unitful.jl

@@ -0,0 +1,36 @@
+import .Unitful: Quantity, FreeUnits
+
+function to_num_str(p::AbstractParticles{T}, d=3) where T <: Quantity
+    s = std(p)
+    if T <: AbstractFloat && s < eps(p)
+        string(mean(p))
+    else
+        string(mean(p), " ± ", s)
+    end
+end
+
+for PT in ParticleSymbols
+
+    @eval begin
+        function Base.promote_rule(::Type{Quantity{S,D,U}}, ::Type{$PT{T,N}}) where {S, D, U, T, N}
+            NT = promote_type(Quantity{S,D,U},T)
+            $PT{NT,N}
+        end
+
+        function Base.convert(::Type{$PT{Quantity{S,D,U},N}}, y::Quantity) where {S, D, U, T, N}
+
+            $PT{Quantity{S,D,U},N}(fill(y, N))
+
+        end
+
+        function Base.:*(p::$PT{T,N}, y::Quantity{S,D,U}) where {S, D, U, T, N}
+            NT = promote_type(S,T)
+            $PT{Quantity{NT,D,U},N}(p.particles .* y)
+        end
+
+        function Base.:*(p::$PT, y::FreeUnits)
+            $PT(p.particles .* y)
+        end
+
+    end
+end

test/runtests.jl

@@ -691,6 +691,7 @@ Random.seed!(0)
         @inferred bymap(sin, 1 ± 2)
     end
 
+    include("test_unitful.jl")
     include("test_forwarddiff.jl")
     include("test_deconstruct.jl")
     include("test_sleefpirates.jl")

test/test_unitful.jl

@@ -0,0 +1,29 @@
+using Unitful
+function unitful_testfunction(Vi)
+    if Vi ≤ 0.0u"V"
+        return 0.0u"V"
+    elseif Vi ≥ 1.0u"V"
+        return 1.0u"V"
+    else
+        return Vi
+    end
+end
+
+register_primitive(unitful_testfunction) # must be outside testset
+
+@testset "Unitful" begin
+    @info "Testing Unitful"
+
+    for PT in (Particles, StaticParticles)
+        p1 = PT(100, Uniform(-0.5,1.5)) * 1u"V"
+        p2 = PT(100, Uniform(-0.5,1.5)) * u"V"
+        @test_nowarn println(p1)
+        @test_nowarn display(p1)
+
+        @test typeof(p1) == typeof(p2)
+
+        p3 = unitful_testfunction(p1)
+        @test extrema(p3) == (0.0u"V", 1.0u"V")
+    end
+
+end