40
$\begingroup$

This generates 75 unitized pairs for a plot that displays physical units:

x = UnitConvert@Quantity["PlanckConstant" "SpeedOfLight" / "BoltzmannConstant"]
x1 = UnitConvert@Quantity[2,"PlanckConstant" ("SpeedOfLight")^2]
B[L_,T_]:= (L^-5) /(Exp[x/(L T)]-1)
p = Table[{Quantity[ll, "Micrometers"], UnitConvert[
 x1 B[Quantity[ll, "Micrometers"], Quantity[1000, "Kelvins"]], 
 "Watts"/("Meters")^3]}, {ll, 0.1, 15, 0.2}]; //Timing

Then plot it

ListLinePlot[p, AxesLabel -> Automatic]

to get nice unitized labels.

However, this naive coding to generate "p" takes 27 seconds on my core-duo 2.4 GHz laptop under windows 7 32-bit, ridiculously slow for 75 pairs! Can anyone accelerate this, say tenfold? Not too tricky please because I need transparency for Mathematica beginners.

$\endgroup$
8
  • $\begingroup$ Confirmed slow: 24 sec on Mac OS X v 10.7.5 on 2.8 GHz Intel Core i7 with 8 GB 1067 MHz DDR3 memory. $\endgroup$
    – DavidC
    Commented Dec 17, 2012 at 2:54
  • $\begingroup$ @DavidCarraher Got about the same in timing on my MacBook. What makes it so slow is calling the function B with the units. Even listing the result p is rather slow. $\endgroup$
    – Matariki
    Commented Dec 17, 2012 at 3:05
  • $\begingroup$ @DavidCarraher Roughly the same here on a similar Mac. Good question (+1) - that looks really weird! $\endgroup$
    – Jens
    Commented Dec 17, 2012 at 3:06
  • 2
    $\begingroup$ The code runs without an internet connection, so connections to Wolfram Alpha are not to blame. $\endgroup$
    – DavidC
    Commented Dec 17, 2012 at 3:17
  • 3
    $\begingroup$ My experience is that the problem lies in using the Quantity function. That is horribly slow. What I do is to create units that will Work in the pre-9 environment, i.e. Meter=Quantity["Meter"] etc. It takes 19 seconds to create the 19 SI-units this way (in a loop) on my computer with 2-Intel 7 (2.8 GHz) and 8GB memory. That is unacceptable. $\endgroup$
    – user6098
    Commented Feb 25, 2013 at 18:15

4 Answers 4

35
$\begingroup$

If we pre-evaluate the expression:

x1 B[Quantity[ll, "Micrometers"], Quantity[1000, "Kelvins"]]

then the run time can be reduced by about factor 10. We can do this by hoisting the expression out of the loop and pre-evaluating using With:

p =
  Module[{ll}
  , With[{v = x1 B[Quantity[ll, "Micrometers"], Quantity[1000, "Kelvins"]]}
      , Table[{Quantity[ll, "Micrometers"], UnitConvert[v, "Watts"/("Meters")^3]}
        , {ll, 0.1, 15, 0.2}
        ]
    ]
  ]; //Timing

(* {2.324415, Null} *)

What's Going On?

I tried this example under Windows 7, 64-bit.

I used TraceScan to count the number of evaluation steps:

Module[{n = 0}
, TraceScan[
    n++ &
  , p = Table[{Quantity[ll, "Micrometers"], UnitConvert[
    x1 B[Quantity[ll, "Micrometers"], Quantity[1000, "Kelvins"]], 
    "Watts"/("Meters")^3]}, {ll, 0.1, 15, 0.2}];
  ]
; n
]

(* 1913468 *)

1,913,468 evaluation steps?! A unitless analog to the expression required only 4,136 steps:

Module[{n = 0}
, TraceScan[
    n++ &
  , p = Table[{1, x1 B[ll, 1]}, {ll, 0.1, 15, 0.2}];
  ]
; n
]

