1 /*++
2 Copyright (c) 2012 Microsoft Corporation
3
4 Module Name:
5
6 trigo.cpp
7
8 Abstract:
9
10 Test trigonometric primitives in the interval class.
11
12 Author:
13
14 Leonardo de Moura (leonardo) 2012-08-20
15
16 Revision History:
17
18 --*/
19 #include<cstdlib>
20 #include "math/interval/interval_def.h"
21 #include "util/dependency.h"
22 #include "util/mpq.h"
23 #include "ast/ast.h"
24 #include "util/debug.h"
25 #include "test/im_float_config.h"
26 #include "util/rlimit.h"
27
28 #define PREC 100000
29
tst_sine_core(std::ostream & out,unsynch_mpq_manager & nm,interval_manager<im_default_config> & im,mpq & a,unsigned k)30 static void tst_sine_core(std::ostream & out, unsynch_mpq_manager & nm, interval_manager<im_default_config> & im, mpq & a, unsigned k) {
31 scoped_mpq lo(nm), hi(nm);
32 im.sine(a, k, lo, hi);
33 nm.display(out, lo);
34 out << " <= Sin["; nm.display(out, a); out << "]\n";
35 out << "Sin["; nm.display(out, a); out << "] <= ";
36 nm.display(out, hi);
37 out << "\n";
38 }
39
tst_sine(std::ostream & out,unsigned N,unsigned k)40 static void tst_sine(std::ostream & out, unsigned N, unsigned k) {
41 unsynch_mpq_manager nm;
42 reslimit rl;
43 interval_manager<im_default_config> im(rl, nm);
44 scoped_mpq a(nm);
45 nm.set(a, 0);
46 tst_sine_core(out, nm, im, a, 1);
47 for (unsigned i = 0; i < N; i++) {
48 nm.set(a, 4 * (rand() % PREC), PREC);
49 if (rand() % 2 == 0)
50 nm.neg(a);
51 tst_sine_core(out, nm, im, a, k);
52 }
53 }
54
55
tst_cosine_core(std::ostream & out,unsynch_mpq_manager & nm,interval_manager<im_default_config> & im,mpq & a,unsigned k)56 static void tst_cosine_core(std::ostream & out, unsynch_mpq_manager & nm, interval_manager<im_default_config> & im, mpq & a, unsigned k) {
57 scoped_mpq lo(nm), hi(nm);
58 im.cosine(a, k, lo, hi);
59 nm.display(out, lo);
60 out << " <= Cos["; nm.display(out, a); out << "]\n";
61 out << "Cos["; nm.display(out, a); out << "] <= ";
62 nm.display(out, hi);
63 out << "\n";
64 }
65
tst_cosine(std::ostream & out,unsigned N,unsigned k)66 static void tst_cosine(std::ostream & out, unsigned N, unsigned k) {
67 reslimit rl;
68 unsynch_mpq_manager nm;
69 interval_manager<im_default_config> im(rl, nm);
70 scoped_mpq a(nm);
71 nm.set(a, 0);
72 tst_cosine_core(out, nm, im, a, 1);
73 for (unsigned i = 0; i < N; i++) {
74 nm.set(a, 4 * (rand() % PREC), PREC);
75 if (rand() % 2 == 0)
76 nm.neg(a);
77 tst_cosine_core(out, nm, im, a, k);
78 }
79 }
80
81
82 template<typename fmanager>
tst_float_sine_core(std::ostream & out,fmanager & fm,interval_manager<im_float_config<fmanager>> & im,typename fmanager::numeral & a,unsigned k)83 static void tst_float_sine_core(std::ostream & out,
84 fmanager & fm,
85 interval_manager<im_float_config<fmanager> > & im,
86 typename fmanager::numeral & a,
87 unsigned k) {
88 _scoped_numeral<fmanager> lo(fm), hi(fm);
89 im.sine(a, k, lo, hi);
90 out << fm.to_rational_string(lo) << " <= Sin[" << fm.to_rational_string(a) << "]\n";
91 out << "Sin[" << fm.to_rational_string(a) << "] <= " << fm.to_rational_string(hi) << "\n";
92 }
93
94 const unsigned EBITS = 11;
95 const unsigned SBITS = 53;
96
97 template<typename fmanager>
tst_float_sine(std::ostream & out,unsigned N,unsigned k)98 static void tst_float_sine(std::ostream & out, unsigned N, unsigned k) {
99 reslimit rl;
100 fmanager fm;
101 interval_manager<im_float_config<fmanager> > im(rl, im_float_config<fmanager>(fm, EBITS, SBITS));
102 _scoped_numeral<fmanager> a(fm);
103 fm.set(a, EBITS, SBITS, static_cast<int>(0));
104 tst_float_sine_core(out, fm, im, a, 1);
105
106 // fm.set(a, EBITS, SBITS, MPF_ROUND_TOWARD_POSITIVE, 25336, 100000);
107 // tst_float_sine_core(out, fm, im, a, k);
108 // return;
109 for (unsigned i = 0; i < N; i++) {
110 unsigned n = 4 * (rand() % PREC);
111 unsigned d = PREC;
112 TRACE("sine", tout << "next-val : " << n << "/" << d << "\n";);
113 fm.set(a, EBITS, SBITS, MPF_ROUND_TOWARD_POSITIVE, n, d);
114 if (rand() % 2 == 0)
115 fm.neg(a);
116 tst_float_sine_core(out, fm, im, a, k);
117 }
118 }
119
120 #if 0
121 static void tst_mpf_bug() {
122 mpf_manager fm;
123 scoped_mpf a(fm), b(fm), c(fm);
124 fm.set(a, EBITS, SBITS, 2);
125 fm.set(b, EBITS, SBITS, 3);
126 std::cout << "a: " << fm.to_double(a) << "\n";
127 std::cout << "b: " << fm.to_double(b) << "\n";
128 fm.mul(MPF_ROUND_TOWARD_NEGATIVE, a, b, c);
129 std::cout << "c: " << fm.to_double(c) << "\n";
130 }
131 #endif
132
tst_e(std::ostream & out)133 static void tst_e(std::ostream & out) {
134 reslimit rl;
135 unsynch_mpq_manager nm;
136 interval_manager<im_default_config> im(rl, nm);
137 im_default_config::interval r;
138 for (unsigned i = 0; i < 64; i++) {
139 im.e(i, r);
140 out << nm.to_string(im.lower(r)) << " <= E\n";
141 out << "E <= " << nm.to_string(im.upper(r)) << "\n";
142 }
143 im.del(r);
144 }
145
tst_e_float(std::ostream & out)146 static void tst_e_float(std::ostream & out) {
147 std::cout << "e float...\n";
148 reslimit rl;
149 unsynch_mpq_manager qm;
150 mpf_manager fm;
151 interval_manager<im_float_config<mpf_manager> > im(rl, fm);
152 scoped_mpq q(qm);
153 im_float_config<mpf_manager>::interval r;
154 for (unsigned i = 0; i < 64; i++) {
155 im.e(i, r);
156 out << fm.to_rational_string(im.lower(r)) << " <= E\n";
157 out << "E <= " << fm.to_rational_string(im.upper(r)) << "\n";
158 }
159 im.del(r);
160 }
161
tst_trigo()162 void tst_trigo() {
163 // enable_trace("sine");
164 // enable_trace("sine_bug");
165 // enable_trace("mpf_mul_bug");
166 std::ofstream out("trigo-lemmas.math");
167 tst_e_float(out);
168 tst_e(out);
169 tst_float_sine<mpf_manager>(out, 100, 5);
170 tst_float_sine<mpf_manager>(out, 100, 7);
171 tst_sine(out, 200, 3);
172 tst_sine(out, 200, 5);
173 tst_sine(out, 200, 9);
174 tst_cosine(out, 200, 3);
175 tst_cosine(out, 200, 5);
176 tst_cosine(out, 200, 9);
177 }
178