Point Cloud Library (PCL)  1.9.1-dev
opennurbs_math.h
1 /* $NoKeywords: $ */
2 /*
3 //
4 // Copyright (c) 1993-2012 Robert McNeel & Associates. All rights reserved.
5 // OpenNURBS, Rhinoceros, and Rhino3D are registered trademarks of Robert
6 // McNeel & Associates.
7 //
8 // THIS SOFTWARE IS PROVIDED "AS IS" WITHOUT EXPRESS OR IMPLIED WARRANTY.
9 // ALL IMPLIED WARRANTIES OF FITNESS FOR ANY PARTICULAR PURPOSE AND OF
10 // MERCHANTABILITY ARE HEREBY DISCLAIMED.
11 //
12 // For complete openNURBS copyright information see <http://www.opennurbs.org>.
13 //
14 ////////////////////////////////////////////////////////////////
15 */
16 
17 #if !defined(ON_MATH_INC_)
18 #define ON_MATH_INC_
19 
20 class ON_3dVector;
21 class ON_Interval;
22 class ON_Line;
23 class ON_Arc;
24 class ON_Plane;
25 
26 /*
27 Description:
28  Class for carefully adding long list of numbers.
29 */
30 class ON_CLASS ON_Sum
31 {
32 public:
33 
34  /*
35  Description:
36  Calls ON_Sum::Begin(x)
37  */
38  void operator=(double x);
39 
40  /*
41  Description:
42  Calls ON_Sum::Plus(x);
43  */
44  void operator+=(double x);
45 
46  /*
47  Description:
48  Calls ON_Sum::Plus(-x);
49  */
50  void operator-=(double x);
51 
52  /*
53  Description:
54  Creates a sum that is ready to be used.
55  */
56  ON_Sum();
57 
58  /*
59  Description:
60  If a sum is being used more than once, call Begin()
61  before starting each sum.
62  Parameters:
63  starting_value - [in] Initial value of sum.
64  */
65  void Begin( double starting_value = 0.0 );
66 
67  /*
68  Description:
69  Add x to the current sum.
70  Parameters:
71  x - [in] value to add to the current sum.
72  */
73  void Plus( double x );
74 
75  /*
76  Description:
77  Calculates the total sum.
78  Parameters:
79  error_estimate - [out] if not NULL, the returned value of
80  *error_estimate is an estimate of the error in the sum.
81  Returns:
82  Total of the sum.
83  Remarks:
84  You can get subtotals by mixing calls to Plus() and Total().
85  In delicate sums, some precision may be lost in the final
86  total if you call Total() to calculate subtotals.
87  */
88  double Total( double* error_estimate = NULL );
89 
90  /*
91  Returns:
92  Number of summands.
93  */
94  int SummandCount() const;
95 
96 private:
97  enum {
98  sum1_max_count=256,
99  sum2_max_count=512,
100  sum3_max_count=1024
101  };
102  double m_sum_err;
103  double m_pos_sum;
104  double m_neg_sum;
105 
106  int m_zero_count; // number of zeros added
107  int m_pos_count; // number of positive numbers added
108  int m_neg_count; // number of negative numbers added
109 
110  int m_pos_sum1_count;
111  int m_pos_sum2_count;
112  int m_pos_sum3_count;
113  double m_pos_sum1[sum1_max_count];
114  double m_pos_sum2[sum2_max_count];
115  double m_pos_sum3[sum3_max_count];
116 
117  int m_neg_sum1_count;
118  int m_neg_sum2_count;
119  int m_neg_sum3_count;
120  double m_neg_sum1[sum1_max_count];
121  double m_neg_sum2[sum2_max_count];
122  double m_neg_sum3[sum3_max_count];
123 
124  double SortAndSum( int, double* );
125 };
126 
127 /*
128 Description:
129  Abstract function with an arbitrary number of parameters
130  and values. ON_Evaluator is used to pass functions to
131  local solvers.
132 */
133 class ON_CLASS ON_Evaluator
134 {
135 public:
136 
137  /*
138  Description:
139  Construction of the class for a function that takes
140  parameter_count input functions and returns
141  value_count values. If the domain is infinite, pass
142  a NULL for the domain[] and periodic[] arrays. If
143  the domain is finite, pass a domain[] array with
144  parameter_count increasing intervals. If one or more of
145  the parameters is periodic, pass the fundamental domain
146  in the domain[] array and a true in the periodic[] array.
147  Parameters:
148  parameter_count - [in] >= 1. Number of input parameters
149  value_count - [in] >= 1. Number of output values.
150  domain - [in] If not NULL, then this is an array
151  of parameter_count increasing intervals
152  that defines the domain of the function.
153  periodic - [in] if not NULL, then this is an array of
154  parameter_count bools where b[i] is true if
155  the i-th parameter is periodic. Valid
156  increasing finite domains must be specificed
157  when this parameter is not NULL.
158  */
159  ON_Evaluator(
160  int parameter_count,
161  int value_count,
162  const ON_Interval* domain,
163  const bool* periodic
164  );
165 
166  virtual ~ON_Evaluator();
167 
168  /*
169  Description:
170  Evaluate the function that takes m_parameter_count parameters
171  and returns a m_value_count dimensional point.
172  Parameters:
173  parameters - [in] array of m_parameter_count evaluation parameters
174  values - [out] array of m_value_count function values
175  jacobian - [out] If NULL, simply evaluate the value of the function.
176  If not NULL, this is the jacobian of the function.
177  jacobian[i][j] = j-th partial of the i-th value
178  0 <= i < m_value_count,
179  0 <= j < m_parameter_count
180  If not NULL, then all the memory for the
181  jacobian is allocated, you just need to fill
182  in the answers.
183  Example:
184  If f(u,v) = square of the distance from a fixed point P to a
185  surface evaluated at (u,v), then
186 
187  values[0] = (S-P)o(S-P)
188  jacobian[0] = ( 2*(Du o (S-P)), 2*(Dv o (S-P)) )
189 
190  where S, Du, Dv = surface point and first partials evaluated
191  at u=parameters[0], v = parameters[1].
192 
193  If the function takes 3 parameters, say (x,y,z), and returns
194  two values, say f(x,y,z) and g(z,y,z), then
195 
196  values[0] = f(x,y,z)
197  values[1] = g(x,y,z)
198 
199  jacobian[0] = (DfDx, DfDy, DfDz)
200  jacobian[1] = (DgDx, DgDy, DgDz)
201 
202  where dfx denotes the first partial of f with respect to x.
203 
204  Returns:
205  0 = unable to evaluate
206  1 = successful evaluation
207  2 = found answer, terminate search
208  */
209  virtual int Evaluate(
210  const double* parameters,
211  double* values,
212  double** jacobian
213  ) = 0;
214 
215  /*
216  Description:
217  OPTIONAL ability to evaluate the hessian in the case when
218  m_value_count is one. If your function has more that
219  one value or it is not feasable to evaluate the hessian,
220  then do not override this function. The default implementation
221  returns -1.
222  Parameters:
223  parameters - [in] array of m_parameter_count evaluation parameters
224  value - [out] value of the function (one double)
225  gradient - [out] The gradient of the function. This is a vector
226  of length m_parameter_count; gradient[i] is
227  the first partial of the function with respect to
228  the i-th parameter.
229  hessian - [out] The hessian of the function. This is an
230  m_parameter_count x m_parameter_count
231  symmetric matrix: hessian[i][j] is the
232  second partial of the function with respect
233  to the i-th and j-th parameters. The evaluator
234  is responsible for filling in both the upper
235  and lower triangles. Since the matrix is
236  symmetrix, you should do something like evaluate
237  the upper triangle and copy the values to the
238  lower tiangle.
239  Returns:
240  -1 = Hessian evaluation not available.
241  0 = unable to evaluate
242  1 = successful evaluation
243  2 = found answer, terminate search
244  */
245  virtual int EvaluateHessian(
246  const double* parameters,
247  double* value,
248  double* gradient,
249  double** hessian
250  );
251 
252  // Number of the function's input parameters. This number
253  // is >= 1 and is specified in the constructor.
254  const int m_parameter_count;
255 
256  // Number of the function's output values. This number
257  // is >= 1 and is specified in the constructor.
258  const int m_value_count;
259 
260  /*
261  Description:
262  Functions can have finite or infinite domains. Finite domains
263  are specified by passing the domain[] array to the constructor
264  or filling in the m_domain[] member variable. If
265  m_domain.Count() == m_parameter_count > 0, then the function
266  has finite domains.
267  Returns:
268  True if the domain of the function is finite.
269  */
270  bool FiniteDomain() const;
271 
272  /*
273  Description:
274  If a function has a periodic parameter, then the m_domain
275  interval for that parameter is the fundamental domain and
276  the m_bPeriodicParameter bool for that parameter is true.
277  A parameter is periodic if, and only if,
278  m_domain.Count() == m_parameter_count, and
279  m_bPeriodicParameter.Count() == m_parameter_count, and
280  m_bPeriodicParameter[parameter_index] is true.
281  Returns:
282  True if the function parameter is periodic.
283  */
284  bool Periodic(
285  int parameter_index
286  ) const;
287 
288  /*
289  Description:
290  If a function has a periodic parameter, then the m_domain
291  interval for that parameter is the fundamental domain and
292  the m_bPeriodicParameter bool for that parameter is true.
293  A parameter is periodic if, and only if,
294  m_domain.Count() == m_parameter_count, and
295  m_bPeriodicParameter.Count() == m_parameter_count, and
296  m_bPeriodicParameter[parameter_index] is true.
297  Returns:
298  The domain of the parameter. If the domain is infinite,
299  the (-1.0e300, +1.0e300) is returned.
300  */
301  ON_Interval Domain(
302  int parameter_index
303  ) const;
304 
305 
306  // If the function has a finite domain or periodic
307  // parameters, then m_domain[] is an array of
308  // m_parameter_count finite increasing intervals.
310 
311  // If the function has periodic parameters, then
312  // m_bPeriodicParameter[] is an array of m_parameter_count
313  // bools. If m_bPeriodicParameter[i] is true, then
314  // the i-th parameter is periodic and m_domain[i] is
315  // the fundamental domain for that parameter.
317 
318 private:
319  ON_Evaluator(); // prohibit default constructor
320  ON_Evaluator& operator=(const ON_Evaluator&); // prohibit operator= (can't copy const members)
321 };
322 
323 /*
324 Description:
325  Test a double to make sure it is a valid number.
326 Returns:
327  True if x != ON_UNSET_VALUE and _finite(x) is true.
328 */
329 ON_DECL
330 bool ON_IsValid( double x );
331 
332 ON_DECL
333 bool ON_IsValidFloat( float x );
334 
335 /*
336 class ON_CLASS ON_TimeLimit
337 {
338  ON_TimeLimit();
339  ON_TimeLimit(ON__UINT64 time_limit_seconds);
340  void SetTimeLimit(ON__UINT64 time_limit_seconds);
341  bool Continue() const;
342  bool IsSet() const;
343 private:
344  ON__UINT64 m_time_limit[2];
345 };
346 */
347 
348 // The ON_IS_FINITE and ON_IS_VALID defines are much faster
349 // than calling ON_IsValid(), but need to be used when
350 // the macro expansion works.
351 
352 #if defined(ON_LITTLE_ENDIAN)
353 
354 // works on little endian CPUs with IEEE doubles
355 #define ON_IS_FINITE(x) (0x7FF0 != (*((unsigned short*)(&x) + 3) & 0x7FF0))
356 #define ON_IS_VALID(x) (x != ON_UNSET_VALUE && 0x7FF0 != (*((unsigned short*)(&x) + 3) & 0x7FF0))
357 #define ON_IS_VALID_FLOAT(x) (x != ON_UNSET_FLOAT)
358 //TODO - ADD FAST ugly bit check#define ON_IS_VALID_FLOAT(x) (x != ON_UNSET_FLOAT && 0x7FF0 != (*((unsigned short*)(&x) + 3) & 0x7FF0))
359 
360 #elif defined(ON_BIG_ENDIAN)
361 
362 // works on big endian CPUs with IEEE doubles
363 #define ON_IS_FINITE(x) (0x7FF0 != (*((unsigned short*)(&x)) & 0x7FF0))
364 #define ON_IS_VALID(x) (x != ON_UNSET_VALUE && 0x7FF0 != (*((unsigned short*)(&x)) & 0x7FF0))
365 #define ON_IS_VALID_FLOAT(x) (x != ON_UNSET_FLOAT)
366 //TODO - ADD FAST ugly bit check#define ON_IS_VALID_FLOAT(x) (x != ON_UNSET_FLOAT && 0x7FF0 != (*((unsigned short*)(&x) + 3) & 0x7FF0))
367 
368 #else
369 
370 // Returns true if x is a finite double. Specifically,
371 // _finite returns a nonzero value (true) if its argument x
372 // is not infinite, that is, if -INF < x < +INF.
373 // It returns 0 (false) if the argument is infinite or a NaN.
374 //
375 // If you are trying to compile opennurbs on a platform
376 // that does not support finite(), then see if you can
377 // use _fpclass(), fpclass(), _isnan(), or isnan(). If
378 // you can't find anything, then just set this
379 // function to return true.
380 
381 #if defined(_GNU_SOURCE)
382 // if you are using an older version of gcc, use finite()
383 //#define ON_IS_FINITE(x) (finite(x)?true:false)
384 #define ON_IS_FINITE(x) (isfinite(x)?true:false)
385 #else
386 #define ON_IS_FINITE(x) (_finite(x)?true:false)
387 #endif
388 
389 #define ON_IS_VALID(x) (x != ON_UNSET_VALUE && ON_IS_FINITE(x))
390 #define ON_IS_VALID_FLOAT(x) (x != ON_UNSET_FLOAT && ON_IS_FINITE(x))
391 
392 #endif
393 
394 
395 ON_DECL
396 float ON_ArrayDotProduct( // returns AoB
397  int, // size of arrays (can be zero)
398  const float*, // A[]
399  const float* // B[]
400  );
401 
402 ON_DECL
403 void ON_ArrayScale(
404  int, // size of arrays (can be zero)
405  float, // a
406  const float*, // A[]
407  float* // returns a*A[]
408  );
409 
410 ON_DECL
411 void ON_Array_aA_plus_B(
412  int, // size of arrays (can be zero)
413  float, // a
414  const float*, // A[]
415  const float*, // B[]
416  float* // returns a*A[] + B[]
417  );
418 
419 ON_DECL
420 double ON_ArrayDotProduct( // returns AoB
421  int, // size of arrays (can be zero)
422  const double*, // A[]
423  const double* // B[]
424  );
425 
426 ON_DECL
427 double ON_ArrayDotDifference( // returns A o ( B - C )
428  int, // size of arrays (can be zero)
429  const double*, // A[]
430  const double*, // B[]
431  const double* // C[]
432  );
433 
434 ON_DECL
435 double ON_ArrayMagnitude( // returns sqrt(AoA)
436  int, // size of arrays (can be zero)
437  const double* // A[]
438  );
439 
440 ON_DECL
441 double ON_ArrayMagnitudeSquared( // returns AoA
442  int, // size of arrays (can be zero)
443  const double* // A[]
444  );
445 
446 ON_DECL
447 double ON_ArrayDistance( // returns sqrt((A-B)o(A-B))
448  int, // size of arrays (can be zero)
449  const double*, // A[]
450  const double* // B[]
451  );
452 
453 ON_DECL
454 double ON_ArrayDistanceSquared( // returns (A-B)o(A-B)
455  int, // size of arrays (can be zero)
456  const double*, // A[]
457  const double* // B[]
458  );
459 
460 ON_DECL
461 void ON_ArrayScale(
462  int, // size of arrays (can be zero)
463  double, // a
464  const double*, // A[]
465  double* // returns a*A[]
466  );
467 
468 ON_DECL
469 void ON_Array_aA_plus_B(
470  int, // size of arrays (can be zero)
471  double, // a
472  const double*, // A[]
473  const double*, // B[]
474  double* // returns a*A[] + B[]
475  );
476 
477 ON_DECL
478 int ON_SearchMonotoneArray( // find a value in an increasing array
479  // returns -1: t < array[0]
480  // i: array[i] <= t < array[i+1] ( 0 <= i < length-1 )
481  // length-1: t == array[length-1]
482  // length: t >= array[length-1]
483  const double*, // array[]
484  int, // length of array
485  double // t = value to search for
486  );
487 
488 
489 /*
490 Description:
491  Compute a binomial coefficient.
492 Parameters:
493  i - [in]
494  j - [in]
495 Returns:
496  (i+j)!/(i!j!), if 0 <= i and 0 <= j, and 0 otherwise.
497 See Also:
498  ON_TrinomialCoefficient()
499 Remarks:
500  If (i+j) <= 52, this function is fast and returns the exact
501  value of the binomial coefficient.
502 
503  For (i+j) > 52, the coefficient is computed recursively using
504  the formula bc(i,j) = bc(i-1,j) + bc(i,j-1).
505  For (i+j) much larger than 60, this is inefficient.
506  If you need binomial coefficients for large i and j, then you
507  should probably be using something like Stirling's Formula.
508  (Look up "Stirling" or "Gamma function" in a calculus book.)
509 */
510 ON_DECL
511 double ON_BinomialCoefficient(
512  int i,
513  int j
514  );
515 
516 
517 /*
518 Description:
519  Compute a trinomial coefficient.
520 Parameters:
521  i - [in]
522  j - [in]
523  k - [in]
524 Returns:
525  (i+j+k)!/(i!j!k!), if 0 <= i, 0 <= j and 0<= k, and 0 otherwise.
526 See Also:
527  ON_BinomialCoefficient()
528 Remarks:
529  The trinomial coefficient is computed using the formula
530 
531  (i+j+k)! (i+j+k)! (j+k)!
532  -------- = -------- * -------
533  i! j! k! i! (j+k)! j! k!
534 
535  = ON_BinomialCoefficient(i,j+k)*ON_BinomialCoefficient(j,k)
536 
537 */
538 ON_DECL
539 double ON_TrinomialCoefficient(
540  int i,
541  int j,
542  int k
543  );
544 
545 
546 ON_DECL
547 ON_BOOL32 ON_GetParameterTolerance(
548  double, double, // domain
549  double, // parameter in domain
550  double*, double* // parameter tolerance (tminus, tplus) returned here
551  );
552 
553 
554 ON_DECL
555 ON_BOOL32 ON_IsValidPointList(
556  int, // dim
557  ON_BOOL32, // true for homogeneous rational points
558  int, // count
559  int, // stride
560  const float*
561  );
562 
563 ON_DECL
564 ON_BOOL32 ON_IsValidPointList(
565  int, // dim
566  ON_BOOL32, // true for homogeneous rational points
567  int, // count
568  int, // stride
569  const double*
570  );
571 
572 /*
573 Description:
574  Determine if a list of points is planar.
575 Parameters:
576  bRational - [in]
577  false if the points are euclidean (x,y,z)
578  true if the points are homogeneous rational (x,y,z,w)
579  point_count - [in]
580  number of points
581  point_stride - [in]
582  number of doubles between point x coordinates
583  first point's x coordinate = points[0],
584  second point's x coordinate = points[point_stride],...
585  points - [in]
586  point coordinates (3d or 4d homogeneous rational)
587  boxMin - [in]
588  boxMax - [in]
589  optional 3d bounding box - pass nulls if not readily available
590  tolerance - [in] >= 0.0
591  plane_equation0 - [in]
592  If you want to test for planarity in a specific plane,
593  pass the plane equation in here. If you want to find
594  a plane containing the points, pass null here.
595  plane_equation - [out]
596  If this point is not null, then the equation of the plane
597  containing the points is retuened here.
598 Returns:
599  0 - points are not coplanar to the specified tolerance
600  1 - points are coplanar to the specified tolerance
601  2 - points are colinear to the specified tolerance
602  (in this case, plane_equation is not a unique answer)
603  3 - points are coincident to the specified tolerance
604  (in this case, plane_equation is not a unique answer)
605 */
606 ON_DECL
607 int ON_IsPointListPlanar(
608  bool bRational,
609  int count,
610  int stride,
611  const double* points,
612  const double* boxMin,
613  const double* boxMax,
614  double tolerance,
615  ON_PlaneEquation* plane_equation
616  );
617 
618 ON_DECL
619 ON_BOOL32 ON_IsValidPointGrid(
620  int, // dim
621  ON_BOOL32, // true for homogeneous rational points
622  int, int, // point_count0, point_count1,
623  int, int, // point_stride0, point_stride1,
624  const double*
625  );
626 
627 ON_DECL
628 bool ON_ReversePointList(
629  int, // dim
630  ON_BOOL32, // true for homogeneous rational points
631  int, // count
632  int, // stride
633  double*
634  );
635 
636 ON_DECL
637 ON_BOOL32 ON_ReversePointGrid(
638  int, // dim
639  ON_BOOL32, // true for homogeneous rational points
640  int, int, // point_count0, point_count1,
641  int, int, // point_stride0, point_stride1,
642  double*,
643  int // dir = 0 or 1
644  );
645 
646 ON_DECL
647 bool ON_SwapPointListCoordinates(
648  int, // count
649  int, // stride
650  float*,
651  int, int // coordinates to swap
652  );
653 
654 ON_DECL
655 bool ON_SwapPointListCoordinates(
656  int, // count
657  int, // stride
658  double*,
659  int, int // coordinates to swap
660  );
661 
662 ON_DECL
663 ON_BOOL32 ON_SwapPointGridCoordinates(
664  int, int, // point_count0, point_count1,
665  int, int, // point_stride0, point_stride1,
666  double*,
667  int, int // coordinates to swap
668  );
669 
670 ON_DECL
671 bool ON_TransformPointList(
672  int, // dim
673  ON_BOOL32, // true for homogeneous rational points
674  int, // count
675  int, // stride
676  float*,
677  const ON_Xform&
678  );
679 
680 ON_DECL
681 bool ON_TransformPointList(
682  int, // dim
683  ON_BOOL32, // true for homogeneous rational points
684  int, // count
685  int, // stride
686  double*,
687  const ON_Xform&
688  );
689 
690 ON_DECL
691 ON_BOOL32 ON_TransformPointGrid(
692  int, // dim
693  ON_BOOL32, // true for homogeneous rational points
694  int, int, // point_count0, point_count1,
695  int, int, // point_stride0, point_stride1,
696  double*,
697  const ON_Xform&
698  );
699 
700 ON_DECL
701 ON_BOOL32 ON_TransformVectorList(
702  int, // dim
703  int, // count
704  int, // stride
705  float*,
706  const ON_Xform&
707  );
708 
709 ON_DECL
710 ON_BOOL32 ON_TransformVectorList(
711  int, // dim
712  int, // count
713  int, // stride
714  double*,
715  const ON_Xform&
716  );
717 
718 /*
719 Parameters:
720  dim - [in]
721  >= 1
722  is_rat - [in]
723  true if the points are rational and points[dim] is the "weight"
724  pointA - [in]
725  pointB - [in]
726  point coordinates
727 Returns:
728  True if the input is valid and for each coordinate pair,
729  |a-b| <= ON_ZERO_TOLERANCE
730  or |a-b| <= (fabs(a)+fabs(b))*ON_RELATIVE_TOLERANCE.
731  False otherwise.
732 */
733 ON_DECL
734 bool ON_PointsAreCoincident(
735  int dim,
736  int is_rat,
737  const double* pointA,
738  const double* pointB
739  );
740 
741 /*
742 Description
743  See ON_PointsAreCoincident() for a description of when opennurbs
744  considers two points to be conincident.
745 Parameters:
746  dim - [in]
747  >= 1
748  is_rat - [in]
749  true if the points are rational and points[dim] is the "weight"
750  point_count - [in]
751  number of points >= 2
752  point_stride - [in]
753  >= (0 != is_rat) ? (dim+1) : dim
754  points - [in]
755  point coordinates
756 Returns:
757  True if the first and last points are coincident and all other
758  points in the list are coincident with the previous point.
759  False if there are points that are not coincident or
760  point_count < 2 or other input parameters are invalid.
761 */
762 ON_DECL
763 bool ON_PointsAreCoincident(
764  int dim,
765  int is_rat,
766  int point_count,
767  int point_stride,
768  const double* points
769  );
770 
771 ON_DECL
772 int ON_ComparePoint( // returns
773  // -1: first < second
774  // 0: first == second
775  // +1: first > second
776  int dim, // dim (>=0)
777  ON_BOOL32 israt, // true for rational CVs
778  const double* cv0, // first CV
779  const double* cv1 // secont CV
780  );
781 
782 ON_DECL
783 int ON_ComparePointList( // returns
784  // -1: first < second
785  // 0: first == second
786  // +1: first > second
787  int, // dim (>=0)
788  ON_BOOL32, // true for rational CVs
789  int, // count
790  // first point list
791  int, // stride
792  const double*, // point
793  // second point list
794  int, // stride
795  const double* // point
796  );
797 
798 ON_DECL
799 ON_BOOL32 ON_IsPointListClosed(
800  int, // dim
801  int, // true for homogeneos rational points
802  int, // count
803  int, // stride
804  const double*
805  );
806 
807 ON_DECL
808 ON_BOOL32 ON_IsPointGridClosed(
809  int, // dim
810  ON_BOOL32, // true for homogeneous rational points
811  int, int, // point_count0, point_count1,
812  int, int, // point_stride0, point_stride1,
813  const double*,
814  int // dir = 0 or 1
815  );
816 
817 ON_DECL
818 int ON_SolveQuadraticEquation( // solve a*X^2 + b*X + c = 0
819  // returns 0: two distinct real roots (r0 < r1)
820  // 1: one real root (r0 = r1)
821  // 2: two complex conjugate roots (r0 +/- (r1)*sqrt(-1))
822  // -1: failure - a = 0, b != 0 (r0 = r1 = -c/b)
823  // -2: failure - a = 0, b = 0 c != 0 (r0 = r1 = 0.0)
824  // -3: failure - a = 0, b = 0 c = 0 (r0 = r1 = 0.0)
825  double, double, double, // a, b, c
826  double*, double* // roots r0 and r1 returned here
827  );
828 
829 ON_DECL
830 ON_BOOL32 ON_SolveTriDiagonal( // solve TriDiagMatrix( a,b,c )*X = d
831  int, // dimension of d and X (>=1)
832  int, // number of equations (>=2)
833  double*, // a[n-1] = sub-diagonal (a is modified)
834  const double*, // b[n] = diagonal
835  double*, // c[n-1] = supra-diagonal
836  const double*, // d[n*dim]
837  double* // X[n*dim] = unknowns
838  );
839 
840 // returns rank - if rank != 2, system is under determined
841 // If rank = 2, then solution to
842 //
843 // a00*x0 + a01*x1 = b0,
844 // a10*x0 + a11*x1 = b1
845 //
846 // is returned
847 ON_DECL
848 int ON_Solve2x2(
849  double, double, // a00 a01 = first row of 2x2 matrix
850  double, double, // a10 a11 = second row of 2x2 matrix
851  double, double, // b0 b1
852  double*, double*, // x0, x1 if not NULL, then solution is returned here
853  double* // if not NULL, then pivot_ratio returned here
854  );
855 
856 // Description:
857 // Solves a system of 3 linear equations and 2 unknowns.
858 //
859 // x*col0[0] + y*col1[0] = d0
860 // x*col0[1] + y*col1[1] = d0
861 // x*col0[2] + y*col1[2] = d0
862 //
863 // Parameters:
864 // col0 - [in] coefficents for "x" unknown
865 // col1 - [in] coefficents for "y" unknown
866 // d0 - [in] constants
867 // d1 - [in]
868 // d2 - [in]
869 // x - [out]
870 // y - [out]
871 // error - [out]
872 // pivot_ratio - [out]
873 //
874 // Returns:
875 // rank of the system.
876 // If rank != 2, system is under determined
877 // If rank = 2, then the solution is
878 //
879 // (*x)*[col0] + (*y)*[col1]
880 // + (*error)*((col0 X col1)/|col0 X col1|)
881 // = (d0,d1,d2).
882 ON_DECL
883 int ON_Solve3x2(
884  const double[3], // col0
885  const double[3], // col1
886  double, // d0
887  double, // d1
888  double, // d2
889  double*, // x
890  double*, // y
891  double*, // error
892  double* // pivot_ratio
893  );
894 
895 /*
896 Description:
897  Use Gauss-Jordan elimination with full pivoting to solve
898  a system of 3 linear equations and 3 unknowns(x,y,z)
899 
900  x*row0[0] + y*row0[1] + z*row0[2] = d0
901  x*row1[0] + y*row1[1] + z*row1[2] = d1
902  x*row2[0] + y*row2[1] + z*row2[2] = d2
903 
904 Parameters:
905  row0 - [in] first row of 3x3 matrix
906  row1 - [in] second row of 3x3 matrix
907  row2 - [in] third row of 3x3 matrix
908  d0 - [in]
909  d1 - [in]
910  d2 - [in] (d0,d1,d2) right hand column of system
911  x_addr - [in] first unknown
912  y_addr - [in] second unknown
913  z_addr - [in] third unknown
914  pivot_ratio - [out] if not NULL, the pivot ration is
915  returned here. If the pivot ratio is "small",
916  then the matrix may be singular or ill
917  conditioned. You should test the results
918  before you use them. "Small" depends on the
919  precision of the input coefficients and the
920  use of the solution. If you can't figure out
921  what "small" means in your case, then you
922  must check the solution before you use it.
923 
924 Returns:
925  The rank of the 3x3 matrix (0,1,2, or 3)
926  If ON_Solve3x3() is successful (returns 3), then
927  the solution is returned in
928  (*x_addr, *y_addr, *z_addr)
929  and *pivot_ratio = min(|pivots|)/max(|pivots|).
930  If the return code is < 3, then (0,0,0) is returned
931  as the "solution".
932 
933 See Also:
934  ON_Solve2x2
935  ON_Solve3x2
936  ON_Solve4x4
937 */
938 ON_DECL
939 int ON_Solve3x3(
940  const double row0[3],
941  const double row1[3],
942  const double row2[3],
943  double d0,
944  double d1,
945  double d2,
946  double* x_addr,
947  double* y_addr,
948  double* z_addr,
949  double* pivot_ratio
950  );
951 
952 /*
953 Description:
954  Use Gauss-Jordan elimination with full pivoting to solve
955  a system of 4 linear equations and 4 unknowns(x,y,z,w)
956 
957  x*row0[0] + y*row0[1] + z*row0[2] + w*row0[3] = d0
958  x*row1[0] + y*row1[1] + z*row1[2] + w*row1[3] = d1
959  x*row2[0] + y*row2[1] + z*row2[2] + w*row2[3] = d2
960  x*row3[0] + y*row3[1] + z*row3[2] + w*row3[2] = d3
961 
962 Parameters:
963  row0 - [in] first row of 4x4 matrix
964  row1 - [in] second row of 4x4 matrix
965  row2 - [in] third row of 4x4 matrix
966  row3 - [in] forth row of 4x4 matrix
967  d0 - [in]
968  d1 - [in]
969  d2 - [in]
970  d3 - [in] (d0,d1,d2,d3) right hand column of system
971  x_addr - [in] first unknown
972  y_addr - [in] second unknown
973  z_addr - [in] third unknown
974  w_addr - [in] forth unknown
975  pivot_ratio - [out] if not NULL, the pivot ration is
976  returned here. If the pivot ratio is "small",
977  then the matrix may be singular or ill
978  conditioned. You should test the results
979  before you use them. "Small" depends on the
980  precision of the input coefficients and the
981  use of the solution. If you can't figure out
982  what "small" means in your case, then you
983  must check the solution before you use it.
984 
985 Returns:
986  The rank of the 4x4 matrix (0,1,2,3, or 4)
987  If ON_Solve4x4() is successful (returns 4), then
988  the solution is returned in
989  (*x_addr, *y_addr, *z_addr, *w_addr)
990  and *pivot_ratio = min(|pivots|)/max(|pivots|).
991  If the return code is < 4, then, it a solution exists,
992  on is returned. However YOU MUST CHECK THE SOLUTION
993  IF THE RETURN CODE IS < 4.
994 
995 See Also:
996  ON_Solve2x2
997  ON_Solve3x2
998  ON_Solve3x3
999 */
1000 ON_DECL
1001 int
1002 ON_Solve4x4(
1003  const double row0[4],
1004  const double row1[4],
1005  const double row2[4],
1006  const double row3[4],
1007  double d0,
1008  double d1,
1009  double d2,
1010  double d3,
1011  double* x_addr,
1012  double* y_addr,
1013  double* z_addr,
1014  double* w_addr,
1015  double* pivot_ratio
1016  );
1017 
1018 /*
1019 Description:
1020  Use Gauss-Jordan elimination to find a numerical
1021  solution to M*X = B where M is a n x n matrix,
1022  B is a known n-dimensional vector and X is
1023  an unknown.
1024 Paramters:
1025  bFullPivot - [in] if true, full pivoting is used,
1026  otherwise partial pivoting is used. In rare
1027  cases full pivoting can produce a more accurate
1028  answer and never produces a less accurate answer.
1029  However full pivoting is slower. If speed is an
1030  issue, then experiement with bFullPivot=false
1031  and see if it makes a difference. Otherwise,
1032  set it to true.
1033  bNormalize - [in]
1034  If bNormalize is true, then the rows of the
1035  matrix are scaled so the sum of their squares
1036  is one. This doesn't make the solution more
1037  accurate but in some cases it makes the pivot
1038  ratio more meaningful. Set bNormalize to
1039  false unless you have a reason for setting it
1040  to true.
1041  n - [in] size of the matrix and vectors.
1042  M - [in] n x n matrix. The values in M are
1043  changed as the solution is calculated.
1044  If you need to preserve M for future use,
1045  pass in a copy.
1046  B - [in] n-dimensional vector. The values in
1047  B are changed as the solution is calculated.
1048  If you need to preserve B for future use,
1049  pass in a copy.
1050  X - [out] solution to M*X = B.
1051 Returns:
1052  If the returned value is <= 0.0, the input matrix
1053  has rank < n and no solution is returned in X.
1054  If the returned value is > 0.0, then a solution is
1055  returned in X and the returned value is the ratio
1056  (minimum pivot)/(maximum pivot). This value is
1057  called the pivot ratio and will be denoted "pr"
1058  the discussion below. If pr <= 1e-15, then
1059  M was nearly degenerate and the solution should be
1060  used with caution. If an accurate solution is
1061  critcial, then check the solution anytime pr <= 1e-10
1062  In general, the difference between M*X and B will be
1063  reasonably small. However, when the pr is small
1064  there tend to be vector E, substantually different
1065  from zero, such that M*(X+E) - B is also reasonably
1066  small.
1067 See Also:
1068  ON_Solve2x2
1069  ON_Solve3x3
1070  ON_Solve4x4
1071  ON_Solve3x2
1072 */
1073 ON_DECL
1074 double ON_SolveNxN(bool bFullPivot, bool bNormalize, int n, double* M[], double B[], double X[]);
1075 
1076 
1077 // return false if determinant is (nearly) singular
1078 ON_DECL
1079 ON_BOOL32 ON_EvJacobian(
1080  double, // ds o ds
1081  double, // ds o dt
1082  double, // dt o dt
1083  double* // jacobian = determinant ( ds_o_ds dt_o_dt / ds_o_dt ds_o_dt )
1084  );
1085 
1086 /*
1087 Description:
1088  Finds scalars x and y so that the component of V in the plane
1089  of A and B is x*A + y*B.
1090 Parameters:
1091  V - [in]
1092  A - [in] nonzero and not parallel to B
1093  B - [in] nonzero and not parallel to A
1094  x - [out]
1095  y - [out]
1096 Returns:
1097  1 - The rank of the problem is 2. The decomposition is unique.
1098  0 - The rank less than 2. Either there is no solution or there
1099  are infinitely many solutions.
1100 
1101 See Also:
1102  ON_Solve2x2
1103 */
1104 ON_DECL
1105 int ON_DecomposeVector(
1106  const ON_3dVector& V,
1107  const ON_3dVector& A,
1108  const ON_3dVector& B,
1109  double* x, double* y
1110  );
1111 
1112 
1113 /*
1114 Description:
1115  Evaluate partial derivatives of surface unit normal
1116 Parameters:
1117  ds - [in]
1118  dt - [in] surface first partial derivatives
1119  dss - [in]
1120  dst - [in]
1121  dtt - [in] surface second partial derivatives
1122  ns - [out]
1123  nt - [out] First partial derivatives of surface unit normal
1124  (If the Jacobian is degenerate, ns and nt are set to zero.)
1125 Returns:
1126  true if Jacobian is nondegenerate
1127  false if Jacobian is degenerate
1128 */
1129 ON_DECL
1130 ON_BOOL32 ON_EvNormalPartials(
1131  const ON_3dVector& ds,
1132  const ON_3dVector& dt,
1133  const ON_3dVector& dss,
1134  const ON_3dVector& dst,
1135  const ON_3dVector& dtt,
1136  ON_3dVector& ns,
1137  ON_3dVector& nt
1138  );
1139 
1140 ON_DECL
1141 ON_BOOL32
1142 ON_Pullback3dVector( // use to pull 3d vector back to surface parameter space
1143  const ON_3dVector&, // 3d vector
1144  double, // signed distance from vector location to closet point on surface
1145  // < 0 if point is below with respect to Du x Dv
1146  const ON_3dVector&, // ds surface first partials
1147  const ON_3dVector&, // dt
1148  const ON_3dVector&, // dss surface 2nd partials
1149  const ON_3dVector&, // dst (used only when dist != 0)
1150  const ON_3dVector&, // dtt
1151  ON_2dVector& // pullback
1152  );
1153 
1154 ON_DECL
1155 ON_BOOL32
1156 ON_GetParameterTolerance(
1157  double, // t0 domain
1158  double, // t1
1159  double, // t parameter in domain
1160  double*, // tminus parameter tolerance (tminus, tplus) returned here
1161  double* // tplus
1162  );
1163 
1164 
1165 ON_DECL
1166 ON_BOOL32 ON_EvNormal(
1167  int, // limit_dir 0=default,1=from quadrant I, 2 = from quadrant II, ...
1168  const ON_3dVector&, const ON_3dVector&, // first partials (Du,Dv)
1169  const ON_3dVector&, const ON_3dVector&, const ON_3dVector&, // optional second partials (Duu, Duv, Dvv)
1170  ON_3dVector& // unit normal returned here
1171  );
1172 
1173 // returns false if first returned tangent is zero
1174 ON_DECL
1175 bool ON_EvTangent(
1176  const ON_3dVector&, // first derivative
1177  const ON_3dVector&, // second derivative
1178  ON_3dVector& // Unit tangent returned here
1179  );
1180 
1181 // returns false if first derivtive is zero
1182 ON_DECL
1183 ON_BOOL32 ON_EvCurvature(
1184  const ON_3dVector&, // first derivative
1185  const ON_3dVector&, // second derivative
1186  ON_3dVector&, // Unit tangent returned here
1187  ON_3dVector& // Curvature returned here
1188  );
1189 
1190 ON_DECL
1191 ON_BOOL32 ON_EvPrincipalCurvatures(
1192  const ON_3dVector&, // Ds,
1193  const ON_3dVector&, // Dt,
1194  const ON_3dVector&, // Dss,
1195  const ON_3dVector&, // Dst,
1196  const ON_3dVector&, // Dtt,
1197  const ON_3dVector&, // N, // unit normal to surface (use ON_EvNormal())
1198  double*, // gauss, // = Gaussian curvature = kappa1*kappa2
1199  double*, // mean, // = mean curvature = (kappa1+kappa2)/2
1200  double*, // kappa1, // = largest principal curvature value (may be negative)
1201  double*, // kappa2, // = smallest principal curvature value (may be negative)
1202  ON_3dVector&, // K1, // kappa1 unit principal curvature direction
1203  ON_3dVector& // K2 // kappa2 unit principal curvature direction
1204  // output K1,K2,N is right handed frame
1205  );
1206 
1207 ON_DECL
1208 ON_BOOL32 ON_EvPrincipalCurvatures(
1209  const ON_3dVector&, // Ds,
1210  const ON_3dVector&, // Dt,
1211  double l, // Dss*N Second fundamental form coefficients
1212  double m, // Dst*N,
1213  double n, // Dtt*N,
1214  const ON_3dVector&, // N, // unit normal to surface (use ON_EvNormal())
1215  double*, // gauss, // = Gaussian curvature = kappa1*kappa2
1216  double*, // mean, // = mean curvature = (kappa1+kappa2)/2
1217  double*, // kappa1, // = largest principal curvature value (may be negative)
1218  double*, // kappa2, // = smallest principal curvature value (may be negative)
1219  ON_3dVector&, // K1, // kappa1 unit principal curvature direction
1220  ON_3dVector& // K2 // kappa2 unit principal curvature direction
1221  // output K1,K2,N is right handed frame
1222  );
1223 
1224 /*
1225 Description:
1226  Evaluate sectional curvature from surface derivatives and
1227  section plane normal.
1228 Parameters:
1229  S10, S01 - [in]
1230  surface 1st partial derivatives
1231  S20, S11, S02 - [in]
1232  surface 2nd partial derivatives
1233  planeNormal - [in]
1234  unit normal to section plane
1235  K - [out] Sectional curvature
1236  Curvature of the intersection curve of the surface
1237  and plane through the surface point where the partial
1238  derivatives were evaluationed.
1239 Returns:
1240  True if successful.
1241  False if first partials are not linearly independent, in
1242  which case the K is set to zero.
1243 */
1244 ON_DECL
1245 bool ON_EvSectionalCurvature(
1246  const ON_3dVector& S10,
1247  const ON_3dVector& S01,
1248  const ON_3dVector& S20,
1249  const ON_3dVector& S11,
1250  const ON_3dVector& S02,
1251  const ON_3dVector& planeNormal,
1252  ON_3dVector& K
1253  );
1254 
1255 
1256 ON_DECL
1257 ON_3dVector ON_NormalCurvature(
1258  const ON_3dVector&, // surface 1rst partial (Ds)
1259  const ON_3dVector&, // surface 1rst partial (Dt)
1260  const ON_3dVector&, // surface 1rst partial (Dss)
1261  const ON_3dVector&, // surface 1rst partial (Dst)
1262  const ON_3dVector&, // surface 1rst partial (Dtt)
1263  const ON_3dVector&, // surface unit normal
1264  const ON_3dVector& // unit tangent direction
1265  );
1266 
1267 /*
1268 Description:
1269  Determing if two curvatrues are different enough
1270  to qualify as a curvature discontinuity.
1271 Parameters:
1272  Km - [in]
1273  Kp - [in]
1274  Km and Kp should be curvatures evaluated at the same
1275  parameters using limits from below (minus) and above (plus).
1276  The assumption is that you have already compared the
1277  points and tangents and consider to curve to be G1 at the
1278  point in question.
1279  cos_angle_tolerance - [in]
1280  If the inut value of cos_angle_tolerance >= -1.0
1281  and cos_angle_tolerance <= 1.0 and
1282  Km o Kp < cos_angle_tolerance*|Km|*|Kp|, then
1283  true is returned. Otherwise it is assumed Km and Kp
1284  are parallel. If the curve being tested is nonplanar,
1285  then use something like cos(2*tangent angle tolerance)
1286  for this parameter. If the curve being tested is planar,
1287  then 0.0 will work fine.
1288  curvature_tolerance - [in]
1289  If |Kp-Km| <= curvature_tolerance,
1290  then false is returned, otherwise other tests are used
1291  to determing continuity.
1292  zero_curvature - [in] (ignored if < 2^-110 = 7.7037197787136e-34)
1293  If |K| <= zero_curvature, then K is treated as zero.
1294  When in doubt, use ON_ZERO_CURVATURE_TOLERANCE.
1295  radius_tolerance - [in]
1296  If radius_tolerance >= 0.0 and the difference between the
1297  radii of curvature is >= radius_tolerance, then true
1298  is returned.
1299  relative_tolerance - [in]
1300  If relative_tolerance > 0 and
1301  |(|Km| - |Kp|)|/max(|Km|,|Kp|) > relative_tolerance,
1302  then true is returned. Note that if the curvatures are
1303  nonzero and rm and rp are the radii of curvature, then
1304  |(|Km| - |Kp|)|/max(|Km|,|Kp|) = |rm-rp|/max(rm,rp).
1305  This means the relative_tolerance insures both the scalar
1306  curvature and the radii of curvature agree to the specified
1307  number of decimal places.
1308  When in double use ON_RELATIVE_CURVATURE_TOLERANCE, which
1309  is currently 0.05.
1310 Returns:
1311  False if the curvatures should be considered G2.
1312  True if the curvatures are different enough that the curve should be
1313  considered not G2.
1314  In addition to the tests described under the curvature_tolerance and
1315  radius_tolerance checks, other hurestic tests are used.
1316 */
1317 ON_DECL
1318 bool ON_IsCurvatureDiscontinuity(
1319  const ON_3dVector Km,
1320  const ON_3dVector Kp,
1321  double cos_angle_tolerance,
1322  double curvature_tolerance,
1323  double zero_curvature,
1324  double radius_tolerance,
1325  double relative_tolerance
1326  );
1327 
1328 ON_DECL
1329 bool ON_IsCurvatureDiscontinuity(
1330  const ON_3dVector Km,
1331  const ON_3dVector Kp,
1332  double cos_angle_tolerance,
1333  double curvature_tolerance,
1334  double zero_curvature,
1335  double radius_tolerance
1336  );
1337 
1338 
1339 /*
1340 Description:
1341  This function is used to test curvature continuity
1342  in IsContinuous and GetNextDiscontinuity functions
1343  when the continuity parameter is ON::G2_continuous.
1344 Parameters:
1345  Km - [in]
1346  Curve's vector curvature evaluated from below
1347  Kp - [in]
1348  Curve's vector curvature evaluated from below
1349 Returns:
1350  True if the change from Km to Kp should be considered
1351  G2 continuous.
1352 */
1353 ON_DECL
1354 bool ON_IsG2CurvatureContinuous(
1355  const ON_3dVector Km,
1356  const ON_3dVector Kp,
1357  double cos_angle_tolerance,
1358  double curvature_tolerance
1359  );
1360 
1361 /*
1362 Description:
1363  This function is used to test curvature continuity
1364  in IsContinuous and GetNextDiscontinuity functions
1365  when the continuity parameter is ON::Gsmooth_continuous.
1366 Parameters:
1367  Km - [in]
1368  Curve's vector curvature evaluated from below
1369  Kp - [in]
1370  Curve's vector curvature evaluated from below
1371 Returns:
1372  True if the change from Km to Kp should be considered
1373  Gsmooth continuous.
1374 */
1375 ON_DECL
1376 bool ON_IsGsmoothCurvatureContinuous(
1377  const ON_3dVector Km,
1378  const ON_3dVector Kp,
1379  double cos_angle_tolerance,
1380  double curvature_tolerance
1381  );
1382 
1383 /*
1384 Description:
1385  Test curve continuity from derivative values.
1386 Parameters:
1387  c - [in] type of continuity to test for. Read ON::continuity
1388  comments for details.
1389  Pa - [in] point on curve A.
1390  D1a - [in] first derviative of curve A.
1391  D2a - [in] second derviative of curve A.
1392  Pb - [in] point on curve B.
1393  D1b - [in] first derviative of curve B.
1394  D3b - [in] second derviative of curve B.
1395  point_tolerance - [in] if the distance between two points is
1396  greater than point_tolerance, then the curve is not C0.
1397  d1_tolerance - [in] if the difference between two first derivatives is
1398  greater than d1_tolerance, then the curve is not C1.
1399  d2_tolerance - [in] if the difference between two second derivatives is
1400  greater than d2_tolerance, then the curve is not C2.
1401  cos_angle_tolerance - [in] default = cos(1 degree) Used only when
1402  c is ON::G1_continuous or ON::G2_continuous. If the cosine
1403  of the angle between two tangent vectors
1404  is <= cos_angle_tolerance, then a G1 discontinuity is reported.
1405  curvature_tolerance - [in] (default = ON_SQRT_EPSILON) Used only when
1406  c is ON::G2_continuous. If K0 and K1 are curvatures evaluated
1407  from above and below and |K0 - K1| > curvature_tolerance,
1408  then a curvature discontinuity is reported.
1409 Returns:
1410  true if the curve has at least the c type continuity at
1411  the parameter t.
1412 */
1413 ON_DECL
1414 ON_BOOL32 ON_IsContinuous(
1415  ON::continuity c,
1416  ON_3dPoint Pa,
1417  ON_3dVector D1a,
1418  ON_3dVector D2a,
1419  ON_3dPoint Pb,
1420  ON_3dVector D1b,
1421  ON_3dVector D2b,
1422  double point_tolerance=ON_ZERO_TOLERANCE,
1423  double d1_tolerance=ON_ZERO_TOLERANCE,
1424  double d2_tolerance=ON_ZERO_TOLERANCE,
1425  double cos_angle_tolerance=ON_DEFAULT_ANGLE_TOLERANCE_COSINE,
1426  double curvature_tolerance=ON_SQRT_EPSILON
1427  );
1428 
1429 
1430 ON_DECL
1431 bool ON_TuneupEvaluationParameter(
1432  int side,
1433  double s0, double s1, // segment domain
1434  double *s // segment parameter
1435  );
1436 
1437 
1438 ON_DECL
1439 int ON_Compare2dex( const ON_2dex* a, const ON_2dex* b);
1440 
1441 ON_DECL
1442 int ON_Compare3dex( const ON_3dex* a, const ON_3dex* b);
1443 
1444 ON_DECL
1445 int ON_Compare4dex( const ON_4dex* a, const ON_4dex* b);
1446 
1447 ON_DECL
1448 const ON_2dex* ON_BinarySearch2dexArray(
1449  int key_i,
1450  const ON_2dex* base,
1451  size_t nel
1452  );
1453 
1454 // These simple intersectors are fast and detect transverse intersections.
1455 // If the intersection is not a simple transverse case, then they
1456 // return false and you will have to use one of the slower but fancier
1457 // models.
1458 
1459 // returns closest points between the two infinite lines
1460 ON_DECL
1461 bool ON_Intersect(
1462  const ON_Line&,
1463  const ON_Line&,
1464  double*, // parameter on first line
1465  double* // parameter on second line
1466  );
1467 
1468 // Returns false unless intersection is a single point
1469 // If returned parameter is < 0 or > 1, then the line
1470 // segment between line.m_point[0] and line.m_point[1]
1471 // does not intersect the plane
1472 ON_DECL
1473 bool ON_Intersect(
1474  const ON_Line&,
1475  const ON_Plane&,
1476  double* // parameter on line
1477  );
1478 
1479 ON_DECL
1480 bool ON_Intersect(
1481  const ON_Plane&,
1482  const ON_Plane&,
1483  ON_Line& // intersection line is returned here
1484  );
1485 
1486 ON_DECL
1487 bool ON_Intersect(
1488  const ON_Plane&,
1489  const ON_Plane&,
1490  const ON_Plane&,
1491  ON_3dPoint& // intersection point is returned here
1492  );
1493 
1494 // returns 0 = no intersections,
1495 // 1 = intersection = single point,
1496 // 2 = intersection = circle
1497 // If 0 is returned, returned circle has radius=0
1498 // and center = point on sphere closest to plane.
1499 // If 1 is returned, intersection is a single
1500 // point and returned circle has radius=0
1501 // and center = intersection point on sphere.
1502 ON_DECL
1503 int ON_Intersect(
1504  const ON_Plane&, const ON_Sphere&, ON_Circle&
1505  );
1506 
1507 // Intersects an infinte line and sphere and returns
1508 // 0 = no intersections,
1509 // 1 = one intersection,
1510 // 2 = 2 intersections
1511 // If 0 is returned, first point is point
1512 // on line closest to sphere and 2nd point is the point
1513 // on the sphere closest to the line.
1514 // If 1 is returned, first point is obtained by evaluating
1515 // the line and the second point is obtained by evaluating
1516 // the sphere.
1517 ON_DECL
1518 int ON_Intersect(
1519  const ON_Line&,
1520  const ON_Sphere&,
1521  ON_3dPoint&,
1522  ON_3dPoint& // intersection point(s) returned here
1523  );
1524 
1525 
1526 // Intersects an infinte line and cylinder and returns
1527 // 0 = no intersections,
1528 // 1 = one intersection,
1529 // 2 = 2 intersections
1530 // 3 = line lies on cylinder
1531 //
1532 // If 0 is returned, first point is point
1533 // on line closest to cylinder and 2nd point is the point
1534 // on the cylinder closest to the line.
1535 // If 1 is returned, first point is obtained by evaluating
1536 // the line and the second point is obtained by evaluating
1537 // the cylinder.
1538 //
1539 // The value of cylinder.IsFinite() determines if the
1540 // intersection is performed on the finite or infinite cylinder.
1541 ON_DECL
1542 int ON_Intersect(
1543  const ON_Line&, // [in]
1544  const ON_Cylinder&, // [in]
1545  ON_3dPoint&, // [out] first intersection point
1546  ON_3dPoint& // [out] second intersection point
1547  );
1548 
1549 // Description:
1550 // Intersect an infinte line and circle.
1551 // Parameters:
1552 // line - [in]
1553 // circle - [in]
1554 // line_t0 - [out] line parameter of first intersection point
1555 // circle_point0 - [out] first intersection point on circle
1556 // line_t1 - [out] line parameter of second intersection point
1557 // circle_point1 - [out] second intersection point on circle
1558 // Returns:
1559 // 0 No intersection
1560 // 1 One intersection at line.PointAt(*line_t0)
1561 // 2 Two intersections at line.PointAt(*line_t0)
1562 // and line.PointAt(*line_t1).
1563 ON_DECL
1564 int ON_Intersect(
1565  const ON_Line& line,
1566  const ON_Circle& circle,
1567  double* line_t0,
1568  ON_3dPoint& circle_point0,
1569  double* line_t1,
1570  ON_3dPoint& circle_point1
1571  );
1572 
1573 
1574 
1575 // Description:
1576 // Intersect a infinte line and arc.
1577 // Parameters:
1578 // line - [in]
1579 // arc - [in]
1580 // line_t0 - [out] line parameter of first intersection point
1581 // arc_point0 - [out] first intersection point on arc
1582 // line_t1 - [out] line parameter of second intersection point
1583 // arc_point1 - [out] second intersection point on arc
1584 // Returns:
1585 // 0 No intersection
1586 // 1 One intersection at line.PointAt(*line_t0)
1587 // 2 Two intersections at line.PointAt(*line_t0)
1588 // and line.PointAt(*line_t1).
1589 ON_DECL
1590 int ON_Intersect(
1591  const ON_Line& line,
1592  const ON_Arc& arc,
1593  double* line_t0,
1594  ON_3dPoint& arc_point0,
1595  double* line_t1,
1596  ON_3dPoint& arc_point1
1597  );
1598 
1599 // Description:
1600 // Intersect a plane and a circle.
1601 // Parameters:
1602 // plane - [in]
1603 // circle - [in]
1604 // point0 - [out] first intersection point
1605 // point1 - [out] second intersection point
1606 // Returns:
1607 // 0 No intersection
1608 // 1 One intersection at point0
1609 // 2 Two intersections at point0
1610 // and point1.
1611 // 3 Circle lies on plane
1612 ON_DECL
1613 int ON_Intersect(
1614  const ON_Plane& plane,
1615  const ON_Circle& circle,
1616  ON_3dPoint& point0,
1617  ON_3dPoint& point1
1618  );
1619 
1620 // Description:
1621 // Intersect a plane and an arc.
1622 // Parameters:
1623 // plane - [in]
1624 // arc - [in]
1625 // point0 - [out] first intersection point
1626 // point1 - [out] second intersection point
1627 // Returns:
1628 // 0 No intersection
1629 // 1 One intersection at point0
1630 // 2 Two intersections at point0
1631 // and point1.
1632 // 3 Arc lies on plane
1633 ON_DECL
1634 int ON_Intersect(
1635  const ON_Plane& plane,
1636  const ON_Arc& arc,
1637  ON_3dPoint& point0,
1638  ON_3dPoint& point1
1639  );
1640 
1641 
1642 // returns 0 = no, 1 = yes, 2 = points are coincident and on line
1643 ON_DECL
1644 int ON_ArePointsOnLine(
1645  int, // dimension of points
1646  int, // is_rat = true if homogeneous rational
1647  int, // count = number of points
1648  int, // stride ( >= is_rat?(dim+1) :dim)
1649  const double*, // point array
1650  const ON_BoundingBox&, // if needed, use ON_GetBoundingBox(dim,is_rat,count,stride,point)
1651  const ON_Line&,
1652  double // tolerance (if 0.0, a tolerance based on bounding box size is used)
1653  );
1654 
1655 // returns 0 = no, 1 = yes, 2 = points are coincident and on line
1656 ON_DECL
1657 int ON_ArePointsOnPlane(
1658  int, // dimension of points
1659  int, // is_rat = true if homogeneous rational
1660  int, // count = number of points
1661  int, // stride ( >= is_rat?(dim+1) :dim)
1662  const double*, // point array
1663  const ON_BoundingBox&, // if needed, use ON_GetBoundingBox(dim,is_rat,count,stride,point)
1664  const ON_Plane&,
1665  double // tolerance (if 0.0, a tolerance based on bounding box size is used)
1666  );
1667 
1668 /*
1669 Description:
1670  Use the quotient rule to compute derivatives of a one parameter
1671  rational function F(t) = X(t)/W(t), where W is a scalar
1672  and F and X are vectors of dimension dim.
1673 Parameters:
1674  dim - [in]
1675  der_count - [in] number of derivative (>=0)
1676  v_stride - [in] (>= dim+1)
1677  v - [in/out]
1678  v[] is an array of length (der_count+1)*v_stride.
1679  The input v[] array contains derivatives of the numerator and
1680  denominator functions in the order (X, W), (Xt, Wt), (Xtt, Wtt), ...
1681  In general, the (dim+1) coordinates of the d-th derivative
1682  are in (v[n],...,v[n+dim]) where n = d*v_stride.
1683  In the output v[] array the derivatives of X are replaced with
1684  the derivatives of F and the derivatives of W are divided by
1685  w = v[dim].
1686 Returns:
1687  True if input is valid; i.e., v[dim] != 0.
1688 See Also:
1689  ON_EvaluateQuotientRule2
1690  ON_EvaluateQuotientRule3
1691 */
1692 ON_DECL
1693 bool ON_EvaluateQuotientRule(
1694  int dim,
1695  int der_count,
1696  int v_stride,
1697  double *v
1698  );
1699 
1700 /*
1701 Description:
1702  Use the quotient rule to compute partial derivatives of a two parameter
1703  rational function F(s,t) = X(s,t)/W(s,t), where W is a scalar
1704  and F and X are vectors of dimension dim.
1705 Parameters:
1706  dim - [in]
1707  der_count - [in] number of derivative (>=0)
1708  v_stride - [in] (>= dim+1)
1709  v - [in/out]
1710  v[] is an array of length (der_count+2)*(der_count+1)*v_stride.
1711  The input array contains derivatives of the numerator and denominator
1712  functions in the order X, W, Xs, Ws, Xt, Wt, Xss, Wss, Xst, Wst, Xtt, Wtt, ...
1713  In general, the (i,j)-th derivatives are in the (dim+1) entries of v[]
1714  v[k], ..., answer[k+dim], where k = ((i+j)*(i+j+1)/2 + j)*v_stride.
1715  In the output v[] array the derivatives of X are replaced with
1716  the derivatives of F and the derivatives of W are divided by
1717  w = v[dim].
1718 Returns:
1719  True if input is valid; i.e., v[dim] != 0.
1720 See Also:
1721  ON_EvaluateQuotientRule
1722  ON_EvaluateQuotientRule3
1723 */
1724 ON_DECL
1725 bool ON_EvaluateQuotientRule2(
1726  int dim,
1727  int der_count,
1728  int v_stride,
1729  double *v
1730  );
1731 
1732 /*
1733 Description:
1734  Use the quotient rule to compute partial derivatives of a 3 parameter
1735  rational function F(r,s,t) = X(r,s,t)/W(r,s,t), where W is a scalar
1736  and F and X are vectors of dimension dim.
1737 Parameters:
1738  dim - [in]
1739  der_count - [in] number of derivative (>=0)
1740  v_stride - [in] (>= dim+1)
1741  v - [in/out]
1742  v[] is an array of length
1743  v_stride*(der_count+1)*(der_count+2)*(der_count+3)/6.
1744  The input v[] array contains derivatives of the numerator and
1745  denominator functions in the order (X, W), (Xr, Wr), (Xs, Ws),
1746  (Xt, Wt), (Xrr, Wrr), (Xrs, Wrs), (Xrt, Wrt), (Xss, Wss),
1747  (Xst, Wst), (Xtt, Wtt), ...
1748  In general, the (dim+1) coordinates of the derivative
1749  (Dr^i Ds^j Dt^k, i+j+k=d) are at v[n], ..., v[n+dim] where
1750  n = v_stride*( d*(d+1)*(d+2)/6 + (d-i)*(d-i+1)/2 + k ).
1751  In the output v[] array the derivatives of X are replaced with
1752  the derivatives of F and the derivatives of W are divided by
1753  w = v[dim].
1754 Returns:
1755  True if input is valid; i.e., v[dim] != 0.
1756 See Also:
1757  ON_EvaluateQuotientRule
1758  ON_EvaluateQuotientRule2
1759 */
1760 ON_DECL
1761 bool ON_EvaluateQuotientRule3(
1762  int dim,
1763  int der_count,
1764  int v_stride,
1765  double *v
1766  );
1767 
1768 ON_DECL
1769 bool ON_GetPolylineLength(
1770  int, // dimension of points
1771  ON_BOOL32, // bIsRational true if points are homogeneous rational
1772  int, // number of points
1773  int, // stride between points
1774  const double*, // points
1775  double* // length returned here
1776  );
1777 
1778 
1779 /*
1780 Description:
1781  Find the index of the point in the point_list that is closest to P.
1782 Parameters:
1783  point_count - [in]
1784  point_list - [in]
1785  P - [in]
1786  closest_point_index - [out]
1787 Returns:
1788  True if successful and *closest_point_index is set.
1789  False if input is not valid, in which case *closest_point_index
1790  is undefined.
1791 */
1792 ON_DECL
1793 bool ON_GetClosestPointInPointList(
1794  int point_count,
1795  const ON_3dPoint* point_list,
1796  ON_3dPoint P,
1797  int* closest_point_index
1798  );
1799 
1800 /*
1801 Description:
1802  Test math library functions.
1803 Parameters:
1804  function_index - [in] Determines which math library function is called.
1805 
1806  1: z = x+y
1807  2: z = x-y
1808  3: z = x*y
1809  4: z = x/y
1810  5: z = fabs(x)
1811  6: z = exp(x)
1812  7: z = log(x)
1813  8: z = logb(x)
1814  9: z = log10(x)
1815  10: z = pow(x,y)
1816  11: z = sqrt(x)
1817  12: z = sin(x)
1818  13: z = cos(x)
1819  14: z = tan(x)
1820  15: z = sinh(x)
1821  16: z = cosh(x)
1822  17: z = tanh(x)
1823  18: z = asin(x)
1824  19: z = acos(x)
1825  20: z = atan(x)
1826  21: z = atan2(y,x)
1827  22: z = fmod(x,y)
1828  23: z = modf(x,&y)
1829  24: z = frexp(x,&y)
1830 
1831  double x - [in]
1832  double y - [in]
1833 Returns:
1834  Returns the "z" value listed in the function_index parameter
1835  description.
1836 Remarks:
1837  This function is used to test the results of class floating
1838  point functions. It is primarily used to see what happens
1839  when opennurbs is used as a DLL and illegal operations are
1840  performed.
1841 */
1842 ON_DECL
1843 double ON_TestMathFunction(
1844  int function_index,
1845  double x,
1846  double y
1847  );
1848 
1849 // If performance is important, then
1850 // you are better off using ((b<a)?a:b)
1851 ON_DECL double ON_Max(double a, double b);
1852 
1853 // If performance is important, then
1854 // you are better off using ((b<a)?a:b)
1855 ON_DECL float ON_Max(float a, float b);
1856 
1857 // If performance is important, then
1858 // you are better off using ((b<a)?a:b)
1859 ON_DECL int ON_Max(int a, int b);
1860 
1861 // If performance is important, then
1862 // you are better off using ((a<b)?a:b)
1863 ON_DECL double ON_Min(double a, double b);
1864 
1865 // If performance is important, then
1866 // you are better off using ((a<b)?a:b)
1867 ON_DECL float ON_Min(float a, float b);
1868 
1869 // If performance is important, then
1870 // you are better off using ((a<b)?a:b)
1871 ON_DECL int ON_Min(int a, int b);
1872 
1873 // Do not call ON_Round() in any opennurbs code, tl code
1874 // or any other code that does critical calculations or
1875 // when there is any possibility that x is invalid or
1876 // fabs(x)>2147483647. Use floor(x+0.5) instead.
1877 ON_DECL int ON_Round(double x);
1878 
1879 
1880 /*
1881 Description:
1882  Find the equation of the parabola, ellipse or hyperbola
1883  (non-degenerate conic) that passes through six distinct points.
1884 Parameters:
1885  stride - [in] (>=2)
1886  points array stride
1887  points2d - [in] (>=2)
1888  i-th point is (points[i*stride],points[i*stride+1])
1889  conic - [out]
1890  Coefficients of the conic equation.
1891  The points on the conic satisfy the equation
1892  0 = conic[0]*x^2 + conic[1]*xy + conic[2]*y^2
1893  + conic[3]*x + conic[4]*y + conic[5]
1894  max_pivot - [out] (can be null)
1895  min_pivot - [out] (can be null)
1896  zero_pivot - [out] (can be null)
1897  If there are some near duplicates in the input point set,
1898  the calculation is not stable. If you want to get an
1899  estimate of the validity of the solution, then inspect
1900  the returned values. max_pivot should around 1,
1901  min_pivot should be > 1e-4 or so, and zero_pivot should
1902  be < 1e-10 or so. If the returned pivots don't satisify
1903  these condtions, then exercise caution when using the
1904  returned solution.
1905 Returns:
1906  True if a there is an ellipse, parabola or hyperbola through the
1907  six points.
1908  False if the input is invalid or the conic degenerate (the
1909  points lie on one or two lines).
1910  If false is returned, then conic[0]=...=conic[5] = 0 and
1911  *min_pivot = *max_pivot = *zero_pivot = 0.
1912 */
1913 ON_DECL bool ON_GetConicEquationThrough6Points(
1914  int stride,
1915  const double* points2d,
1916  double conic[6],
1917  double* max_pivot,
1918  double* min_pivot,
1919  double* zero_pivot
1920  );
1921 
1922 /*
1923 Description:
1924  Test a conic equation to see if it defines and ellipse. If so,
1925  return the center and axes of the ellipse.
1926 Parameters:
1927  conic - [in]
1928  Coefficients of the conic equation.
1929  The points on the conic satisfy the equation
1930  0 = conic[0]*x^2 + conic[1]*xy + conic[2]*y^2
1931  + conic[3]*x + conic[4]*y + conic[5]
1932  center - [out]
1933  major_axis - [out]
1934  minor_axis - [out]
1935  major_radius - [out]
1936  minor_radius - [out]
1937 Returns:
1938  True if the conic is an ellipse and the center and axes were found.
1939  False if the conic is not an ellipse, in which case the input values
1940  of center, major_axis, minor_axis, major_radius, and minor_radius
1941  are not changed.
1942 */
1943 ON_DECL bool ON_IsConicEquationAnEllipse(
1944  const double conic[6],
1945  ON_2dPoint& center,
1946  ON_2dVector& major_axis,
1947  ON_2dVector& minor_axis,
1948  double* major_radius,
1949  double* minor_radius
1950  );
1951 
1952 /*
1953 Description:
1954  Get the conic equation of an ellipse.
1955 Parameters:
1956  a - [in] (a>0)
1957  b - [in] (b>0)
1958  a and b are the lengths of the axes. Either one
1959  may be largest and they can be equal.
1960  x0 - [in]
1961  y0 - [in]
1962  (x0,y0) is the enter of the ellipse.
1963  alpha - [in] (angle in radians)
1964  When alpha is 0, a corresponds to the x-axis and
1965  b corresponds to the y axis. If alpha is non-zero
1966  it specifies the rotation of these axes.
1967  conic - [out]
1968  Coefficients of the conic equation.
1969  The points on the conic satisfy the equation
1970  0 = conic[0]*x^2 + conic[1]*xy + conic[2]*y^2
1971  + conic[3]*x + conic[4]*y + conic[5]
1972  center - [out]
1973  major_axis - [out]
1974  minor_axis - [out]
1975  major_radius - [out]
1976  minor_radius - [out]
1977 Remarks:
1978  Here is the way to evaluate a point on the ellipse:
1979 
1980 
1981  double t = ellipse paramter in radians;
1982  double x = a*cos(t);
1983  double y = b*sin(t);
1984  ON_2dPoint ellipse_point;
1985  ellipse_point.x = x0 + x*cos(alpha) + y*sin(alpha);
1986  ellipse_point.y = y0 - x*sin(alpha) + y*cos(alpha);
1987 
1988 Returns:
1989  True if the input is valid and conic[] was filled in.
1990  Falis if the input is not valid. In this case the values in conic[]
1991  are not changed.
1992 */
1993 ON_DECL bool ON_GetEllipseConicEquation(
1994  double a, double b,
1995  double x0, double y0,
1996  double alpha,
1997  double conic[6]
1998  );
1999 
2000 /*
2001 Descripton:
2002  Return the length of a 2d vector (x,y)
2003 Returns:
2004  sqrt(x^2 + y^2) calculated in as precisely and safely as possible.
2005 */
2006 ON_DECL double ON_Length2d( double x, double y );
2007 
2008 /*
2009 Descripton:
2010  Return the length of a 3d vector (x,y,z)
2011 Returns:
2012  sqrt(x^2 + y^2 + z^2) calculated in as precisely and safely as possible.
2013 */
2014 ON_DECL double ON_Length3d( double x, double y, double z );
2015 
2016 
2017 /*
2018 Description:
2019  Convert a double x to the largest float f such that
2020  the mathematical value of f is at most the value of x.
2021 Parameters:
2022  x - [in]
2023 Returns
2024  The largest float f such that the mathematical value
2025  of f is at most the value of x.
2026 */
2027 ON_DECL float ON_FloatFloor(double x);
2028 
2029 /*
2030 Description:
2031  Convert a double x to the smallest float f such that
2032  the mathematical value of f is at least the value of x.
2033 Parameters:
2034  x - [in]
2035 Returns
2036  The smallest float f such that the mathematical value
2037  of f is at least the value of x.
2038 */
2039 ON_DECL float ON_FloatCeil(double x);
2040 
2041 #endif
ON_SimpleArray< bool > m_bPeriodicParameter
const int m_parameter_count
ON_SimpleArray< ON_Interval > m_domain
Definition: norms.h:54
Definition: norms.h:54
const int m_value_count