(* 4136 *)

TracePrint applied to the expression produces a phenomenally large output (incidentally generating a few INTERNAL SELF-TEST ERRORs along the way). Even a conversion from a single loop iteration produces lots of output:

(* WARNING: this expression produces SELF-TEST errors *)
TracePrint @
  UnitConvert[
    x1 B[Quantity[0.1, "Micrometers"], Quantity[1000, "Kelvins"]]
  , "Watts"/("Meters")^3
  ]

A cursory examination of the TracePrint output suggests that all unit operations are being performed symbolically from first principles each time. In particular, I did not see any sign of "the usual" unit calculation optimizations (like representing quantities as scaled vectors using a standard unit basis). More importantly, the implementation appears to operate upon large lists of all defined units instead of operating only upon the referenced units. For example, a large Dispatch table containing all defined units is built repeatedly.

It would seem that there is plenty of room for optimization :)

$\endgroup$
6
  • 10
    $\begingroup$ "If you optimize on V1, you're spoiling the next contract" (one of my first managers, eons ago) $\endgroup$ Commented Dec 17, 2012 at 4:26
  • $\begingroup$ The "transparency for MMA beginners" will probably mean the punctuation will end up where you don't like it. +1 $\endgroup$
    – Rojo
    Commented Dec 17, 2012 at 5:16
  • $\begingroup$ @Rojo ... but beginners have no punctuation prejudice to unlearn :D Seriously though, I don't think my answer is truly suitable for use by beginners. The Module / With / symbolic evaluation combination is pretty subtle and a lot to ask of someone unfamiliar with Mathematica. I believe that this kind of problem argues for a "definition-time evaluation" facility in Mathematica: an attribute on a symbol that triggers evaluation within the body of a SetDelayed expression. In this case, it could be used to predigest a unit specification into a more efficient representation. $\endgroup$
    – WReach
    Commented Dec 17, 2012 at 15:57
  • $\begingroup$ Considerable improvement running my code at top of this thread in v10.0 ... 9.3 seconds on the same laptop hardware, down from 27 seconds in v9 and 12 seconds (for a functionally similar code, see my comment below) in v8. $\endgroup$ Commented Jul 14, 2014 at 22:44
  • $\begingroup$ Units are still very slow in v11. Is there a fundamentally algorithmic impediment to make this fast that we are missing? Otherwise it seems that optimizing this should be straightforward. $\endgroup$
    – a06e
    Commented Apr 18, 2017 at 0:03
17
$\begingroup$

Just some analysis to try to find where the slow down. On my PC, it took 25 seconds to build the table. ps. I never used Units before.

Your main loop:

x = UnitConvert@Quantity["PlanckConstant" "SpeedOfLight"/"BoltzmannConstant"]
x1 = UnitConvert@Quantity[2, "PlanckConstant" ("SpeedOfLight")^2]
B[L_, T_] := (L^(-5))/(Exp[x/(L T)] - 1)
c = Quantity[1000, "Kelvins"]

Table[{Quantity[ll, "Micrometers"],
     UnitConvert[x1 B[Quantity[ll, "Micrometers"], c],"Watts"/("Meters")^3]}, 
     {ll, 0.1, 15, 0.2}
]; // Timing

 (*   {25.693365, Null}  *)

CPU Analysis

Raising to -5 powers is fast

 Table[Quantity[ll, "Micrometers"]^(-5) , {ll, 0.1, 15,0.2}]; // Timing
 (*   {0.764405, Null}  *)

The Exp use is also fast

   Table[Exp[Quantity[ll, "Micrometers"]] , {ll, 0.1, 15, 0.2}]; // Timing
   (*  {0.015600, Null}  *)

Products of Quantities is little slow, but it depends on which Quantities:

 c = Quantity[1000, "Kelvins"]
 Table[Quantity[ll, "Micrometers"] c , {ll, 0.1, 15,0.2}]; // Timing
 (*  {2.886019, Null}  (1) *)

