-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathmydouble.h
202 lines (193 loc) · 5.63 KB
/
mydouble.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
#ifndef _MY_DOUBLE_H_
#define _MY_DOUBLE_H_
#include <limits>
#include <math.h>
#include <myerror.h>
/* This class behaves to the user like a non-negative double, but
is stored internally as the natural logarithm. Standard mathematical
operations are performed on the logarithm of the number so that it
should not underflow or overflow like a double. */
class mydouble {
protected:
double _log;
bool zero;
public:
/*Default constructor*/
mydouble() {
zero = false;
};
/*Copy constructor*/
mydouble(const double &_doub) {
zero = false;
if(_doub<0.0) myutils::error("mydouble::mydouble(const double&): cannot initialize with negative number");
if(_doub==0.0) setzero();
else _log = log(_doub);
};
/*Copy constructor*/
mydouble(const mydouble &_mydoub) {
zero = _mydoub.zero;
_log = _mydoub._log;
}
/*Conversion operator*/
operator double const() {
return (zero) ? 0.0 : exp(_log);
};
/*Assignment operator*/
mydouble& operator=(const double &_doub) {
zero = false;
if(_doub<0.0) myutils::error("mydouble::operator=(const double&): cannot assign a negative number");
if(_doub==0.0) setzero();
else _log = log(_doub);
return *this;
}
/*Assignment operator*/
mydouble& operator=(const mydouble &_mydoub) {
zero = _mydoub.zero;
_log = _mydoub._log;
return *this;
}
void setlog(const double &log) {
zero = false;
_log = log;
}
mydouble& setzero() {
zero = true;
_log = -std::numeric_limits<double>::max();
return *this;
}
bool iszero() {
return zero;
}
/*** MULTIPLICATION ***/
mydouble operator*(const double &dbl) const {
return operator*(mydouble(dbl));
}
mydouble operator*(const mydouble &mydbl) const {
mydouble a;
if(zero || mydbl.zero) a.setzero();
else a.setlog(_log + mydbl._log);
return a;
}
mydouble& operator*=(const double &dbl) {
if(zero || dbl==0.0) setzero();
else _log += mydouble(dbl)._log;
return *this;
}
mydouble& operator*=(const mydouble &mydbl) {
if(zero || mydbl.zero) setzero();
else _log += mydbl._log;
return *this;
}
/*** DIVISION ***/
mydouble operator/(const double &dbl) const {
return operator/(mydouble(dbl));
}
mydouble operator/(const mydouble &mydbl) const {
mydouble a;
if(zero) a.setzero();
else if(mydbl.zero) error("mydouble::operator/(const mydouble&): division by zero");
else a.setlog(_log - mydbl._log);
return a;
}
mydouble& operator/=(const double &dbl) {
if(dbl==0.0) error("mydouble::operator/=(const double&): division by zero");
else if(!zero) _log -= mydouble(dbl)._log;
return *this;
}
mydouble& operator/=(const mydouble &mydbl) {
if(mydbl.zero) error("mydouble::operator/=(const mydouble&): division by zero");
else if(!zero) _log -= mydbl._log;
return *this;
}
/*** ADDITION ***/
mydouble operator+(const double &dbl) const {
if(dbl==0.0) return mydouble(*this);
return operator+(mydouble(dbl));
}
mydouble operator+(const mydouble &mydbl) const {
mydouble a;
if(zero) a = mydouble(mydbl);
else if(mydbl.zero) a = mydouble(*this);
else {
double diff = _log - mydbl._log;
if(diff==0.0) a.setlog(log(2.0) + _log);
else if(diff<0.0) a.setlog(mydbl._log + log(1.0 + exp(diff)));
else a.setlog(_log + log(1.0 + exp(-diff)));
}
return a;
}
mydouble& operator+=(const double &dbl) {
if(dbl==0.0) return *this;
return operator+=(mydouble(dbl));
}
mydouble& operator+=(const mydouble &mydbl) {
if(zero) *this = mydbl;
else if(!mydbl.zero) {
double diff = _log - mydbl._log;
if(diff==0.0) _log += log(2.0);
else if(diff<0.0) _log = mydbl._log + log(1.0 + exp(diff));
else _log += log(1.0 + exp(-diff));
}
return *this;
}
/*** SUBTRACTION - warning cannot have negative numbers ***/
mydouble operator-(const double &dbl) const {
if(dbl==0.0) return mydouble(*this);
return operator-(mydouble(dbl));
}
mydouble operator-(const mydouble &mydbl) const {
mydouble a;
if(mydbl.zero) a = mydouble(*this);
else if(zero) error("mydouble::operator-(const mydouble&): subtracting a positive number from zero");
else {
/* diff must always be positive */
double diff = _log - mydbl._log;
if(diff==0.0) a.setzero();
else if(diff<0.0) myutils::error("mydouble::operator-(const mydouble&) cannot handle negative numbers");
else a.setlog(_log + log(1.0 - exp(-diff)));
}
return a;
}
mydouble& operator-=(const double &dbl) {
if(dbl==0.0) return *this;
return operator-=(mydouble(dbl));
}
mydouble& operator-=(const mydouble &mydbl) {
if(!mydbl.zero) {
if(zero) error("mydouble::operator-=(const mydouble&): subtracting a positive number from zero");
/* diff must always be positive */
double diff = _log - mydbl._log;
if(diff==0.0) setzero();
else if(diff<0.0) myutils::error("mydouble::operator-=(const mydouble&) cannot handle negative numbers");
else _log += log(1.0 - exp(-diff));
}
return *this;
}
/*** SPECIAL OPERATIONS ***/
double LOG() const {
return _log;
}
/* Caution: ^ has lower precedence than /*+- */
mydouble operator^(const double &dbl) const {
mydouble a;
if(zero) a.setzero();
else a.setlog(_log * dbl);
return a;
}
/* Caution: ^ has lower precedence than /*+- */
mydouble operator^(const mydouble &mydbl) const {
mydouble a;
if(zero) a.setzero();
else a.setlog(_log * exp(mydbl._log));
return a;
}
mydouble& operator^=(const double &dbl) {
if(!zero) _log *= dbl;
return *this;
}
mydouble& operator^=(const mydouble &mydbl) {
if(!zero) _log *= exp(mydbl._log);
return *this;
}
};
#endif//_MY_DOUBLE_H_