Upgrade to Pro
— share decks privately, control downloads, hide ads and more …
Speaker Deck
Speaker Deck
PRO
Sign in
Sign up for free
Speed, Correctness, or Simplicity: Choose 3
Tom Switzer
January 30, 2015
Programming
1
190
Speed, Correctness, or Simplicity: Choose 3
This talk introduces the floating point filter implementation in Spire (spire.math.FpFilter).
Tom Switzer
January 30, 2015
Tweet
Share
Other Decks in Programming
See All in Programming
"What's new in Swift"の要約 / swift_5_7_summary
uhooi
1
100
ドメインモデル方式のクラス設計 座談会
masuda220
PRO
3
1k
A Philosophy of Software Design 後半
yosuke_furukawa
PRO
8
2.4k
Migrating to Kotlin State & Shared Flows
heyitsmohit
1
180
dotdotdot_in_predict_function
bk_18
1
170
engineer
spacemarket
0
460
個人開発でReact Native + Expo製アプリを作った話
ryonakae
1
430
Amazon ECSのネットワーク関連コストの話
msato
0
610
はじめてのプルリク - BLEA 編
watany
0
140
One does not simply: migrating to Android 12 🤯
oleur
1
120
Jakarta EE 10 - Feature by Feature with My Open Source Journey
ivargrimstad
0
1.8k
Angular‘s Future without NgModules: Architectures with Standalone Components @enterJS
manfredsteyer
PRO
0
170
Featured
See All Featured
What's new in Ruby 2.0
geeforr
336
30k
Agile that works and the tools we love
rasmusluckow
319
19k
Bash Introduction
62gerente
597
210k
I Don’t Have Time: Getting Over the Fear to Launch Your Podcast
jcasabona
12
920
Thoughts on Productivity
jonyablonski
43
2.3k
VelocityConf: Rendering Performance Case Studies
addyosmani
316
22k
The Mythical Team-Month
searls
209
39k
CoffeeScript is Beautiful & I Never Want to Write Plain JavaScript Again
sstephenson
151
13k
Designing for humans not robots
tammielis
241
23k
Java REST API Framework Comparison - PWX 2021
mraible
PRO
11
4.7k
The Web Native Designer (August 2011)
paulrobertlloyd
74
1.9k
Designing on Purpose - Digital PM Summit 2013
jponch
106
5.6k
Transcript
Speed, Correctness, or Simplicity: Choose 3 Tom
Switzer @9xxit h;ps://github.com/9xxit/fpfilter-‐talk
Overview Floa9ng point is “good enough”…
most of the 9me.
Op9ons Use Double, live with the errors.
Use higher precision type, live with performance loss. But, there is a 3rd op9on…
Floa9ng Point Filters Use floa9ng point when you
can. Use higher precision when you can’t.
Err… Not So Simple Solve problem using floa9ng point
approxima9on… Maintain an error bound on approxima9on. Re-‐evaluate with exact type if error too large.
The Catch
What is the determinant of my matrix?
Not Good For: Minimizing Errors in Floa9ng Point Arithme9c
What is the sign of the determinant of
my matrix?
Good For: Making a Decision
FpFilter[A] Simple wrapper: FpFilter[Rational] Standard Opera2ons +, -‐,
*, /, .sqrt, etc Fast predictes signum, compare, isWhole, etc.
FpFilter[A] class FpFilter[A]( apx: Double, mes: Double, ind: Int,
exact: => A ) { … } floa9ng point approxima9on error bounds
FpFilter[A] class FpFilter[A]( apx: Double, mes: Double, ind: Int,
exact: => A ) { … } error bounds “Exact Geometric Computa2on Using Cascading.” Burnikel, Funke & Seel.
FpFilter[A] class FpFilter[A]( apx: Double, mes: Double, ind: Int,
exact: => A ) { … } error bounds thunk for higher precision Welcome to …
… Macro City def abs(implicit ev: Signed[A]): FpFilter[A] =
macro FpFilter.absImpl[A] def unary_- (implicit ev: Rng[A]) : FpFilter[A] = macro FpFilter.negateImpl[A] def +(rhs: FpFilter[A])(implicit ev: Semiring[A]): FpFilter[A] = macro FpFilter.plusImpl[A] def -(rhs: FpFilter[A])(implicit ev: Rng[A]): FpFilter[A] = macro FpFilter.minusImpl[A] def *(rhs: FpFilter[A])(implicit ev: Semiring[A]): FpFilter[A] = macro FpFilter.timesImpl[A] def /(rhs: FpFilter[A])(implicit ev: Field[A]): FpFilter[A] = macro FpFilter.divideImpl[A] def sqrt(implicit ev: NRoot[A]): FpFilter[A] = macro FpFilter.sqrtImpl[A] def <(rhs: FpFilter[A])(implicit ev0: Signed[A], ev1: Rng[A]): Boolean = macro FpFilter.ltImpl[A] def >(rhs: FpFilter[A])(implicit ev0: Signed[A], ev1: Rng[A]): Boolean = macro FpFilter.gtImpl[A] def <=(rhs: FpFilter[A])(implicit ev0: Signed[A], ev1: Rng[A]): Boolean = macro FpFilter.ltEqImpl[A] def >=(rhs: FpFilter[A])(implicit ev0: Signed[A], ev1: Rng[A]): Boolean = macro FpFilter.gtEqImpl[A] def ===(rhs: FpFilter[A])(implicit ev0: Signed[A], ev1: Rng[A]): Boolean = macro FpFilter.eqImpl[A] def signum(implicit ev: Signed[A]): Int = macro FpFilter.signImpl[A]
… Macro City • Operator fusion – No intermediate
alloca9ons • In-‐line approxima9on + error bounds – Fast, Double arithme9c • Thunk becomes inner defs – Compile down to private methods
Turn this… (x + y).signum
… into this. val fpf$tmp$macro$38 = x.value; val fpf$apx$macro$39
= fpf$tmp$macro$38; val fpf$mes$macro$40 = java.lang.Math.abs(fpf$tmp$macro$38); def fpf$exact$macro$42 = spire.algebra.Field.apply[spire.math.Algebraic](Algebraic.AlgebraicAlgebra).fromDouble(fpf $tmp$macro$38); val fpf$tmp$macro$43 = y.value; val fpf$apx$macro$44 = fpf$tmp$macro$43; val fpf$mes$macro$45 = java.lang.Math.abs(fpf$tmp$macro$43); def fpf$exact$macro$47 = spire.algebra.Field.apply[Algebraic](Algebraic.AlgebraicAlgebra).fromDouble(fpf$tmp$macro $43); val fpf$apx$macro$48 = fpf$apx$macro$39.+(fpf$apx$macro$44); val fpf$mes$macro$49 = fpf$mes$macro$40.+(fpf$mes$macro$45); def fpf$exact$macro$51 = Algebraic.AlgebraicAlgebra.plus( fpf$exact$macro$42, fpf$exact$macro$47); val fpf$err$macro$52 = fpf$mes$macro$49.$times(1).$times(2.220446049250313E-16); if (fpf$apx$macro$48 > fpf$err$macro$52 && fpf$apx$macro$48 < Double.POSITIVE_INFINITY) 1 else if (fpf$apx$macro$48 < fpf$err$macro$52.unary_$minus && fpf$apx$macro$48 > Double.NEGATIVE_INFINITY) -1 else if (fpf$err$macro$52 == 0.0) 0 else Algebraic.AlgebraicAlgebra.signum(fpf$exact$macro$51)
Examples
2D Orienta2on
p q r
p q r
p q r RIGHT
p r q
p r q LEFT
p r q
p r q NO TURN
trait Turn[@spec A] { def apply( px: A, py: A,
qx: A, qy: A, rx: A, ry: A ): Int }
object FastTurn extends Turn[Double] { def apply( px: Double, py:
Double, qx: Double, qy: Double, rx: Double, ry: Double ): Int = signum { (qx - px) * (ry - py) - (rx - px) * (qy - py) } }
Accuracy of Fast Turn
object ExactTurn extends Turn[Double] { def apply( px: Double, py:
Double, qx: Double, qy: Double, rx: Double, ry: Double ): Int = { val pxa = Algebraic(px) val pya = Algebraic(py) val qxa = Algebraic(qx) val qya = Algebraic(qy) val rxa = Algebraic(rx) val rya = Algebraic(ry) ((qxa - pxa) * (rya - pya) – (rxa - pxa) * (qya - pya)).signum } }
10,000x Slower!
Let’s try again…
object FilteredTurn extends Turn[Double] { def apply( px: Double, py:
Double, qx: Double, qy: Double, rx: Double, ry: Double ): Int = { val pxf = FpFilter.exact[Algebraic](px) val pyf = FpFilter.exact[Algebraic](py) val qxf = FpFilter.exact[Algebraic](qx) val qyf = FpFilter.exact[Algebraic](qy) val rxf = FpFilter.exact[Algebraic](rx) val ryf = FpFilter.exact[Algebraic](ry) ((qxf - pxf) * (ryf - pyf) – (rxf - pxf) * (qyf - pyf)).signum } }
FilteredTurn Speed Rela9ve to FastTurn
Polynomial Root Finding
Polynomial[A]
Interval[A] Root
“Quadra2c Interval Refinement for Real Roots.” John AbboT.
QIR for short.
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=8)
QIR (N=16)
QIR (N=16)
QIR (N=16)
QIR • Requires 2 polynomial evalua9ons – High precision
generally required • Very fast convergence (quadra9c) – Under some assump9ons • Occasionally fails when assump9ons not met – Fallsback to bisec9on!
QIR (N=8)
QIR (N=8)
QIR (N=8) FAILED!
Falls back to bisec9on…
Bisec9on
Bisec9on
Bisec9on
Bisec9on
Bisec9on
Bisec9on
Bisec9on
Bisec9on Requires only sign tests Converges
slowly, 1 bit at a 9me!
trait SignTest[A] { def apply( poly: Polynomial[A], x: A ):
Sign }
final class FilteredSignTest[@sp A: Semiring]( implicit A: IsAlgebraic[A] ) extends
SignTest[A] { def apply(poly: Polynomial[A], x: A): Sign = { val x0 = FpFilter(A.toDouble(x), A.toAlgebraic(x)) @tailrec def loop(acc: FpFilter[Algebraic], i: Int): Sign = if (i >= 0) { val c = poly.nth(i) val cftr = FpFilter(A.toDouble(c), A.toAlgebraic(c)) loop(cftr + acc * x0, i - 1) } else { Sign(acc.signum) } loop(FpFilter.approx(Algebraic.Zero), poly.degree) } }
Accuracy Using Double
Speed Up from Exact Sign Test Fast (d=8)
Fast (d=16) Fast (d=32) Filtered (d=8) Filtered (d=16) Filtered (d=16)
Summary • Works like any other number type
– Operator fusion + inlining within expressions • Speeds up predicates – Sign tests, comparisons, etc. • Near-‐Double performance – 2-‐4x in most cases h;p://github.com/9xxit/fpfilter-‐talk
Thanks! h;p://github.com/non/spire Tom Switzer @9xxit