Skip to content

Commit cb06375

Browse files
committed
major redesign of floating point comparison
1 parent a08bb5a commit cb06375

File tree

5 files changed

+997
-243
lines changed

5 files changed

+997
-243
lines changed

math.hpp

Lines changed: 279 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -20,16 +20,36 @@ namespace unlib {
2020

2121
namespace detail {
2222

23-
template<typename U, typename S, typename V, typename T, std::intmax_t N, std::intmax_t D>
24-
auto pow(const quantity<U,S,T,V>& q, const std::ratio<N,D>) {
25-
return pow_quantity_t<quantity<U,S,T,V>, std::ratio<N,D>>{ std::pow( static_cast<double>(q.get()), static_cast<double>(N)/D ) };
23+
template<typename Float1, typename T>
24+
using if_float1_pt_t = typename std::enable_if< is_floating_point<Float1>::value, T >::type;
25+
template<typename Float1, typename Float2, typename T>
26+
using if_float2_pt_t = typename std::enable_if< is_floating_point<Float1>::value
27+
and is_floating_point<Float2>::value, T >::type;
28+
29+
template<typename D, typename S, typename V, typename T, std::intmax_t Nom, std::intmax_t Den>
30+
auto pow(const quantity<D,S,T,V>& q, const std::ratio<Nom,Den>) {
31+
return pow_quantity_t<quantity<D,S,T,V>, std::ratio<Nom,Den>>{ std::pow( static_cast<double>(q.get()), static_cast<double>(Nom)/Den ) };
2632
}
2733

2834
template<typename U, typename S, typename V, typename T>
2935
auto pow(const quantity<U,S,T,V>& q, const std::ratio<1,2>) {
3036
return pow_quantity_t<quantity<U,S,T,V>,std::ratio<1,2>>(std::sqrt(q.get()));
3137
}
3238

39+
template<typename TolTag, typename F, typename>
40+
struct tolerance_aux {
41+
using tolerance_tag_type = TolTag;
42+
using value_type = F;
43+
value_type v;
44+
};
45+
46+
struct tolerance_val_tag;
47+
struct tolerance_nom_tag;
48+
struct tolerance_frc_tag;
49+
50+
template<typename V, typename S=scale_t<1>>
51+
using fraction = quantity<dimensionless, S, V, no_tag>;
52+
3353
}
3454

3555
/**
@@ -113,6 +133,262 @@ template<typename U, typename S, typename V, typename T>
113133
auto cbrt(quantity<U,S,T,V> q) {return pow<std::ratio<1,3>>(q);}
114134
/** @} */
115135

136+
template<typename D, typename F, typename S1, typename S2, typename T>
137+
auto min(const unlib::quantity<D,S1,F,T> lhs, const unlib::quantity<D,S2,F,T> rhs) {
138+
using std::min;
139+
return unlib::quantity<D,S1,F,T>{min(lhs.get(), rhs.template get_scaled<S1>())};
140+
}
141+
142+
template<typename D, typename F, typename S1, typename S2, typename T>
143+
auto max(const unlib::quantity<D,S1,F,T> lhs, const unlib::quantity<D,S2,F,T> rhs) {
144+
using std::max;
145+
return unlib::quantity<D,S1,F,T>{max(lhs.get(), rhs.template get_scaled<S1>())};
146+
}
147+
148+
149+
/**
150+
* @{
151+
*
152+
* @brief Floating point comparisons
153+
*
154+
* Due to the limitations of floating point precision, more often then
155+
* never values which, mathematically, ought to be equal, are in fact not
156+
* totally equal. Therefore, floating point comparisons should be done with
157+
* a certain tolerance.
158+
*
159+
* This provides a framework for comparing floating point types and
160+
* quantities of floating point types. The tolerance can either be an
161+
* absolute value, or a fraction of a nominal value. For example, in order to
162+
* test whether two masses are equal, the tolerance could be given either as
163+
* an absolute value (10mg) or as a fraction (0.1%) of a nominal weight (1t).
164+
* The library provides default values for the nominal and the weight, which
165+
* can be used, but this must be indicated explicitly.
166+
*
167+
* The defaults are accessed via traits which can be specialized by users in
168+
* order to provide their own defaults.
169+
*/
170+
171+
constexpr double tolerance_default_nominal = 1.; /**< default tolerance nominal for is_near() etc. */
172+
constexpr double tolerance_default_fraction = 0.0001; /**< default tolerance fraction for is_near() etc. */
173+
174+
/**
175+
* @{
176+
*
177+
* @brief wrapper for tolerance values
178+
*
179+
* In order to distinguish between tolerance values, nominals, and fractions
180+
* in floating point comparisons, those are wrapped in these types.
181+
*
182+
* @note These are the result types of @ref tolerance_value(), @ref tolerance_nominal(),
183+
* and @ref tolerance_fraction(). Users should not have to deal with them
184+
* directly.
185+
*/
186+
template<typename V> using tolerance_val = detail::tolerance_aux<detail::tolerance_val_tag,V,V>;
187+
template<typename V> using tolerance_nom = detail::tolerance_aux<detail::tolerance_nom_tag,V,V>;
188+
template<typename V, typename Q> using tolerance_frc = detail::tolerance_aux<detail::tolerance_frc_tag,V,Q>;
189+
/** @} */
190+
191+
/**
192+
* @{
193+
*
194+
* @brief tolerance value calculations
195+
*
196+
* A tolerance value can be calculated by multiplying a nominal and a fraction.
197+
*
198+
* @param n tolerance nominal value
199+
* @param f tolerance fraction value
200+
*
201+
* @return tolerance value
202+
*/
203+
template<typename U, typename S, typename SF, typename V, typename T, typename Q>
204+
inline constexpr
205+
tolerance_val<quantity<U,S,V,T>> operator*(tolerance_nom<quantity<U,S,V,T>> n, tolerance_frc<detail::fraction<V,SF>,Q> f)
206+
{return tolerance_val<quantity<U,S,V,T>>{n.v * f.v.template get_scaled<no_scaling>()};}
207+
208+
template<typename F, typename SF, typename Q>
209+
inline constexpr
210+
tolerance_val<F> operator*(tolerance_nom<F> n, tolerance_frc<detail::fraction<F,SF>,Q> f)
211+
{return tolerance_val<F>{n.v * f.v.template get_scaled<no_scaling>()};}
212+
/** @} */
213+
214+
/**
215+
* @{
216+
* @brief tolerance traits
217+
*
218+
* These provide the necessary functions to obtain default tolerance values,
219+
* nominals, and fractions.
220+
*
221+
* @tparam F Floating point type.
222+
*/
223+
template<typename F>
224+
struct tolerance_traits {
225+
using float_t = detail::if_float1_pt_t<F,F>;
226+
using fract_t = detail::fraction<float_t>;
227+
228+
static constexpr auto value () {return nominal() * fraction();}
229+
static constexpr auto nominal () {return tolerance_nom<float_t >{ static_cast<float_t>(tolerance_default_nominal ) };}
230+
static constexpr auto fraction() {return tolerance_frc<fract_t,float_t>{fract_t{static_cast<float_t>(tolerance_default_fraction)}};}
231+
};
232+
233+
/** specialization for quantities of floating point types */
234+
template<typename U, typename S, typename F, typename T>
235+
struct tolerance_traits<quantity<U,S,F,T>> {
236+
using float_t = detail::if_float1_pt_t<F,F>;
237+
using quant_t = quantity<U,S,float_t,T>;
238+
using fract_t = detail::fraction<float_t>;
239+
240+
static constexpr auto value () {return nominal() * fraction();}
241+
static constexpr auto nominal () {return tolerance_nom<quant_t >{quant_t{tolerance_traits<float_t>::nominal ().v}};}
242+
static constexpr auto fraction() {return tolerance_frc<fract_t,quant_t>{ tolerance_traits<float_t>::fraction().v };}
243+
};
244+
/** @} */
245+
246+
/**
247+
* @{
248+
*
249+
* @brief create a tolerance value, nominal, or fraction
250+
*
251+
* These function turn a value into a tolerance value, nominal, or fraction.
252+
*
253+
* @param val value to turn into a tolerance value, nominal, or fraction
254+
*
255+
* @return tolerance value, nominal, or fraction
256+
*/
257+
template<typename T> constexpr auto tolerance_value (T val) {return tolerance_val<T >{val};}
258+
template<typename T> constexpr auto tolerance_nominal (T nom) {return tolerance_nom<T >{nom};}
259+
template<typename T, typename F> constexpr auto tolerance_fraction(F frc) {return tolerance_frc<F,T>{frc};}
260+
/** @} */
261+
262+
/**
263+
* @{
264+
*
265+
* @brief get the default tolerance value, nominal, or fraction
266+
*
267+
* These function return the tolerance value, nominal, or fraction for a
268+
* given type.
269+
*
270+
* @return default tolerance value, nominal, or fraction
271+
*/
272+
template<typename T> constexpr auto tolerance_value () {return tolerance_traits<T>::value ();}
273+
template<typename T> constexpr auto tolerance_nominal () {return tolerance_traits<T>::nominal ();}
274+
template<typename T> constexpr auto tolerance_fraction() {return tolerance_traits<T>::fraction();}
275+
/** @} */
276+
277+
namespace detail {
278+
279+
template<typename V, typename F, typename Q> constexpr auto get_tol_val(tolerance_nom<V> n, tolerance_frc<F,Q> f) {return n * f;}
280+
template<typename V, typename F, typename Q> constexpr auto get_tol_val(tolerance_frc<F,Q> f, tolerance_nom<V> n) {return get_tol_val(n,f);}
281+
template<typename V > constexpr auto get_tol_val(tolerance_aux<tolerance_val_tag,V,V> v) {return v;}
282+
template<typename V > constexpr auto get_tol_val(tolerance_aux<tolerance_nom_tag,V,V> n) {return get_tol_val(n, tolerance_traits<V>::fraction());}
283+
template<typename F, typename Q> constexpr auto get_tol_val(tolerance_aux<tolerance_frc_tag,F,Q> f) {return get_tol_val(f, tolerance_traits<Q>::nominal ());}
284+
285+
template<typename T1, typename T2, typename Tol> constexpr bool is_near_impl (T1 lval, T2 rval, tolerance_val<Tol> tol) {using std::abs; return abs(lval-rval) <= abs(tol.v);}
286+
template<typename T1, typename T2, typename Tol> constexpr bool is_smaller_impl(T1 lval, T2 rval, tolerance_val<Tol> tol) {using std::abs; return not is_near(lval,rval,tol) and lval-rval < abs(tol.v);}
287+
template<typename T1, typename T2, typename Tol> constexpr bool is_greater_impl(T1 lval, T2 rval, tolerance_val<Tol> tol) {using std::abs; return not is_near(lval,rval,tol) and lval-rval > abs(tol.v);}
288+
289+
template<typename T1, typename T2, typename Aux> constexpr detail::if_float2_pt_t<T1,T2,bool> is_near (T1 lval, T2 rval, Aux tol) {return is_near_impl (lval, rval, get_tol_val(tol));}
290+
template<typename T1, typename T2, typename Aux> constexpr detail::if_float2_pt_t<T1,T2,bool> is_smaller(T1 lval, T2 rval, Aux tol) {return is_smaller_impl(lval, rval, get_tol_val(tol));}
291+
template<typename T1, typename T2, typename Aux> constexpr detail::if_float2_pt_t<T1,T2,bool> is_greater(T1 lval, T2 rval, Aux tol) {return is_greater_impl(lval, rval, get_tol_val(tol));}
292+
template<typename T1, typename T2, typename Aux1, typename Aux2> constexpr detail::if_float2_pt_t<T1,T2,bool> is_near (T1 lval, T2 rval, Aux1 tol1, Aux2 tol2) {return is_near_impl (lval, rval, get_tol_val(tol1,tol2));}
293+
template<typename T1, typename T2, typename Aux1, typename Aux2> constexpr detail::if_float2_pt_t<T1,T2,bool> is_smaller(T1 lval, T2 rval, Aux1 tol1, Aux2 tol2) {return is_smaller_impl(lval, rval, get_tol_val(tol1,tol2));}
294+
template<typename T1, typename T2, typename Aux1, typename Aux2> constexpr detail::if_float2_pt_t<T1,T2,bool> is_greater(T1 lval, T2 rval, Aux1 tol1, Aux2 tol2) {return is_greater_impl(lval, rval, get_tol_val(tol1,tol2));}
295+
296+
}
297+
298+
/** @{
299+
* @brief compare two floating point values
300+
*
301+
* This compares two floating point values, or quantities of floating point
302+
* values.
303+
*
304+
* @param lval left value to compare
305+
* @param rval right value to compare
306+
* @param tol, tol1, tol2 tolerance
307+
*
308+
* @note As tolerance, either an absolute tolerance value (created via
309+
* @ref tolerance_value()) or a tolerance nominal (created via @ref
310+
* tolerance_nominal()) and a tolerance fraction (created via @ref
311+
* tolerance_fraction()) can be passed. Either the nominal or the
312+
* fraction can also be omitted, in which case the default is obtained
313+
* via the @ref tolerance_traits.
314+
*
315+
* @note The quantities do not need to be of the same scale, but they need to
316+
* have the same dimensions, value types, and tags.
317+
*
318+
* @return true if both quantities are equal/greater/smaller within the
319+
* tolerance given
320+
*/
321+
template<typename F1, typename F2, typename TT, typename TF, typename X>
322+
bool is_near(F1 lval,F2 rval, detail::tolerance_aux<TT,TF,X> tol)
323+
{return detail::is_near(lval, rval, tol);}
324+
325+
template<typename F1, typename F2, typename TT, typename TF, typename X>
326+
bool is_greater(F1 lval,F2 rval, detail::tolerance_aux<TT,TF,X> tol)
327+
{return detail::is_greater(lval, rval, tol);}
328+
329+
template<typename F1, typename F2, typename TT, typename TF, typename X>
330+
bool is_smaller(F1 lval,F2 rval, detail::tolerance_aux<TT,TF,X> tol)
331+
{return detail::is_smaller(lval, rval, tol);}
332+
333+
template<typename F1, typename F2, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
334+
bool is_near(F1 lval,F2 rval, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
335+
{return detail::is_near(lval, rval, tol1, tol2);}
336+
337+
template<typename F1, typename F2, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
338+
bool is_smaller(F1 lval,F2 rval, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
339+
{return detail::is_smaller(lval, rval, tol1, tol2);}
340+
341+
template<typename F1, typename F2, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
342+
bool is_greater(F1 lval,F2 rval, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
343+
{return detail::is_greater(lval, rval, tol1, tol2);}
344+
/** @} */
345+
346+
/** @{
347+
* @brief compare a floating point value with zero
348+
*
349+
* This compares a floating point value, or a quantity of a floating point
350+
* value, with zero.
351+
*
352+
* @param val value to compare
353+
* @param tol, tol1, tol2 tolerance
354+
*
355+
* @note As tolerance, either an absolute tolerance value (created via
356+
* @ref tolerance_value()) or a tolerance nominal (created via @ref
357+
* tolerance_nominal()) and a tolerance fraction (created via @ref
358+
* tolerance_fraction()) can be passed. Either the nominal or the
359+
* fraction can also be omitted, in which case the default is obtained
360+
* via the @ref tolerance_traits.
361+
*
362+
* @return true if the quantity is equal/greater/smaller zero within the
363+
* tolerance given
364+
*/
365+
template<typename F, typename TT, typename TF, typename X>
366+
bool is_near_zero(F val, detail::tolerance_aux<TT,TF,X> tol)
367+
{return detail::is_near(val,F{}, tol);}
368+
369+
template<typename F, typename TT, typename TF, typename X>
370+
bool is_smaller_zero(F val, detail::tolerance_aux<TT,TF,X> tol)
371+
{return detail::is_smaller(val,F{}, tol);}
372+
373+
template<typename F, typename TT, typename TF, typename X>
374+
bool is_greater_zero(F val, detail::tolerance_aux<TT,TF,X> tol)
375+
{return detail::is_greater(val,F{}, tol);}
376+
377+
template<typename F, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
378+
bool is_near_zero(F val, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
379+
{return detail::is_near(val,F{}, tol1, tol2);}
380+
381+
template<typename F, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
382+
bool is_smaller_zero(F val, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
383+
{return detail::is_smaller(val,F{}, tol1, tol2);}
384+
385+
template<typename F, typename TT1, typename TF1, typename TT2, typename TF2, typename X1, typename X2>
386+
bool is_greater_zero(F val, detail::tolerance_aux<TT1,TF1,X1> tol1, detail::tolerance_aux<TT2,TF2,X2> tol2)
387+
{return detail::is_greater(val,F{}, tol1, tol2);}
388+
/** @} */
389+
390+
/** @} */
391+
116392
}
117393

118394
#endif //UNLIB_MATH_HPP

quantity.hpp

Lines changed: 14 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -48,10 +48,22 @@ class quantity;
4848
template<typename T> struct is_quantity : std::false_type {};
4949
template<typename U, typename S, typename V, typename T> struct is_quantity<quantity<U,S,V,T>> : std::true_type {};
5050

51-
template<typename T> constexpr bool is_quantity_v = is_quantity<T>::value;
52-
5351
namespace detail {
5452

53+
template<typename TFloat>
54+
struct is_floating_point : std::is_floating_point<TFloat> {};
55+
template<typename TU, typename TS, typename TV, typename TT>
56+
struct is_floating_point<quantity<TU,TS,TV,TT>> : is_floating_point<TV> {};
57+
58+
template<typename TFloat>
59+
struct is_integral: std::is_integral<TFloat> {};
60+
template<typename TU, typename TS, typename TV, typename TT>
61+
struct is_integral<quantity<TU,TS,TV,TT>> : is_integral<TV> {};
62+
63+
template<typename T> constexpr bool is_quantity_v = is_quantity<T>::value;
64+
template<typename T> constexpr bool is_floating_point_v = is_floating_point<T>::value;
65+
template<typename T> constexpr bool is_integral_v = is_integral<T>::value;
66+
5567
template<typename NewScale, typename OldScale>
5668
struct value_rescaler {
5769
using conversion_scale = std::ratio_divide<OldScale, NewScale>;
@@ -73,16 +85,6 @@ struct value_rescaler<Scale,Scale> {
7385
template<typename NewScale, typename OldScale, typename ValueType>
7486
constexpr ValueType rescale_value(ValueType v) {return value_rescaler<NewScale,OldScale>::rescale_value(v);}
7587

76-
template<typename T> constexpr bool is_floating_point_v = std::is_floating_point<T>::value;
77-
template<typename T> constexpr bool is_integral_v = std::is_integral<T>::value;
78-
79-
template<typename T>
80-
constexpr std::enable_if_t<is_integral_v<T>,T> get_tolerance() {return T{};}
81-
template<typename F>
82-
constexpr std::enable_if_t<is_floating_point_v<F>,F> get_tolerance() {return 100*std::numeric_limits<float>::epsilon();}
83-
template<typename Q>
84-
constexpr std::enable_if_t<is_quantity_v<Q>, Q> get_tolerance() {return Q{get_tolerance<typename Q::value_type>()};}
85-
8688
}
8789

8890
/**
@@ -379,14 +381,6 @@ class quantity< unit< std::ratio< TimeNum, TimeDen>
379381
*
380382
* @note There is a free function of the same name.
381383
*/
382-
template<typename U, typename S, typename V, typename T, typename D = quantity>
383-
constexpr bool is_near(const quantity<U,S,V,T>& rhs, const D tolerance = detail::get_tolerance<D>()) const
384-
{
385-
static_assert(are_units_compatible_v<unit_type , U>, "fundamentally incompatible units" );
386-
static_assert(detail::is_same_v <value_type, V>, "different value types (use value_cast)");
387-
static_assert(detail::is_same_v <tag_type , T>, "incompatible tags (use tag_cast)" );
388-
return (*this-rhs).is_near_zero(tolerance);
389-
}
390384

391385
/**
392386
* @brief Check whether a floating point quantity is near zero
@@ -400,16 +394,8 @@ class quantity< unit< std::ratio< TimeNum, TimeDen>
400394
*
401395
* @note There is a free function of the same name.
402396
*/
403-
template<typename D = quantity>
404-
constexpr bool is_near_zero(const D tolerance = detail::get_tolerance<D>()) const
405-
{return is_near_zero_impl(tolerance);}
406397

407398
private:
408-
template<typename D>
409-
constexpr bool is_near_zero_impl(const D tolerance) const {return std::abs(value) <= tolerance;}
410-
template<typename UD, typename SD, typename VD, typename TD>
411-
constexpr bool is_near_zero_impl(const quantity<UD,SD,VD,TD>& tolerance) const
412-
{return quantity{std::abs(value)} <= tolerance;}
413399

414400
value_type value;
415401
};
@@ -424,15 +410,6 @@ template<typename NewScale, typename U, typename S, typename V, typename T>
424410
constexpr typename quantity<U,S,V,T>::value_type get_scaled(const quantity<U,S,V,T>& q, NewScale = NewScale{})
425411
{return q.template get_scaled<NewScale>();}
426412

427-
template< typename U1, typename S1, typename V1, typename T1
428-
, typename U2, typename S2, typename V2, typename T2
429-
, typename D = quantity<U1,S1,V1,T1> >
430-
constexpr bool is_near(quantity<U1,S1,V1,T1> lhs, const quantity<U2,S2,V2,T2>& rhs, const D tolerance = detail::get_tolerance<D>())
431-
{return lhs.is_near(rhs,tolerance);}
432-
433-
template< typename U, typename S, typename V, typename T, typename D = quantity<U,S,V,T>>
434-
constexpr bool is_near_zero(const quantity<U,S,V,T>& q, const D tolerance = detail::get_tolerance<D>())
435-
{return q.is_near_zero(tolerance);}
436413

437414
/**
438415
* @{

0 commit comments

Comments
 (0)