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
|
// miscfp.cc -- misc floating-point functions
//
// jtutil::signum - sign of a floating point number
// jtutil::hypot3 - 3D Euclidean distance
// jtutil::arctan_xy - 4-quadrant arc tangent
// jtutil::modulo_reduce - reduce an angle modulo 2*pi radians (360 degrees)
// jtutil::zero_C_array - set a C-style array to all zeros
//
// ***** template instantiations *****
//
#include <math.h>
#include <stdlib.h>
// we want to instantiate templates with CCTK_* types
#ifdef STANDALONE_TEST
#include "fake_cctk.h"
#else
#include "cctk.h"
#endif
#include "stdc.h"
#include "util.hh"
// everything in this file is inside these namespaces
namespace AHFinderDirect
{
namespace jtutil
{
//******************************************************************************
//******************************************************************************
//******************************************************************************
//
// This function computes the floating point "signum" function (as in APL),
// signum(x) = -1.0, 0.0, or +1.0, according to the sign of x.
//
double signum(double x)
{
if (x == 0.0)
then return 0.0;
else return (x > 0.0) ? 1.0 : -1.0;
}
//******************************************************************************
//
// This function computes the 3-dimensional Euclidean distance,
// analogously to the standard math library function hypot(2) .
//
// Arguments:
// (x,y,z) = (in) The rectangular coordinates of a point in $\Re^3$.
//
// Result:
// The function returns the Euclidean distance of (x,y,z) from the origin.
//
// Bugs:
// Unlike hypot(2), this function takes no special care to avoid
// unwarranted IEEE exceptions if any of |x|, |y|, or |z| is close to
// the overflow and/or underflow threshold.
//
double hypot3(double x, double y, double z)
{
return sqrt(x*x + y*y + z*z);
}
//******************************************************************************
//
// This function computes a four-quadrant arctangent, using what I think
// is the "right" ordering of the arguments, and returning 0.0 if both
// arguments are 0.0.
//
// Arguments:
// (x,y) = (in) The rectangular coordinates of a point in $\Re^2$.
//
// Result:
// The function returns a value $\theta$ such that $x + iy = R e^{i\theta}$ for
// some real $R$, i.e. it returns the angle between the positive $x$ axis and
// the line joining the origin and the point $(x,y)$.
//
double arctan_xy(double x, double y)
{
// note reversed argument order (y,x) in std::atan2() function
return ((x == 0.0) && (y == 0.0)) ? 0.0 : atan2(y,x);
}
//******************************************************************************
//
// This function reduces x modulo xmod to be (fuzzily) in the range
// [xmin, xmax] , or does an error_exit() if no such value exists.
//
double modulo_reduce(double x, double xmod, double xmin, double xmax)
{
double xx = x;
while (fuzzy<double>::LT(xx, xmin))
{
xx += xmod;
}
while (fuzzy<double>::GT(xx, xmax))
{
xx -= xmod;
}
if (! (fuzzy<double>::GE(xx, xmin) && fuzzy<double>::LE(xx, xmax)) )
then error_exit(ERROR_EXIT,
"***** modulo_reduce(): no modulo value is fuzzily within specified range!\n"
" x = %g xmod = %g\n"
" [xmin,xmax] = [%g,%g]\n"
" ==> xx = %g\n"
,
x, xmod,
xmin, xmax,
xx); /*NOTREACHED*/
return xx;
}
//******************************************************************************
//
// This function sets a C-style array to all zeros.
//
template <typename fp_t>
void zero_C_array(int N, fp_t array[])
{
for (int i = 0 ; i < N ; ++i)
{
array[i] = 0;
}
}
//******************************************************************************
//******************************************************************************
//******************************************************************************
//
// ***** template instantiations *****
//
template void zero_C_array<CCTK_REAL>(int, CCTK_REAL[]);
//******************************************************************************
//******************************************************************************
//******************************************************************************
} // namespace jtutil
} // namespace AHFinderDirect
|