Lets go back and add the division now, to try to build the B function as it was. Dividing quantities now:

  x = UnitConvert@Quantity["PlanckConstant" "SpeedOfLight"/"BoltzmannConstant"]
  Table[x/Quantity[ll, "Micrometers"], {ll, 0.1, 15, 0.2}]; // Timing
  (*  {7.347647, Null}   (2) *)

Lets add the c back, now will do products of Quantities, notice CPU doubles *)

  Table[x/(Quantity[ll, "Micrometers"] c), {ll, 0.1, 15, 0.2}]; // Timing
  (*  {14.040090, Null}   (3) *)

Add the Exp back. CPU does not change much as expected

  Table[Exp[x/(Quantity[ll, "Micrometers"] c)], {ll, 0.1, 15, 0.2}]; // Timing
  (*  {14.008890, Null}  *)

Add the second division of Quantity you had. Notice no change in CPU since Units are doing division by a number now (no units in result of Exp[...], just a number)

   Table[ Quantity[ll, "Micrometers"]^(-5)/Exp[x/(Quantity[ll, "Micrometers"] c)], 
   {ll, 0.1, 15,0.2}]; // Timing
   (* {14.882495, Null}     *)

Now add x1 back. This is product of Quantity again, Notice the CPU doubles now

   x1 = UnitConvert@Quantity[2, "PlanckConstant" ("SpeedOfLight")^2]
   Table[ x1 (Quantity[ll,"Micrometers"]^(-5)/Exp[x/(Quantity[ll,"Micrometers"]c)]), 
   {ll, 0.1, 15,0.2}]; // Timing
   (*  {24.679358, Null}   (4)  *)

here we go.

Summary

  1. Product of "Kelvins" with "Micrometers" about 3 second. step (1)
  2. division of meter kelvin by "Micrometers" is about 4 seconds. step (2)
  3. division of "Micrometers" by micro meter kelvin about 7 seconds step (3)
  4. product of kg meter to the fourth per second cube with micrometer to the fifth 14 seconds. step (4).

Tried to Compile these, but can't on Units.

Conclusion

product and division between units seems to be slow (compared to just numbers!), it also seems to depend on which units are involved.

More thorough analysis is needed to find why that is and optimize these operations.

$\endgroup$
5
  • 13
    $\begingroup$ My conclusion: just set $\hbar = c = 1$. $\endgroup$
    – Jens
    Commented Dec 17, 2012 at 3:37
  • 5
    $\begingroup$ @Jens you forgot to include $\pi = e = i = 1$ just for convenience. :) $\endgroup$
    – rcollyer
    Commented Dec 17, 2012 at 4:46
  • $\begingroup$ @rcollyer That comment would work really well with this question... $\endgroup$
    – Jens
    Commented Dec 17, 2012 at 5:12
  • 2
    $\begingroup$ @rcollyer And 2=1 also helps sometimes. $\endgroup$ Commented Dec 17, 2012 at 6:31
  • $\begingroup$ @Jens you're right. It seems I forgot to include a few things $\mathbb{R} = 1$ should do it. $\endgroup$
    – rcollyer
    Commented Dec 17, 2012 at 6:32
7
$\begingroup$

Mathematica 10 got out just a couple of days ago and I have compared speed between v 9.0.0 and v 10.0.0 on the same computer (64-bit Win7) needed to generate p table:

Version 9: 17.5 s

Version 10: 5 s

Over 3 times faster, but still slow. Anyway - if you use units move to v 10 seems like a good idea.

