this repo has no description
dekker.one
1/* -*- mode: C++; c-basic-offset: 2; indent-tabs-mode: nil -*- */
2/*
3 * Main authors:
4 * Vincent Barichard <Vincent.Barichard@univ-angers.fr>
5 *
6 * Copyright:
7 * Vincent Barichard, 2012
8 *
9 * This file is part of Gecode, the generic constraint
10 * development environment:
11 * http://www.gecode.org
12 *
13 * Permission is hereby granted, free of charge, to any person obtaining
14 * a copy of this software and associated documentation files (the
15 * "Software"), to deal in the Software without restriction, including
16 * without limitation the rights to use, copy, modify, merge, publish,
17 * distribute, sublicense, and/or sell copies of the Software, and to
18 * permit persons to whom the Software is furnished to do so, subject to
19 * the following conditions:
20 *
21 * The above copyright notice and this permission notice shall be
22 * included in all copies or substantial portions of the Software.
23 *
24 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
25 * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
26 * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
27 * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
28 * LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
29 * OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
30 * WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
31 *
32 */
33
34namespace Gecode { namespace Float { namespace Trigonometric {
35 /*
36 * ATan projection function
37 *
38 */
39 template<class V>
40 void
41 aTanProject(Rounding& r, const V& aTanIv, FloatNum& iv_min, FloatNum& iv_max, int& n_min, int& n_max) {
42 #define I0_PI_2I FloatVal(0,pi_half_upper())
43 #define POS(X) ((I0_PI_2I.in(X))?0:1)
44 #define ATANINF_DOWN r.atan_down(aTanIv.min())
45 #define ATANSUP_UP r.atan_up(aTanIv.max())
46
47 // 0 <=> in [0;PI/2]
48 // 1 <=> in [PI/2;PI]
49 switch ( POS(iv_min) )
50 {
51 case 0:
52 if (r.tan_down(iv_min) > aTanIv.max()) { n_min++; iv_min = ATANINF_DOWN; }
53 else if (r.tan_up(iv_min) < aTanIv.min()) { iv_min = ATANINF_DOWN; }
54 break;
55 case 1:
56 if (r.tan_down(iv_min) > aTanIv.max()) { n_min+=2; iv_min = ATANINF_DOWN; }
57 else if (r.tan_up(iv_min) < aTanIv.min()) { n_min++; iv_min = ATANINF_DOWN; }
58 break;
59 default:
60 GECODE_NEVER;
61 break;
62 }
63
64 // 0 <=> in [0;PI/2]
65 // 1 <=> in [PI/2;PI]
66 switch ( POS(iv_max) )
67 {
68 case 0:
69 if (r.tan_down(iv_max) > aTanIv.max()) { iv_max = ATANSUP_UP; }
70 else if (r.tan_up(iv_max) < aTanIv.min()) { n_max--; iv_max = ATANSUP_UP; }
71 break;
72 case 1:
73 if (r.tan_down(iv_max) > aTanIv.max()) { n_max++; iv_max = ATANSUP_UP; }
74 else if (r.tan_up(iv_max) < aTanIv.min()) { iv_max = ATANSUP_UP; }
75 break;
76 default:
77 GECODE_NEVER;
78 break;
79 }
80 #undef ATANINF_DOWN
81 #undef ATANSUP_UP
82 #undef POS
83 #undef I0_PI_2I
84 }
85
86 /*
87 * Bounds consistent tangent operator
88 *
89 */
90
91 template<class A, class B>
92 ExecStatus
93 Tan<A,B>::dopropagate(Space& home, A x0, B x1) {
94 Rounding r;
95 int n_min = static_cast<int>(r.div_up(x0.min() + pi_half_upper(), pi_upper()));
96 int n_max = static_cast<int>(r.div_up(x0.max() + pi_half_upper(), pi_upper()));
97
98 if (x0 == x1) {
99 #define I0_PI_2I FloatVal(0,pi_half_upper())
100 if (I0_PI_2I.in(x0.max())) GECODE_ME_CHECK(x0.lq(home, 0));
101 if (I0_PI_2I.in(-x0.min())) GECODE_ME_CHECK(x0.gq(home, 0));
102 #undef I0_PI_2I
103
104 n_min = static_cast<int>(r.div_up(x0.min(), pi_upper()));
105 n_max = static_cast<int>(r.div_up(x0.max(), pi_upper()));
106
107 FloatNum x0_min;
108 FloatNum x0_max;
109 FloatNum t = x0.min();
110 do {
111 x0_min = t;
112 if (r.tan_down(x0_min) > x0_min) n_min++;
113 t = r.add_down(r.mul_up(n_min,pi_upper()),r.tan_down(x0_min));
114 } while (t > x0_min);
115 t = r.sub_down(r.mul_up(2*n_max,pi_upper()),x0.max());
116 do {
117 x0_max = t;
118 if (r.tan_down(x0_max) < x0_max) n_max--;
119 t = r.add_up(r.mul_up(n_max,pi_upper()),r.tan_up(x0_max));
120 } while (t > x0_max);
121 x0_max = r.sub_up(r.mul_up(2*n_max,pi_upper()),x0_max);
122
123 if (x0_min > x0_max) return ES_FAILED;
124 GECODE_ME_CHECK(x0.eq(home,FloatVal(x0_min,x0_max)));
125 } else {
126 GECODE_ME_CHECK(x1.eq(home,tan(x0.val())));
127 n_min = static_cast<int>(r.div_up(x0.min(), pi_upper()));
128 n_max = static_cast<int>(r.div_up(x0.max(), pi_upper()));
129 if (x0.min() < 0) n_min--;
130 if (x0.max() < 0) n_max--;
131 FloatNum iv_min = r.sub_down(x0.min(),r.mul_down(n_min, pi_upper()));
132 FloatNum iv_max = r.sub_up (x0.max(),r.mul_down(n_max, pi_upper()));
133 aTanProject(r,x1,iv_min,iv_max,n_min,n_max);
134 FloatNum n_iv_min = r.add_down(iv_min,r.mul_down(n_min, pi_upper()));
135 FloatNum n_iv_max = r.add_up (iv_max,r.mul_down(n_max, pi_upper()));
136 if (n_iv_min > n_iv_max) return ES_FAILED;
137 GECODE_ME_CHECK(x0.eq(home,FloatVal(n_iv_min,n_iv_max)));
138 GECODE_ME_CHECK(x1.eq(home,tan(x0.val()))); // Redo tan because with x0 reduction, sin may be more accurate
139 }
140
141 return ES_OK;
142 }
143
144 template<class A, class B>
145 forceinline
146 Tan<A,B>::Tan(Home home, A x0, B x1)
147 : MixBinaryPropagator<A,PC_FLOAT_BND,B,PC_FLOAT_BND>(home,x0,x1) {}
148
149 template<class A, class B>
150 ExecStatus
151 Tan<A,B>::post(Home home, A x0, B x1) {
152 if (x0 == x1) {
153 #define I0_PI_2I FloatVal(0,pi_half_upper())
154 if (I0_PI_2I.in(x0.max())) GECODE_ME_CHECK(x0.lq(home, 0));
155 if (I0_PI_2I.in(-x0.min())) GECODE_ME_CHECK(x0.gq(home, 0));
156 #undef I0_PI_2I
157 }
158 GECODE_ES_CHECK(dopropagate(home,x0,x1));
159 (void) new (home) Tan<A,B>(home,x0,x1);
160 return ES_OK;
161 }
162
163 template<class A, class B>
164 forceinline
165 Tan<A,B>::Tan(Space& home, Tan<A,B>& p)
166 : MixBinaryPropagator<A,PC_FLOAT_BND,B,PC_FLOAT_BND>(home,p) {}
167
168 template<class A, class B>
169 Actor*
170 Tan<A,B>::copy(Space& home) {
171 return new (home) Tan<A,B>(home,*this);
172 }
173
174 template<class A, class B>
175 ExecStatus
176 Tan<A,B>::propagate(Space& home, const ModEventDelta&) {
177 GECODE_ES_CHECK(dopropagate(home,x0,x1));
178 return (x0.assigned()) ? home.ES_SUBSUMED(*this) : ES_FIX;
179 }
180
181 /*
182 * Bounds consistent arc tangent operator
183 *
184 */
185
186 template<class A, class B>
187 forceinline
188 ATan<A,B>::ATan(Home home, A x0, B x1)
189 : MixBinaryPropagator<A,PC_FLOAT_BND,B,PC_FLOAT_BND>(home,x0,x1) {}
190
191 template<class A, class B>
192 ExecStatus
193 ATan<A,B>::post(Home home, A x0, B x1) {
194 if (x0 == x1) {
195 GECODE_ME_CHECK(x0.eq(home,0.0));
196 } else {
197 GECODE_ME_CHECK(x1.eq(home,atan(x0.domain())));
198 GECODE_ME_CHECK(x0.eq(home,tan(x1.domain())));
199 (void) new (home) ATan<A,B>(home,x0,x1);
200 }
201 return ES_OK;
202 }
203
204
205 template<class A, class B>
206 forceinline
207 ATan<A,B>::ATan(Space& home, ATan<A,B>& p)
208 : MixBinaryPropagator<A,PC_FLOAT_BND,B,PC_FLOAT_BND>(home,p) {}
209
210 template<class A, class B>
211 Actor*
212 ATan<A,B>::copy(Space& home) {
213 return new (home) ATan<A,B>(home,*this);
214 }
215
216 template<class A, class B>
217 ExecStatus
218 ATan<A,B>::propagate(Space& home, const ModEventDelta&) {
219 GECODE_ME_CHECK(x1.eq(home,atan(x0.domain())));
220 GECODE_ME_CHECK(x0.eq(home,tan(x1.domain())));
221 return (x0.assigned() && x1.assigned()) ? home.ES_SUBSUMED(*this) : ES_FIX;
222 }
223
224}}}
225
226// STATISTICS: float-prop
227