1# fuzzy.tcl --
2#
3#    Script to define tolerant floating-point comparisons
4#    (Tcl-only version)
5#
6#    version 0.2: improved and extended, march 2002
7#    version 0.2.1: fix bug #2933130, january 2010
8
9package provide math::fuzzy 0.2.1
10
11namespace eval ::math::fuzzy {
12   variable eps3 2.2e-16
13
14   namespace export teq tne tge tgt tle tlt tfloor tceil tround troundn
15
16# DetermineTolerance
17#    Determine the epsilon value
18#
19# Arguments:
20#    None
21#
22# Result:
23#    None
24#
25# Side effects:
26#    Sets variable eps3
27#
28proc DetermineTolerance { } {
29   variable eps3
30   set eps 1.0
31   while { [expr {1.0+$eps}] != 1.0 } {
32      set eps3 [expr 3.0*$eps]
33      set eps  [expr 0.5*$eps]
34   }
35   #set check [expr {1.0+2.0*$eps}]
36   #puts "Eps3: $eps3 ($eps) ([expr {1.0-$check}] [expr 1.0-$check]"
37}
38
39# Absmax --
40#    Return the absolute maximum of two numbers
41#
42# Arguments:
43#    first      First number
44#    second     Second number
45#
46# Result:
47#    Maximum of the absolute values
48#
49proc Absmax { first second } {
50   return [expr {abs($first) > abs($second)? abs($first) : abs($second)}]
51}
52
53# teq, tne, tge, tgt, tle, tlt --
54#    Compare two floating-point numbers and return the logical result
55#
56# Arguments:
57#    first      First number
58#    second     Second number
59#
60# Result:
61#    1 if the condition holds, 0 if not.
62#
63proc teq { first second } {
64   variable eps3
65   set scale [Absmax $first $second]
66   return [expr {abs($first-$second) <= $eps3 * $scale}]
67}
68
69proc tne { first second } {
70   variable eps3
71
72   return [expr {![teq $first $second]}]
73}
74
75proc tgt { first second } {
76   variable eps3
77   set scale [Absmax $first $second]
78   return [expr {($first-$second) > $eps3 * $scale}]
79}
80
81proc tle { first second } {
82   return [expr {![tgt $first $second]}]
83}
84
85proc tlt { first second } {
86   expr { [tle $first $second] && [tne $first $second] }
87}
88
89proc tge { first second } {
90   if { [tgt $first $second] } {
91      return 1
92   } else {
93      return [teq $first $second]
94   }
95}
96
97# tfloor --
98#    Determine the "floor" of a number and return the result
99#
100# Arguments:
101#    number     Number in question
102#
103# Result:
104#    Largest integer number that is tolerantly smaller than the given
105#    value
106#
107proc tfloor { number } {
108   variable eps3
109
110   set q      [expr {($number < 0.0)? (1.0-$eps3) : 1.0 }]
111   set rmax   [expr {$q / (2.0 - $eps3)}]
112   set eps5   [expr {$eps3/$q}]
113   set vmin1  [expr {$eps5*abs(1.0+floor($number))}]
114   set vmin2  [expr {($rmax < $vmin1)? $rmax : $vmin1}]
115   set vmax   [expr {($eps3 > $vmin2)? $eps3 : $vmin2}]
116   set result [expr {floor($number+$vmax)}]
117   if { $number <= 0.0 || ($result-$number) < $rmax } {
118      return $result
119   } else {
120      return [expr {$result-1.0}]
121   }
122}
123
124# tceil --
125#    Determine the "ceil" of a number and return the result
126#
127# Arguments:
128#    number     Number in question
129#
130# Result:
131#    Smallest integer number that is tolerantly greater than the given
132#    value
133#
134proc tceil { number } {
135   expr {-[tfloor [expr {-$number}]]}
136}
137
138# tround --
139#    Round off a number and return the result
140#
141# Arguments:
142#    number     Number in question
143#
144# Result:
145#    Nearest integer number
146#
147proc tround { number } {
148   tfloor [expr {$number+0.5}]
149}
150
151# troundn --
152#    Round off a number to a given precision and return the result
153#
154# Arguments:
155#    number     Number in question
156#    ndec       Number of decimals to keep
157#
158# Result:
159#    Nearest number with given precision
160#
161proc troundn { number ndec } {
162   set scale   [expr {pow(10.0,$ndec)}]
163   set rounded [tfloor [expr {$number*$scale+0.5}]]
164   expr {$rounded/$scale}
165}
166
167#
168# Determine the tolerance once and for all
169#
170DetermineTolerance
171rename DetermineTolerance {}
172
173} ;# End of namespace
174