$\endgroup$
4
  • $\begingroup$ On my 64bit Windows 7 PC, using a slightly different calculation, Mathematica 10 is 30% slower with units than Mathematica 9. I moved back to Mathematica 9! $\endgroup$ Commented Sep 23, 2014 at 11:57
  • $\begingroup$ On my computer 10.2 still needs 4 seconds while 9.0 needs 9 seconds. Notice that in a fresh kernel AbsoluteTiming@Quantity[1, "Years"] needs 0.8 seconds in 10.2 compared to 0.4 seconds in 9.0. This is all way too slow. $\endgroup$ Commented Jul 20, 2015 at 18:37
  • $\begingroup$ @RolfMertig I get 7 seconds with 9, 5 seconds with 10.0 and 1.7 seconds with 10.4. There are clear improvements, but to make this actually usable for real work, it would need to improve by a couple of orders of magnitude. The problem is that new functions are increasingly adopting units, and thus sooner or later we will be forced to use them and put up with the huge slowdowns. $\endgroup$
    – Szabolcs
    Commented Jul 6, 2016 at 14:18
  • $\begingroup$ @Szabolcs Certainly somebody at Wolfram must care about performance. If not, Mathematica will loose even more to Python and alternatives. I have no time now to compare to timings in Pint. Maybe in the long run we will all use Python from Mathematica and translate finally to Python. $\endgroup$ Commented Jul 6, 2016 at 20:33
5
$\begingroup$

To explore my question a bit I tried the Mathematica 8 AutomaticUnits add-on from http://blog.wolfram.com/2010/12/09/automatic-physical-units-in-mathematica/ You must modify built-in PhysicalConstants.m (benign) as per comments in that blog to play properly with AutomaticUnits. Once done, you can write clearly in Mathematica 8

<<PhysicalConstants`
<<AutomaticUnits`
B[ll_, T_] := (2 PlanckConstant (SpeedOfLight^2) (ll^-5) )/(
  Exp[(PlanckConstant SpeedOfLight)/(ll BoltzmannConstant T)] - 1)
Plot[B[ll Meter, 1000 Kelvin], {ll, 10^-7, 10^-5}, PlotRange -> All, 
  AxesLabel -> {"Meter", "W/\!\(\*SuperscriptBox[\(m\), \(3\)]\)"}]

Which produces

a Planck plot http://www.physics.unc.edu/~cecil/ftp/mathplot.png

in 12 seconds on my 2.4 GHz Core-Duo running Mathematica 8.0.4 under Fedora 17 64-bit, half the time of Mathematica 9.0 units. However, you can't convert units in axes labels with this package, they must be converted manually, i.e. AxesLabel->Automatic gives $J/s/m^3$.

Nonetheless, I find the syntax of this older package to be cleaner for newcomers such as my students than the cumbersome Quantity[...] 9.0 syntax of the new units. And one can use Plot instead of ListLinePlot to benefit from automatic function sampling without need to tabulate a plotting list. Perhaps that's what accelerates?

Luckily AutomaticUnits still works in Mathematica 9 if you globally edit AutomaticUnits.m after installation to replace all occurrences of CompatibleUnitQ to DimensionCompatibleUnitQ as described in the Wolfram blog entry.

$\endgroup$
3
  • 3
    $\begingroup$ Anything over a few milliseconds for such a plot is just wrong. But the comparison to the earlier package is probably useful. Using that as a yardstick, version 9 is only three times slower... But in reality, that's not a yardstick that inspires a lot of confidence. Just one detail: how did you do the timing? Using AbsoluteTiming? $\endgroup$
    – Jens
    Commented Dec 21, 2012 at 1:06
  • $\begingroup$ Timing says MMA 8 / AutomaticUnits runs in 12 seconds vs 27 for MMA 9 units on my machine so 2.3x speedup not 3x. $\endgroup$ Commented Dec 21, 2012 at 13:03
  • 1
    $\begingroup$ Please use *italics* to style and not $LaTeX$ as the latter requires significant rendering overhead, and moreover is a misapplication of the feature. $\endgroup$
    – Mr.Wizard
    Commented Dec 21, 2012 at 23:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.