com_croftsoft_core/math/math_lib/mod.rs
1// =============================================================================
2//! - A collection of math functions
3//!
4//! # Metadata
5//! - Copyright: © 2022-2023 [`CroftSoft Inc`]
6//! - Author: [`David Wallace Croft`]
7//! - Java created: 1998-12-27
8//! - Java updated: 2008-08-09
9//! - Rust created: 2022-08-24
10//! - Rust updated: 2023-03-20
11//!
12//! # History
13//! - Adapted from the Java class com.croftsoft.core.math.MathLib
14//! - In the Java-based [`CroftSoft Core Library`]
15//!
16//! [`CroftSoft Core Library`]: https://www.croftsoft.com/library/code/
17//! [`CroftSoft Inc`]: https://www.croftsoft.com/
18//! [`David Wallace Croft`]: https://www.croftsoft.com/people/david/
19// =============================================================================
20
21// -----------------------------------------------------------------------------
22/// Clips the value to a minimum and maximum range
23///
24/// # Alternative
25/// <https://doc.rust-lang.org/std/primitive.f64.html#method.clamp>
26///
27/// # Examples
28/// ```
29/// use com_croftsoft_core::math::math_lib::*;
30/// assert_eq!(
31/// Clip {
32/// maximum: 1.0,
33/// minimum: -1.0,
34/// value: 0.0,
35/// }.calculate().unwrap(),
36/// 0.0);
37/// assert_eq!(
38/// Clip {
39/// maximum: 1.0,
40/// minimum: -1.0,
41/// value: -2.0,
42/// }.calculate().unwrap(),
43/// -1.0);
44/// assert_eq!(
45/// Clip {
46/// maximum: 1.0,
47/// minimum: -1.0,
48/// value: 2.0,
49/// }.calculate().unwrap(),
50/// 1.0);
51/// assert_eq!(
52/// Clip {
53/// maximum: f64::NAN,
54/// minimum: -1.0,
55/// value: 0.0,
56/// }.calculate().unwrap_err(),
57/// ClipError::MaximumIsNotANumber);
58/// assert_eq!(
59/// Clip {
60/// maximum: f64::INFINITY,
61/// minimum: -1.0,
62/// value: 0.0,
63/// }.calculate().unwrap_err(),
64/// ClipError::MaximumIsInfinite(f64::INFINITY));
65/// assert_eq!(
66/// Clip {
67/// maximum: f64::NEG_INFINITY,
68/// minimum: -1.0,
69/// value: 0.0,
70/// }.calculate().unwrap_err(),
71/// ClipError::MaximumIsInfinite(f64::NEG_INFINITY));
72/// assert_eq!(
73/// Clip {
74/// maximum: -1.0,
75/// minimum: 1.0,
76/// value: 0.0,
77/// }.calculate().unwrap_err(),
78/// ClipError::MinimumIsGreaterThanMaximum);
79/// assert_eq!(
80/// Clip {
81/// maximum: 1.0,
82/// minimum: f64::INFINITY,
83/// value: 0.0,
84/// }.calculate().unwrap_err(),
85/// ClipError::MinimumIsInfinite(f64::INFINITY));
86/// assert_eq!(
87/// Clip {
88/// maximum: 1.0,
89/// minimum: f64::NEG_INFINITY,
90/// value: 0.0,
91/// }.calculate().unwrap_err(),
92/// ClipError::MinimumIsInfinite(f64::NEG_INFINITY));
93/// assert_eq!(
94/// Clip {
95/// maximum: 1.0,
96/// minimum: f64::NAN,
97/// value: 0.0,
98/// }.calculate().unwrap_err(),
99/// ClipError::MinimumIsNotANumber);
100/// assert_eq!(
101/// Clip {
102/// maximum: 1.0,
103/// minimum: -1.0,
104/// value: f64::INFINITY,
105/// }.calculate().unwrap_err(),
106/// ClipError::ValueIsInfinite(f64::INFINITY));
107/// assert_eq!(
108/// Clip {
109/// maximum: 1.0,
110/// minimum: -1.0,
111/// value: f64::NEG_INFINITY,
112/// }.calculate().unwrap_err(),
113/// ClipError::ValueIsInfinite(f64::NEG_INFINITY));
114/// assert_eq!(
115/// Clip {
116/// maximum: 1.0,
117/// minimum: -1.0,
118/// value: f64::NAN,
119/// }.calculate().unwrap_err(),
120/// ClipError::ValueIsNotANumber);
121/// ```
122// -----------------------------------------------------------------------------
123#[derive(Clone, Debug, PartialEq)]
124pub struct Clip {
125 pub maximum: f64,
126 pub minimum: f64,
127 pub value: f64,
128}
129
130#[derive(Debug, PartialEq)]
131pub enum ClipError {
132 MaximumIsInfinite(f64),
133 MaximumIsNotANumber,
134 MinimumIsGreaterThanMaximum,
135 MinimumIsInfinite(f64),
136 MinimumIsNotANumber,
137 ValueIsInfinite(f64),
138 ValueIsNotANumber,
139}
140
141impl Clip {
142 pub fn calculate(&self) -> Result<f64, ClipError> {
143 let max: f64 = self.maximum;
144 let min: f64 = self.minimum;
145 let val: f64 = self.value;
146 if max.is_infinite() {
147 return Err(ClipError::MaximumIsInfinite(max));
148 }
149 if min.is_infinite() {
150 return Err(ClipError::MinimumIsInfinite(min));
151 }
152 if val.is_infinite() {
153 return Err(ClipError::ValueIsInfinite(val));
154 }
155 if max.is_nan() {
156 return Err(ClipError::MaximumIsNotANumber);
157 }
158 if min.is_nan() {
159 return Err(ClipError::MinimumIsNotANumber);
160 }
161 if val.is_nan() {
162 return Err(ClipError::ValueIsNotANumber);
163 }
164 if min > max {
165 return Err(ClipError::MinimumIsGreaterThanMaximum);
166 }
167 Ok(val.clamp(min, max))
168 }
169}
170
171// -----------------------------------------------------------------------------
172/// Cumulative Distribution Function (CDF)
173///
174/// # Links
175/// <https://en.wikipedia.org/wiki/Cumulative_distribution_function>
176// -----------------------------------------------------------------------------
177#[derive(Clone, Debug, PartialEq)]
178pub struct CumulativeDistributionFunction {
179 pub x: f64,
180 pub lambda: f64,
181}
182
183impl CumulativeDistributionFunction {
184 pub fn calculate(&self) -> f64 {
185 if self.x <= 0.0 {
186 return 0.0;
187 }
188 1.0 - (-self.lambda * self.x).exp()
189 }
190}
191
192// -----------------------------------------------------------------------------
193/// Coordinates specified as angle and radius from the origin
194// -----------------------------------------------------------------------------
195#[derive(Clone, Debug, PartialEq)]
196pub struct PolarCoordinates {
197 pub angle: f64,
198 pub radius: f64,
199}
200
201impl PolarCoordinates {
202 // ---------------------------------------------------------------------------
203 /// Converts from polar to rectangular coordinates
204 ///
205 /// # Examples
206 /// ```
207 /// use com_croftsoft_core::math::math_lib::*;
208 /// assert_eq!(
209 /// PolarCoordinates {
210 /// angle: 0.0,
211 /// radius: 1.0,
212 /// }.to_rectangular_coordinates(),
213 /// RectangularCoordinates {
214 /// x: 1.0,
215 /// y: 0.0,
216 /// });
217 /// assert_eq!(
218 /// PolarCoordinates {
219 /// angle: core::f64::consts::FRAC_PI_2,
220 /// radius: 1.0,
221 /// }.to_rectangular_coordinates(),
222 /// RectangularCoordinates {
223 /// x: 6.123233995736766e-17,
224 /// y: 1.0,
225 /// });
226 /// assert_eq!(
227 /// PolarCoordinates {
228 /// angle: core::f64::consts::PI,
229 /// radius: 1.0,
230 /// }.to_rectangular_coordinates(),
231 /// RectangularCoordinates {
232 /// x: -1.0,
233 /// y: 1.2246467991473532e-16,
234 /// });
235 /// assert_eq!(
236 /// PolarCoordinates {
237 /// angle: 3.0 * core::f64::consts::FRAC_PI_2,
238 /// radius: 2.0,
239 /// }.to_rectangular_coordinates(),
240 /// RectangularCoordinates {
241 /// x: -3.6739403974420594e-16,
242 /// y: -2.0,
243 /// });
244 /// ```
245 // ---------------------------------------------------------------------------
246 pub fn to_rectangular_coordinates(&self) -> RectangularCoordinates {
247 let angle = self.angle;
248 let radius = self.radius;
249 RectangularCoordinates {
250 x: radius * angle.cos(),
251 y: radius * angle.sin(),
252 }
253 }
254}
255
256// -----------------------------------------------------------------------------
257/// Cartesian (x, y) coordinates
258// -----------------------------------------------------------------------------
259#[derive(Clone, Debug, PartialEq)]
260pub struct RectangularCoordinates {
261 pub x: f64,
262 pub y: f64,
263}
264
265// -----------------------------------------------------------------------------
266/// Wraps the value to [minimum, minimum + range)
267///
268/// # Examples
269/// ```
270/// use com_croftsoft_core::math::math_lib::*;
271/// assert_eq!(
272/// Wrap {
273/// minimum: -180.0,
274/// range: 360.0,
275/// value: -190.0,
276/// }.calculate().unwrap(),
277/// 170.0);
278/// assert_eq!(
279/// Wrap {
280/// minimum: -180.0,
281/// range: 360.0,
282/// value: 190.0,
283/// }.calculate().unwrap(),
284/// -170.0);
285/// assert_eq!(
286/// Wrap {
287/// minimum: -180.0,
288/// range: 360.0,
289/// value: 180.0,
290/// }.calculate().unwrap(),
291/// -180.0);
292/// assert_eq!(
293/// Wrap {
294/// minimum: f64::MAX,
295/// range: 360.0,
296/// value: 190.0,
297/// }.calculate().unwrap_err(),
298/// WrapError::FloatResolution(
299/// WrapErrorFloatResolution::DeltaIsNegativeMinimum));
300/// assert_eq!(
301/// Wrap {
302/// minimum: -180.0,
303/// range: 360.0,
304/// value: f64::MAX,
305/// }.calculate().unwrap_err(),
306/// WrapError::FloatResolution(
307/// WrapErrorFloatResolution::DeltaIsValue));
308/// assert_eq!(
309/// Wrap {
310/// minimum: f64::INFINITY,
311/// range: 360.0,
312/// value: 190.0,
313/// }.calculate().unwrap_err(),
314/// WrapError::InvalidArgument(
315/// WrapErrorInvalidArgument::MinimumIsInfinite(f64::INFINITY)));
316/// assert_eq!(
317/// Wrap {
318/// minimum: f64::NEG_INFINITY,
319/// range: 360.0,
320/// value: 190.0,
321/// }.calculate().unwrap_err(),
322/// WrapError::InvalidArgument(
323/// WrapErrorInvalidArgument::MinimumIsInfinite(f64::NEG_INFINITY)));
324/// assert_eq!(
325/// Wrap {
326/// minimum: f64::NAN,
327/// range: 360.0,
328/// value: 190.0,
329/// }.calculate().unwrap_err(),
330/// WrapError::InvalidArgument(
331/// WrapErrorInvalidArgument::MinimumIsNotANumber));
332/// assert_eq!(
333/// Wrap {
334/// minimum: -180.0,
335/// range: f64::INFINITY,
336/// value: 190.0,
337/// }.calculate().unwrap_err(),
338/// WrapError::InvalidArgument(
339/// WrapErrorInvalidArgument::RangeIsInfinite(f64::INFINITY)));
340/// assert_eq!(
341/// Wrap {
342/// minimum: -180.0,
343/// range: f64::NEG_INFINITY,
344/// value: 190.0,
345/// }.calculate().unwrap_err(),
346/// WrapError::InvalidArgument(
347/// WrapErrorInvalidArgument::RangeIsInfinite(f64::NEG_INFINITY)));
348/// assert_eq!(
349/// Wrap {
350/// minimum: -180.0,
351/// range: -360.0,
352/// value: 180.0,
353/// }.calculate().unwrap_err(),
354/// WrapError::InvalidArgument(
355/// WrapErrorInvalidArgument::RangeIsNonPositive(-360.0)));
356/// assert_eq!(
357/// Wrap {
358/// minimum: -180.0,
359/// range: f64::NAN,
360/// value: 190.0,
361/// }.calculate().unwrap_err(),
362/// WrapError::InvalidArgument(
363/// WrapErrorInvalidArgument::RangeIsNotANumber));
364/// assert_eq!(
365/// Wrap {
366/// minimum: -180.0,
367/// range: 360.0,
368/// value: f64::INFINITY,
369/// }.calculate().unwrap_err(),
370/// WrapError::InvalidArgument(
371/// WrapErrorInvalidArgument::ValueIsInfinite(f64::INFINITY)));
372/// assert_eq!(
373/// Wrap {
374/// minimum: -180.0,
375/// range: 360.0,
376/// value: f64::NEG_INFINITY,
377/// }.calculate().unwrap_err(),
378/// WrapError::InvalidArgument(
379/// WrapErrorInvalidArgument::ValueIsInfinite(f64::NEG_INFINITY)));
380/// assert_eq!(
381/// Wrap {
382/// minimum: -180.0,
383/// range: 360.0,
384/// value: f64::NAN,
385/// }.calculate().unwrap_err(),
386/// WrapError::InvalidArgument(
387/// WrapErrorInvalidArgument::ValueIsNotANumber));
388/// ```
389// -----------------------------------------------------------------------------
390#[derive(Clone, Debug, PartialEq)]
391pub struct Wrap {
392 pub minimum: f64,
393 pub range: f64,
394 pub value: f64,
395}
396
397#[derive(Debug, PartialEq)]
398pub enum WrapError {
399 FloatResolution(WrapErrorFloatResolution),
400 InvalidArgument(WrapErrorInvalidArgument),
401}
402
403#[derive(Debug, PartialEq)]
404pub enum WrapErrorFloatResolution {
405 CalculatedIsLessThanMinimum(f64),
406 CalculatedIsNotLessThanMinimumPlusRange(f64),
407 CyclesIsZero,
408 DeltaIsNegativeMinimum,
409 DeltaIsValue,
410 OffsetIsZero,
411}
412
413#[derive(Debug, PartialEq)]
414pub enum WrapErrorInvalidArgument {
415 MinimumIsInfinite(f64),
416 MinimumIsNotANumber,
417 RangeIsInfinite(f64),
418 RangeIsNonPositive(f64),
419 RangeIsNotANumber,
420 ValueIsInfinite(f64),
421 ValueIsNotANumber,
422}
423
424impl Wrap {
425 pub fn calculate(&self) -> Result<f64, WrapError> {
426 let min = self.minimum;
427 let rng = self.range;
428 let val = self.value;
429 if min.is_infinite() {
430 return Err(WrapError::InvalidArgument(
431 WrapErrorInvalidArgument::MinimumIsInfinite(min),
432 ));
433 }
434 if rng.is_infinite() {
435 return Err(WrapError::InvalidArgument(
436 WrapErrorInvalidArgument::RangeIsInfinite(rng),
437 ));
438 }
439 if val.is_infinite() {
440 return Err(WrapError::InvalidArgument(
441 WrapErrorInvalidArgument::ValueIsInfinite(val),
442 ));
443 }
444 if min.is_nan() {
445 return Err(WrapError::InvalidArgument(
446 WrapErrorInvalidArgument::MinimumIsNotANumber,
447 ));
448 }
449 if rng.is_nan() {
450 return Err(WrapError::InvalidArgument(
451 WrapErrorInvalidArgument::RangeIsNotANumber,
452 ));
453 }
454 if val.is_nan() {
455 return Err(WrapError::InvalidArgument(
456 WrapErrorInvalidArgument::ValueIsNotANumber,
457 ));
458 }
459 if rng <= 0.0 {
460 return Err(WrapError::InvalidArgument(
461 WrapErrorInvalidArgument::RangeIsNonPositive(rng),
462 ));
463 }
464 let max = min + rng;
465 if min <= val && val < max {
466 return Ok(val);
467 }
468 let delta = val - min;
469 if delta == -min {
470 return Err(WrapError::FloatResolution(
471 WrapErrorFloatResolution::DeltaIsNegativeMinimum,
472 ));
473 }
474 if delta == val {
475 return Err(WrapError::FloatResolution(
476 WrapErrorFloatResolution::DeltaIsValue,
477 ));
478 }
479 let cycles = (delta / rng).floor();
480 if cycles == 0.0 {
481 return Err(WrapError::FloatResolution(
482 WrapErrorFloatResolution::CyclesIsZero,
483 ));
484 }
485 let offset = cycles * rng;
486 if offset == 0.0 {
487 return Err(WrapError::FloatResolution(
488 WrapErrorFloatResolution::OffsetIsZero,
489 ));
490 }
491 let calculated = val - offset;
492 if calculated < min {
493 return Err(WrapError::FloatResolution(
494 WrapErrorFloatResolution::CalculatedIsLessThanMinimum(calculated),
495 ));
496 }
497 if calculated >= max {
498 return Err(WrapError::FloatResolution(
499 WrapErrorFloatResolution::CalculatedIsNotLessThanMinimumPlusRange(
500 calculated,
501 ),
502 ));
503 }
504 Ok(calculated)
505 }
506}
507
508#[derive(Debug, Eq, PartialEq)]
509pub enum FactorError {
510 ArgumentIsZeroOrOne(u64),
511}
512
513// -----------------------------------------------------------------------------
514/// Factors a number into its primes
515/// ```
516/// use com_croftsoft_core::math::math_lib::*;
517/// assert_eq!(factor(0).unwrap_err(), FactorError::ArgumentIsZeroOrOne(0));
518/// assert_eq!(factor(1).unwrap_err(), FactorError::ArgumentIsZeroOrOne(1));
519/// assert_eq!(factor(2).unwrap(), vec!(2));
520/// assert_eq!(factor(3).unwrap(), vec!(3));
521/// assert_eq!(factor(4).unwrap(), vec!(2, 2));
522/// assert_eq!(factor(5).unwrap(), vec!(5));
523/// assert_eq!(factor(6).unwrap(), vec!(2, 3));
524/// assert_eq!(factor(7).unwrap(), vec!(7));
525/// assert_eq!(factor(8).unwrap(), vec!(2, 2, 2));
526/// assert_eq!(factor(9).unwrap(), vec!(3, 3));
527/// ```
528// -----------------------------------------------------------------------------
529pub fn factor(n: u64) -> Result<Vec<u64>, FactorError> {
530 if n < 2 {
531 return Err(FactorError::ArgumentIsZeroOrOne(n));
532 }
533 let mut prime_vec = Vec::new();
534 let mut dividend = n;
535 let mut divisor = 2;
536 loop {
537 if dividend % divisor == 0 {
538 prime_vec.push(divisor);
539 dividend /= divisor;
540 if dividend == 1 {
541 break;
542 }
543 } else {
544 divisor += 1;
545 }
546 }
547 Ok(prime_vec)
548}
549
550#[derive(Debug, Eq, PartialEq)]
551pub enum GreatestCommonFactorError {
552 ArgumentIsZero,
553}
554
555// -----------------------------------------------------------------------------
556/// Computes the greatest common factor for two positive integers
557///
558/// ```
559/// use com_croftsoft_core::math::math_lib::*;
560/// assert_eq!(
561/// greatest_common_factor(0, 1).unwrap_err(),
562/// GreatestCommonFactorError::ArgumentIsZero);
563/// assert_eq!(greatest_common_factor(1, 2).unwrap(), 1);
564/// assert_eq!(greatest_common_factor(2, 3).unwrap(), 1);
565/// assert_eq!(greatest_common_factor(2, 4).unwrap(), 2);
566/// assert_eq!(greatest_common_factor(3, 6).unwrap(), 3);
567/// ```
568// -----------------------------------------------------------------------------
569pub fn greatest_common_factor(
570 n0: u64,
571 n1: u64,
572) -> Result<u64, GreatestCommonFactorError> {
573 if n0 == 0 || n1 == 0 {
574 return Err(GreatestCommonFactorError::ArgumentIsZero);
575 }
576 if n0 == 1 || n1 == 1 {
577 return Ok(1);
578 }
579 let mut gcf: u64 = 1;
580 let factor_vec_0 = factor(n0).unwrap();
581 let mut factor_vec_1 = factor(n1).unwrap();
582 for (index, factor_0) in factor_vec_0.iter().enumerate() {
583 if factor_vec_1.contains(factor_0) {
584 gcf *= factor_0;
585 factor_vec_1.remove(index);
586 }
587 }
588 Ok(gcf)
589}
590
591// -----------------------------------------------------------------------------
592/// The sigmoid or logistic function
593///
594/// # Examples
595/// ```
596/// use com_croftsoft_core::math::math_lib::sigmoid;
597/// assert_eq!(
598/// sigmoid(f64::NEG_INFINITY),
599/// 0.0);
600/// assert_eq!(
601/// sigmoid(-1.0),
602/// 0.2689414213699951);
603/// assert_eq!(
604/// sigmoid(0.0),
605/// 0.5);
606/// assert_eq!(
607/// sigmoid(1.0),
608/// 0.7310585786300049);
609/// assert_eq!(
610/// sigmoid(f64::INFINITY),
611/// 1.0);
612/// ```
613// -----------------------------------------------------------------------------
614pub fn sigmoid(a: f64) -> f64 {
615 1.0 / (1.0 + (-a).exp())
616}
617
618// -----------------------------------------------------------------------------
619/// The derivative of the sigmoid function with respect to the argument
620///
621/// # Examples
622/// ```
623/// use com_croftsoft_core::math::math_lib::sigmoid_derivative;
624/// assert_eq!(
625/// sigmoid_derivative(f64::NEG_INFINITY),
626/// 0.0);
627/// assert_eq!(
628/// sigmoid_derivative(-1.0),
629/// 0.19661193324148185);
630/// assert_eq!(
631/// sigmoid_derivative(0.0),
632/// 0.25);
633/// assert_eq!(
634/// sigmoid_derivative(1.0),
635/// 0.19661193324148185);
636/// assert_eq!(
637/// sigmoid_derivative(f64::INFINITY),
638/// 0.0);
639/// ```
640// -----------------------------------------------------------------------------
641pub fn sigmoid_derivative(a: f64) -> f64 {
642 let s = sigmoid(a);
643 s * (1.0 - s)
644}