polysolver-inl.h
Go to the documentation of this file.
1 /*
2  * Software License Agreement (BSD License)
3  *
4  * Copyright (c) 2011-2014, Willow Garage, Inc.
5  * Copyright (c) 2014-2016, Open Source Robotics Foundation
6  * All rights reserved.
7  *
8  * Redistribution and use in source and binary forms, with or without
9  * modification, are permitted provided that the following conditions
10  * are met:
11  *
12  * * Redistributions of source code must retain the above copyright
13  * notice, this list of conditions and the following disclaimer.
14  * * Redistributions in binary form must reproduce the above
15  * copyright notice, this list of conditions and the following
16  * disclaimer in the documentation and/or other materials provided
17  * with the distribution.
18  * * Neither the name of Open Source Robotics Foundation nor the names of its
19  * contributors may be used to endorse or promote products derived
20  * from this software without specific prior written permission.
21  *
22  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
23  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
24  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
25  * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
26  * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
27  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
28  * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
29  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
30  * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
32  * ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
33  * POSSIBILITY OF SUCH DAMAGE.
34  */
35 
38 #ifndef FCL_NARROWPHASE_DETAIL_POLYSOLVER_INL_H
39 #define FCL_NARROWPHASE_DETAIL_POLYSOLVER_INL_H
40 
42 
43 #include <cmath>
44 #include "fcl/common/types.h"
45 
46 namespace fcl
47 {
48 
49 namespace detail {
50 
51 //==============================================================================
52 extern template
53 class FCL_EXPORT PolySolver<double>;
54 
55 //==============================================================================
56 template <typename S>
57 int PolySolver<S>::solveLinear(S c[2], S s[1])
58 {
59  if(isZero(c[1]))
60  return 0;
61  s[0] = - c[0] / c[1];
62  return 1;
63 }
64 
65 //==============================================================================
66 template <typename S>
67 int PolySolver<S>::solveQuadric(S c[3], S s[2])
68 {
69  S p, q, D;
70 
71  // make sure we have a d2 equation
72 
73  if(isZero(c[2]))
74  return solveLinear(c, s);
75 
76  // normal for: x^2 + px + q
77  p = c[1] / (2.0 * c[2]);
78  q = c[0] / c[2];
79  D = p * p - q;
80 
81  if(isZero(D))
82  {
83  // one S root
84  s[0] = s[1] = -p;
85  return 1;
86  }
87 
88  if(D < 0.0)
89  // no real root
90  return 0;
91  else
92  {
93  // two real roots
94  S sqrt_D = sqrt(D);
95  s[0] = sqrt_D - p;
96  s[1] = -sqrt_D - p;
97  return 2;
98  }
99 }
100 
101 //==============================================================================
102 template <typename S>
103 int PolySolver<S>::solveCubic(S c[4], S s[3])
104 {
105  int i, num;
106  S sub, A, B, C, sq_A, p, q, cb_p, D;
107  const S ONE_OVER_THREE = 1 / 3.0;
108  const S PI = 3.14159265358979323846;
109 
110  // make sure we have a d2 equation
111  if(isZero(c[3]))
112  return solveQuadric(c, s);
113 
114  // normalize the equation:x ^ 3 + Ax ^ 2 + Bx + C = 0
115  A = c[2] / c[3];
116  B = c[1] / c[3];
117  C = c[0] / c[3];
118 
119  // substitute x = y - A / 3 to eliminate the quadratic term: x^3 + px + q = 0
120  sq_A = A * A;
121  p = (-ONE_OVER_THREE * sq_A + B) * ONE_OVER_THREE;
122  q = 0.5 * (2.0 / 27.0 * A * sq_A - ONE_OVER_THREE * A * B + C);
123 
124  // use Cardano's formula
125  cb_p = p * p * p;
126  D = q * q + cb_p;
127 
128  if(isZero(D))
129  {
130  if(isZero(q))
131  {
132  // one triple solution
133  s[0] = 0.0;
134  num = 1;
135  }
136  else
137  {
138  // one single and one S solution
139  S u = cbrt(-q);
140  s[0] = 2.0 * u;
141  s[1] = -u;
142  num = 2;
143  }
144  }
145  else
146  {
147  if(D < 0.0)
148  {
149  // three real solutions
150  S phi = ONE_OVER_THREE * acos(-q / sqrt(-cb_p));
151  S t = 2.0 * sqrt(-p);
152  s[0] = t * cos(phi);
153  s[1] = -t * cos(phi + PI / 3.0);
154  s[2] = -t * cos(phi - PI / 3.0);
155  num = 3;
156  }
157  else
158  {
159  // one real solution
160  S sqrt_D = sqrt(D);
161  S u = cbrt(sqrt_D + fabs(q));
162  if(q > 0.0)
163  s[0] = - u + p / u ;
164  else
165  s[0] = u - p / u;
166  num = 1;
167  }
168  }
169 
170  // re-substitute
171  sub = ONE_OVER_THREE * A;
172  for(i = 0; i < num; i++)
173  s[i] -= sub;
174  return num;
175 }
176 
177 //==============================================================================
178 template <typename S>
180 {
181  return (v < getNearZeroThreshold()) && (v > -getNearZeroThreshold());
182 }
183 
184 //==============================================================================
185 template <typename S>
187 {
188  return std::pow(v, 1.0 / 3.0);
189 }
190 
191 //==============================================================================
192 template <typename S>
194 {
195  return 1e-9;
196 }
197 
198 } // namespace detail
199 } // namespace fcl
200 
201 #endif
fcl::detail::PolySolver::solveQuadric
static int solveQuadric(S c[3], S s[2])
Solve a quadratic function with coefficients c, return roots s and number of roots.
Definition: polysolver-inl.h:67
fcl::detail::PolySolver::solveLinear
static int solveLinear(S c[2], S s[1])
Solve a linear equation with coefficients c, return roots s and number of roots.
Definition: polysolver-inl.h:57
types.h
polysolver.h
fcl::detail::PolySolver::getNearZeroThreshold
static constexpr S getNearZeroThreshold()
Definition: polysolver-inl.h:193
fcl::detail::PolySolver< double >
template class FCL_EXPORT PolySolver< double >
fcl::detail::PolySolver::isZero
static bool isZero(S v)
Check whether v is zero.
Definition: polysolver-inl.h:179
fcl::detail::PolySolver::cbrt
static bool cbrt(S v)
Compute v^{1/3}.
Definition: polysolver-inl.h:186
fcl
Main namespace.
Definition: broadphase_bruteforce-inl.h:45
fcl::detail::PolySolver::solveCubic
static int solveCubic(S c[4], S s[3])
Solve a cubic function with coefficients c, return roots s and number of roots.
Definition: polysolver-inl.h:103


fcl
Author(s):
autogenerated on Tue Dec 5 2023 03:40:48