1#![allow(clippy::needless_lifetimes, clippy::suspicious_arithmetic_impl)]
46
47use crate::passive::Passive;
48use std::fmt;
49use std::ops::{Add, AddAssign, Div, DivAssign, Mul, MulAssign, Neg, Sub, SubAssign};
50
51#[derive(Clone, Debug)]
66pub struct Jet1Vec<T: Passive = f64> {
67 pub real: T,
68 pub dual: Vec<T>,
69}
70
71#[inline]
72fn len_check<T: Passive>(a: &Jet1Vec<T>, b: &Jet1Vec<T>) {
73 debug_assert_eq!(a.dual.len(), b.dual.len(), "Jet1Vec tangent length mismatch");
74}
75
76impl<T: Passive> Jet1Vec<T> {
77 #[inline]
79 pub fn new(real: T, dual: Vec<T>) -> Self {
80 Jet1Vec { real, dual }
81 }
82
83 #[inline]
85 pub fn constant(real: T, n: usize) -> Self {
86 Jet1Vec {
87 real,
88 dual: vec![T::zero(); n],
89 }
90 }
91
92 #[inline]
98 pub fn variable(value: T, i: usize, n: usize) -> Self {
99 let mut dual = vec![T::zero(); n];
100 dual[i] = T::one();
101 Jet1Vec { real: value, dual }
102 }
103
104 #[inline]
106 pub fn real(&self) -> T {
107 self.real
108 }
109
110 #[inline]
112 pub fn dual(&self) -> &[T] {
113 &self.dual
114 }
115
116 #[inline]
118 pub fn partial(&self, i: usize) -> T {
119 self.dual[i]
120 }
121
122 #[inline]
124 pub fn num_vars(&self) -> usize {
125 self.dual.len()
126 }
127
128 #[inline]
135 fn scaled(&self, k: T) -> Vec<T> {
136 self.dual.iter().map(|&d| d * k).collect()
137 }
138
139 #[inline]
144 fn scaled_into(mut self, k: T) -> Vec<T> {
145 for d in self.dual.iter_mut() {
146 *d *= k;
147 }
148 self.dual
149 }
150
151 #[inline]
153 pub fn powf(&self, n: T) -> Jet1Vec<T> {
154 let vn = self.real.powf(n);
155 let gp = n * self.real.powf(n - T::one());
156 Jet1Vec {
157 real: vn,
158 dual: self.scaled(gp),
159 }
160 }
161
162 #[inline]
164 pub fn powi(&self, n: i32) -> Jet1Vec<T> {
165 let vn = self.real.powi(n);
166 let gp = T::from(n).unwrap() * self.real.powi(n - 1);
167 Jet1Vec {
168 real: vn,
169 dual: self.scaled(gp),
170 }
171 }
172
173 #[inline]
175 pub fn exp(&self) -> Jet1Vec<T> {
176 let e = self.real.exp();
177 Jet1Vec {
178 real: e,
179 dual: self.scaled(e),
180 }
181 }
182
183 #[inline]
185 pub fn ln(&self) -> Jet1Vec<T> {
186 let inv = T::one() / self.real;
187 Jet1Vec {
188 real: self.real.ln(),
189 dual: self.scaled(inv),
190 }
191 }
192
193 #[inline]
195 pub fn sqrt(&self) -> Jet1Vec<T> {
196 let s = self.real.sqrt();
197 Jet1Vec {
198 real: s,
199 dual: self.scaled(T::from(0.5).unwrap() / s),
200 }
201 }
202
203 #[inline]
205 pub fn sin(&self) -> Jet1Vec<T> {
206 Jet1Vec {
207 real: self.real.sin(),
208 dual: self.scaled(self.real.cos()),
209 }
210 }
211
212 #[inline]
214 pub fn cos(&self) -> Jet1Vec<T> {
215 Jet1Vec {
216 real: self.real.cos(),
217 dual: self.scaled(-self.real.sin()),
218 }
219 }
220
221 #[inline]
223 pub fn tan(&self) -> Jet1Vec<T> {
224 let t = self.real.tan();
225 let sec2 = T::one() + t * t;
226 Jet1Vec {
227 real: t,
228 dual: self.scaled(sec2),
229 }
230 }
231
232 #[inline]
234 pub fn tanh(&self) -> Jet1Vec<T> {
235 let t = self.real.tanh();
236 Jet1Vec {
237 real: t,
238 dual: self.scaled(T::one() - t * t),
239 }
240 }
241
242 #[inline]
244 pub fn abs(&self) -> Jet1Vec<T> {
245 let sign = if self.real >= T::zero() {
246 T::one()
247 } else {
248 -T::one()
249 };
250 Jet1Vec {
251 real: self.real.abs(),
252 dual: self.scaled(sign),
253 }
254 }
255
256 #[inline]
260 pub fn erf(&self) -> Jet1Vec<T> {
261 let v = self.real;
262 let r = crate::math::erf::<T>(v);
263 let g = T::from(std::f64::consts::FRAC_2_SQRT_PI).unwrap() * (-v * v).exp();
267 Jet1Vec {
268 real: r,
269 dual: self.scaled(g),
270 }
271 }
272
273 #[inline]
277 pub fn norm_cdf(&self) -> Jet1Vec<T> {
278 let v = self.real;
279 let r = crate::math::norm_cdf::<T>(v);
280 let g = crate::math::norm_pdf::<T>(v);
281 Jet1Vec {
282 real: r,
283 dual: self.scaled(g),
284 }
285 }
286
287 #[inline]
291 pub fn inv_norm_cdf(&self) -> Jet1Vec<T> {
292 let v = self.real;
293 let r = crate::math::inv_norm_cdf::<T>(v);
294 let g = T::one() / crate::math::norm_pdf::<T>(r);
295 Jet1Vec {
296 real: r,
297 dual: self.scaled(g),
298 }
299 }
300
301 #[inline]
314 pub fn into_powf(self, n: T) -> Jet1Vec<T> {
315 let real = self.real;
316 let vn = real.powf(n);
317 let gp = n * real.powf(n - T::one());
318 Jet1Vec {
319 real: vn,
320 dual: self.scaled_into(gp),
321 }
322 }
323
324 #[inline]
326 pub fn into_powi(self, n: i32) -> Jet1Vec<T> {
327 let real = self.real;
328 let vn = real.powi(n);
329 let gp = T::from(n).unwrap() * real.powi(n - 1);
330 Jet1Vec {
331 real: vn,
332 dual: self.scaled_into(gp),
333 }
334 }
335
336 #[inline]
338 pub fn into_exp(self) -> Jet1Vec<T> {
339 let e = self.real.exp();
340 Jet1Vec {
341 real: e,
342 dual: self.scaled_into(e),
343 }
344 }
345
346 #[inline]
348 pub fn into_ln(self) -> Jet1Vec<T> {
349 let real = self.real;
350 Jet1Vec {
351 real: real.ln(),
352 dual: self.scaled_into(T::one() / real),
353 }
354 }
355
356 #[inline]
358 pub fn into_sqrt(self) -> Jet1Vec<T> {
359 let s = self.real.sqrt();
360 Jet1Vec {
361 real: s,
362 dual: self.scaled_into(T::from(0.5).unwrap() / s),
363 }
364 }
365
366 #[inline]
368 pub fn into_sin(self) -> Jet1Vec<T> {
369 let real = self.real;
370 Jet1Vec {
371 real: real.sin(),
372 dual: self.scaled_into(real.cos()),
373 }
374 }
375
376 #[inline]
378 pub fn into_cos(self) -> Jet1Vec<T> {
379 let real = self.real;
380 Jet1Vec {
381 real: real.cos(),
382 dual: self.scaled_into(-real.sin()),
383 }
384 }
385
386 #[inline]
388 pub fn into_tan(self) -> Jet1Vec<T> {
389 let t = self.real.tan();
390 let sec2 = T::one() + t * t;
391 Jet1Vec {
392 real: t,
393 dual: self.scaled_into(sec2),
394 }
395 }
396
397 #[inline]
399 pub fn into_tanh(self) -> Jet1Vec<T> {
400 let t = self.real.tanh();
401 Jet1Vec {
402 real: t,
403 dual: self.scaled_into(T::one() - t * t),
404 }
405 }
406
407 #[inline]
409 pub fn into_abs(self) -> Jet1Vec<T> {
410 let sign = if self.real >= T::zero() {
411 T::one()
412 } else {
413 -T::one()
414 };
415 Jet1Vec {
416 real: self.real.abs(),
417 dual: self.scaled_into(sign),
418 }
419 }
420
421 #[inline]
423 pub fn into_erf(self) -> Jet1Vec<T> {
424 let v = self.real;
425 let r = crate::math::erf::<T>(v);
426 let g = T::from(std::f64::consts::FRAC_2_SQRT_PI).unwrap() * (-v * v).exp();
427 Jet1Vec {
428 real: r,
429 dual: self.scaled_into(g),
430 }
431 }
432
433 #[inline]
435 pub fn into_norm_cdf(self) -> Jet1Vec<T> {
436 let v = self.real;
437 let r = crate::math::norm_cdf::<T>(v);
438 let g = crate::math::norm_pdf::<T>(v);
439 Jet1Vec {
440 real: r,
441 dual: self.scaled_into(g),
442 }
443 }
444
445 #[inline]
447 pub fn into_inv_norm_cdf(self) -> Jet1Vec<T> {
448 let v = self.real;
449 let r = crate::math::inv_norm_cdf::<T>(v);
450 let g = T::one() / crate::math::norm_pdf::<T>(r);
451 Jet1Vec {
452 real: r,
453 dual: self.scaled_into(g),
454 }
455 }
456}
457
458impl<'a, 'b, T: Passive> Add<&'b Jet1Vec<T>> for &'a Jet1Vec<T> {
468 type Output = Jet1Vec<T>;
469 #[inline]
470 fn add(self, rhs: &'b Jet1Vec<T>) -> Jet1Vec<T> {
471 len_check(self, rhs);
472 Jet1Vec {
473 real: self.real + rhs.real,
474 dual: self
475 .dual
476 .iter()
477 .zip(&rhs.dual)
478 .map(|(a, b)| *a + *b)
479 .collect(),
480 }
481 }
482}
483
484impl<T: Passive> Add for Jet1Vec<T> {
485 type Output = Jet1Vec<T>;
486 #[inline]
487 fn add(mut self, rhs: Jet1Vec<T>) -> Jet1Vec<T> {
488 self += &rhs;
489 self
490 }
491}
492
493impl<T: Passive> Add<&Jet1Vec<T>> for Jet1Vec<T> {
494 type Output = Jet1Vec<T>;
495 #[inline]
496 fn add(mut self, rhs: &Jet1Vec<T>) -> Jet1Vec<T> {
497 self += rhs;
498 self
499 }
500}
501
502impl<T: Passive> Add<Jet1Vec<T>> for &Jet1Vec<T> {
503 type Output = Jet1Vec<T>;
504 #[inline]
505 fn add(self, mut rhs: Jet1Vec<T>) -> Jet1Vec<T> {
506 rhs += self;
507 rhs
508 }
509}
510
511impl<T: Passive> Add<T> for &Jet1Vec<T> {
512 type Output = Jet1Vec<T>;
513 #[inline]
514 fn add(self, rhs: T) -> Jet1Vec<T> {
515 Jet1Vec {
516 real: self.real + rhs,
517 dual: self.dual.clone(),
518 }
519 }
520}
521
522impl<T: Passive> Add<T> for Jet1Vec<T> {
523 type Output = Jet1Vec<T>;
524 #[inline]
525 fn add(mut self, rhs: T) -> Jet1Vec<T> {
526 self.real += rhs;
527 self
528 }
529}
530
531impl<'a, 'b, T: Passive> Sub<&'b Jet1Vec<T>> for &'a Jet1Vec<T> {
533 type Output = Jet1Vec<T>;
534 #[inline]
535 fn sub(self, rhs: &'b Jet1Vec<T>) -> Jet1Vec<T> {
536 len_check(self, rhs);
537 Jet1Vec {
538 real: self.real - rhs.real,
539 dual: self
540 .dual
541 .iter()
542 .zip(&rhs.dual)
543 .map(|(a, b)| *a - *b)
544 .collect(),
545 }
546 }
547}
548
549impl<T: Passive> Sub for Jet1Vec<T> {
550 type Output = Jet1Vec<T>;
551 #[inline]
552 fn sub(mut self, rhs: Jet1Vec<T>) -> Jet1Vec<T> {
553 self -= &rhs;
554 self
555 }
556}
557
558impl<T: Passive> Sub<&Jet1Vec<T>> for Jet1Vec<T> {
559 type Output = Jet1Vec<T>;
560 #[inline]
561 fn sub(mut self, rhs: &Jet1Vec<T>) -> Jet1Vec<T> {
562 self -= rhs;
563 self
564 }
565}
566
567impl<T: Passive> Sub<Jet1Vec<T>> for &Jet1Vec<T> {
568 type Output = Jet1Vec<T>;
569 #[inline]
570 fn sub(self, mut rhs: Jet1Vec<T>) -> Jet1Vec<T> {
571 len_check(self, &rhs);
573 rhs.real = self.real - rhs.real;
574 for (b, &a) in rhs.dual.iter_mut().zip(&self.dual) {
575 *b = a - *b;
576 }
577 rhs
578 }
579}
580
581impl<T: Passive> Sub<T> for &Jet1Vec<T> {
582 type Output = Jet1Vec<T>;
583 #[inline]
584 fn sub(self, rhs: T) -> Jet1Vec<T> {
585 Jet1Vec {
586 real: self.real - rhs,
587 dual: self.dual.clone(),
588 }
589 }
590}
591
592impl<T: Passive> Sub<T> for Jet1Vec<T> {
593 type Output = Jet1Vec<T>;
594 #[inline]
595 fn sub(mut self, rhs: T) -> Jet1Vec<T> {
596 self.real -= rhs;
597 self
598 }
599}
600
601impl<'a, 'b, T: Passive> Mul<&'b Jet1Vec<T>> for &'a Jet1Vec<T> {
604 type Output = Jet1Vec<T>;
605 #[inline]
606 fn mul(self, rhs: &'b Jet1Vec<T>) -> Jet1Vec<T> {
607 len_check(self, rhs);
608 let a = self.real;
609 let b = rhs.real;
610 Jet1Vec {
611 real: a * b,
612 dual: self
613 .dual
614 .iter()
615 .zip(&rhs.dual)
616 .map(|(&da, &db)| da * b + db * a)
617 .collect(),
618 }
619 }
620}
621
622impl<T: Passive> Mul for Jet1Vec<T> {
623 type Output = Jet1Vec<T>;
624 #[inline]
625 fn mul(mut self, rhs: Jet1Vec<T>) -> Jet1Vec<T> {
626 len_check(&self, &rhs);
627 let a = self.real;
628 let b = rhs.real;
629 self.real = a * b;
630 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
631 *da = *da * b + db * a;
632 }
633 self
634 }
635}
636
637impl<T: Passive> Mul<&Jet1Vec<T>> for Jet1Vec<T> {
638 type Output = Jet1Vec<T>;
639 #[inline]
640 fn mul(mut self, rhs: &Jet1Vec<T>) -> Jet1Vec<T> {
641 len_check(&self, rhs);
642 let a = self.real;
643 let b = rhs.real;
644 self.real = a * b;
645 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
646 *da = *da * b + db * a;
647 }
648 self
649 }
650}
651
652impl<T: Passive> Mul<Jet1Vec<T>> for &Jet1Vec<T> {
653 type Output = Jet1Vec<T>;
654 #[inline]
655 fn mul(self, mut rhs: Jet1Vec<T>) -> Jet1Vec<T> {
656 len_check(self, &rhs);
657 let a = self.real;
658 let b = rhs.real;
659 rhs.real = a * b;
660 for (db, &da) in rhs.dual.iter_mut().zip(&self.dual) {
661 *db = da * b + *db * a;
662 }
663 rhs
664 }
665}
666
667impl<T: Passive> Mul<T> for &Jet1Vec<T> {
668 type Output = Jet1Vec<T>;
669 #[inline]
670 fn mul(self, rhs: T) -> Jet1Vec<T> {
671 Jet1Vec {
672 real: self.real * rhs,
673 dual: self.dual.iter().map(|&d| d * rhs).collect(),
674 }
675 }
676}
677
678impl<T: Passive> Mul<T> for Jet1Vec<T> {
679 type Output = Jet1Vec<T>;
680 #[inline]
681 fn mul(mut self, rhs: T) -> Jet1Vec<T> {
682 self.real *= rhs;
683 for d in self.dual.iter_mut() {
684 *d *= rhs;
685 }
686 self
687 }
688}
689
690impl<'a, 'b, T: Passive> Div<&'b Jet1Vec<T>> for &'a Jet1Vec<T> {
693 type Output = Jet1Vec<T>;
694 #[inline]
695 fn div(self, rhs: &'b Jet1Vec<T>) -> Jet1Vec<T> {
696 len_check(self, rhs);
697 let inv = T::one() / rhs.real;
698 let inv2 = inv * inv;
699 let a_inv2 = self.real * inv2;
700 Jet1Vec {
701 real: self.real * inv,
702 dual: self
703 .dual
704 .iter()
705 .zip(&rhs.dual)
706 .map(|(&da, &db)| da * inv - db * a_inv2)
707 .collect(),
708 }
709 }
710}
711
712impl<T: Passive> Div for Jet1Vec<T> {
713 type Output = Jet1Vec<T>;
714 #[inline]
715 fn div(mut self, rhs: Jet1Vec<T>) -> Jet1Vec<T> {
716 len_check(&self, &rhs);
717 let inv = T::one() / rhs.real;
718 let inv2 = inv * inv;
719 let a_inv2 = self.real * inv2;
720 self.real *= inv;
721 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
722 *da = *da * inv - db * a_inv2;
723 }
724 self
725 }
726}
727
728impl<T: Passive> Div<&Jet1Vec<T>> for Jet1Vec<T> {
729 type Output = Jet1Vec<T>;
730 #[inline]
731 fn div(mut self, rhs: &Jet1Vec<T>) -> Jet1Vec<T> {
732 len_check(&self, rhs);
733 let inv = T::one() / rhs.real;
734 let inv2 = inv * inv;
735 let a_inv2 = self.real * inv2;
736 self.real *= inv;
737 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
738 *da = *da * inv - db * a_inv2;
739 }
740 self
741 }
742}
743
744impl<T: Passive> Div<Jet1Vec<T>> for &Jet1Vec<T> {
745 type Output = Jet1Vec<T>;
746 #[inline]
747 fn div(self, rhs: Jet1Vec<T>) -> Jet1Vec<T> {
748 self / &rhs
752 }
753}
754
755impl<T: Passive> Div<T> for &Jet1Vec<T> {
756 type Output = Jet1Vec<T>;
757 #[inline]
758 fn div(self, rhs: T) -> Jet1Vec<T> {
759 let inv = T::one() / rhs;
760 Jet1Vec {
761 real: self.real * inv,
762 dual: self.dual.iter().map(|&d| d * inv).collect(),
763 }
764 }
765}
766
767impl<T: Passive> Div<T> for Jet1Vec<T> {
768 type Output = Jet1Vec<T>;
769 #[inline]
770 fn div(mut self, rhs: T) -> Jet1Vec<T> {
771 let inv = T::one() / rhs;
772 self.real *= inv;
773 for d in self.dual.iter_mut() {
774 *d *= inv;
775 }
776 self
777 }
778}
779
780impl<T: Passive> Neg for &Jet1Vec<T> {
782 type Output = Jet1Vec<T>;
783 #[inline]
784 fn neg(self) -> Jet1Vec<T> {
785 Jet1Vec {
786 real: -self.real,
787 dual: self.dual.iter().map(|&d| -d).collect(),
788 }
789 }
790}
791
792impl<T: Passive> Neg for Jet1Vec<T> {
793 type Output = Jet1Vec<T>;
794 #[inline]
795 fn neg(mut self) -> Jet1Vec<T> {
796 self.real = -self.real;
797 for d in self.dual.iter_mut() {
798 *d = -*d;
799 }
800 self
801 }
802}
803
804impl<T: Passive> AddAssign<&Jet1Vec<T>> for Jet1Vec<T> {
806 #[inline]
807 fn add_assign(&mut self, rhs: &Jet1Vec<T>) {
808 len_check(self, rhs);
809 self.real += rhs.real;
810 for (a, &b) in self.dual.iter_mut().zip(&rhs.dual) {
811 *a += b;
812 }
813 }
814}
815
816impl<T: Passive> AddAssign for Jet1Vec<T> {
817 #[inline]
818 fn add_assign(&mut self, rhs: Jet1Vec<T>) {
819 len_check(self, &rhs);
820 self.real += rhs.real;
821 for (a, &b) in self.dual.iter_mut().zip(&rhs.dual) {
822 *a += b;
823 }
824 }
825}
826
827impl<T: Passive> AddAssign<T> for Jet1Vec<T> {
828 #[inline]
829 fn add_assign(&mut self, rhs: T) {
830 self.real += rhs;
831 }
832}
833
834impl<T: Passive> SubAssign<&Jet1Vec<T>> for Jet1Vec<T> {
835 #[inline]
836 fn sub_assign(&mut self, rhs: &Jet1Vec<T>) {
837 len_check(self, rhs);
838 self.real -= rhs.real;
839 for (a, &b) in self.dual.iter_mut().zip(&rhs.dual) {
840 *a -= b;
841 }
842 }
843}
844
845impl<T: Passive> SubAssign for Jet1Vec<T> {
846 #[inline]
847 fn sub_assign(&mut self, rhs: Jet1Vec<T>) {
848 len_check(self, &rhs);
849 self.real -= rhs.real;
850 for (a, &b) in self.dual.iter_mut().zip(&rhs.dual) {
851 *a -= b;
852 }
853 }
854}
855
856impl<T: Passive> SubAssign<T> for Jet1Vec<T> {
857 #[inline]
858 fn sub_assign(&mut self, rhs: T) {
859 self.real -= rhs;
860 }
861}
862
863impl<T: Passive> MulAssign<&Jet1Vec<T>> for Jet1Vec<T> {
864 #[inline]
865 fn mul_assign(&mut self, rhs: &Jet1Vec<T>) {
866 len_check(self, rhs);
867 let a = self.real;
868 let b = rhs.real;
869 self.real = a * b;
870 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
871 *da = *da * b + db * a;
872 }
873 }
874}
875
876impl<T: Passive> MulAssign<T> for Jet1Vec<T> {
877 #[inline]
878 fn mul_assign(&mut self, rhs: T) {
879 self.real *= rhs;
880 for d in self.dual.iter_mut() {
881 *d *= rhs;
882 }
883 }
884}
885
886impl<T: Passive> DivAssign<&Jet1Vec<T>> for Jet1Vec<T> {
887 #[inline]
888 fn div_assign(&mut self, rhs: &Jet1Vec<T>) {
889 len_check(self, rhs);
890 let inv = T::one() / rhs.real;
891 let inv2 = inv * inv;
892 let a_inv2 = self.real * inv2;
893 self.real *= inv;
894 for (da, &db) in self.dual.iter_mut().zip(&rhs.dual) {
895 *da = *da * inv - db * a_inv2;
896 }
897 }
898}
899
900impl<T: Passive> DivAssign<T> for Jet1Vec<T> {
901 #[inline]
902 fn div_assign(&mut self, rhs: T) {
903 let inv = T::one() / rhs;
904 self.real *= inv;
905 for d in self.dual.iter_mut() {
906 *d *= inv;
907 }
908 }
909}
910
911impl<T: Passive> PartialEq for Jet1Vec<T> {
913 fn eq(&self, other: &Self) -> bool {
914 self.real == other.real
915 }
916}
917
918impl<T: Passive> PartialEq<T> for Jet1Vec<T> {
919 fn eq(&self, other: &T) -> bool {
920 self.real == *other
921 }
922}
923
924impl<T: Passive> PartialOrd for Jet1Vec<T> {
925 fn partial_cmp(&self, other: &Self) -> Option<std::cmp::Ordering> {
926 self.real.partial_cmp(&other.real)
927 }
928}
929
930impl<T: Passive> PartialOrd<T> for Jet1Vec<T> {
931 fn partial_cmp(&self, other: &T) -> Option<std::cmp::Ordering> {
932 self.real.partial_cmp(other)
933 }
934}
935
936impl<T: Passive> fmt::Display for Jet1Vec<T> {
938 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
939 write!(f, "{}", self.real)
940 }
941}
942
943macro_rules! impl_scalar_lhs_jet1vec_binop {
953 ($scalar:ty) => {
954 impl Add<Jet1Vec<$scalar>> for $scalar {
956 type Output = Jet1Vec<$scalar>;
957 #[inline]
958 fn add(self, rhs: Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
959 rhs + self
960 }
961 }
962 impl Add<&Jet1Vec<$scalar>> for $scalar {
963 type Output = Jet1Vec<$scalar>;
964 #[inline]
965 fn add(self, rhs: &Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
966 rhs + self
967 }
968 }
969 impl Sub<Jet1Vec<$scalar>> for $scalar {
971 type Output = Jet1Vec<$scalar>;
972 #[inline]
973 fn sub(self, mut rhs: Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
974 rhs.real = self - rhs.real;
975 for d in rhs.dual.iter_mut() {
976 *d = -*d;
977 }
978 rhs
979 }
980 }
981 impl Sub<&Jet1Vec<$scalar>> for $scalar {
982 type Output = Jet1Vec<$scalar>;
983 #[inline]
984 fn sub(self, rhs: &Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
985 Jet1Vec {
986 real: self - rhs.real,
987 dual: rhs.dual.iter().map(|&d| -d).collect(),
988 }
989 }
990 }
991 impl Mul<Jet1Vec<$scalar>> for $scalar {
993 type Output = Jet1Vec<$scalar>;
994 #[inline]
995 fn mul(self, rhs: Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
996 rhs * self
997 }
998 }
999 impl Mul<&Jet1Vec<$scalar>> for $scalar {
1000 type Output = Jet1Vec<$scalar>;
1001 #[inline]
1002 fn mul(self, rhs: &Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
1003 rhs * self
1004 }
1005 }
1006 impl Div<Jet1Vec<$scalar>> for $scalar {
1008 type Output = Jet1Vec<$scalar>;
1009 #[inline]
1010 fn div(self, mut rhs: Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
1011 let inv = 1.0 as $scalar / rhs.real;
1012 let neg_scale = -self * inv * inv;
1013 rhs.real = self * inv;
1014 for d in rhs.dual.iter_mut() {
1015 *d *= neg_scale;
1016 }
1017 rhs
1018 }
1019 }
1020 impl Div<&Jet1Vec<$scalar>> for $scalar {
1021 type Output = Jet1Vec<$scalar>;
1022 #[inline]
1023 fn div(self, rhs: &Jet1Vec<$scalar>) -> Jet1Vec<$scalar> {
1024 let inv = 1.0 as $scalar / rhs.real;
1025 let neg_scale = -self * inv * inv;
1026 Jet1Vec {
1027 real: self * inv,
1028 dual: rhs.dual.iter().map(|&d| d * neg_scale).collect(),
1029 }
1030 }
1031 }
1032 };
1033}
1034
1035impl_scalar_lhs_jet1vec_binop!(f64);
1036impl_scalar_lhs_jet1vec_binop!(f32);
1037
1038#[derive(Clone, Debug)]
1077pub struct NamedJet1Vec<T: Passive = f64> {
1078 pub(crate) inner: Jet1Vec<T>,
1079 #[cfg(debug_assertions)]
1083 pub(crate) gen_id: u64,
1084}
1085
1086impl<T: Passive> NamedJet1Vec<T> {
1087 #[inline]
1091 pub(crate) fn __from_inner(inner: Jet1Vec<T>) -> Self {
1092 Self {
1093 inner,
1094 #[cfg(debug_assertions)]
1095 gen_id: crate::forward_tape::current_gen(),
1096 }
1097 }
1098
1099 #[inline]
1101 pub fn real(&self) -> T {
1102 self.inner.real
1103 }
1104
1105 pub fn partial(&self, name: &str) -> T {
1113 let idx = crate::forward_tape::with_active_registry(|r| {
1114 let r =
1115 r.expect("NamedJet1Vec::partial called outside a frozen NamedForwardTape scope");
1116 r.index_of(name).unwrap_or_else(|| {
1117 panic!(
1118 "NamedJet1Vec::partial: name {:?} not present in registry",
1119 name
1120 )
1121 })
1122 });
1123 self.inner.partial(idx)
1124 }
1125
1126 pub fn gradient(&self) -> Vec<(String, T)> {
1132 crate::forward_tape::with_active_registry(|r| {
1133 let r =
1134 r.expect("NamedJet1Vec::gradient called outside a frozen NamedForwardTape scope");
1135 let n = r.len();
1136 let mut out = Vec::with_capacity(n);
1137 for (i, name) in r.iter().enumerate() {
1138 out.push((name.to_string(), self.inner.partial(i)));
1139 }
1140 out
1141 })
1142 }
1143
1144 #[inline]
1146 pub fn inner(&self) -> &Jet1Vec<T> {
1147 &self.inner
1148 }
1149
1150 #[inline]
1157 pub fn exp(&self) -> Self {
1158 Self {
1159 inner: self.inner.exp(),
1160 #[cfg(debug_assertions)]
1161 gen_id: self.gen_id,
1162 }
1163 }
1164
1165 #[inline]
1167 pub fn ln(&self) -> Self {
1168 Self {
1169 inner: self.inner.ln(),
1170 #[cfg(debug_assertions)]
1171 gen_id: self.gen_id,
1172 }
1173 }
1174
1175 #[inline]
1177 pub fn sqrt(&self) -> Self {
1178 Self {
1179 inner: self.inner.sqrt(),
1180 #[cfg(debug_assertions)]
1181 gen_id: self.gen_id,
1182 }
1183 }
1184
1185 #[inline]
1187 pub fn sin(&self) -> Self {
1188 Self {
1189 inner: self.inner.sin(),
1190 #[cfg(debug_assertions)]
1191 gen_id: self.gen_id,
1192 }
1193 }
1194
1195 #[inline]
1197 pub fn cos(&self) -> Self {
1198 Self {
1199 inner: self.inner.cos(),
1200 #[cfg(debug_assertions)]
1201 gen_id: self.gen_id,
1202 }
1203 }
1204
1205 #[inline]
1207 pub fn tan(&self) -> Self {
1208 Self {
1209 inner: self.inner.tan(),
1210 #[cfg(debug_assertions)]
1211 gen_id: self.gen_id,
1212 }
1213 }
1214
1215 #[inline]
1217 pub fn norm_cdf(&self) -> Self {
1218 Self {
1219 inner: self.inner.norm_cdf(),
1220 #[cfg(debug_assertions)]
1221 gen_id: self.gen_id,
1222 }
1223 }
1224
1225 #[inline]
1227 pub fn inv_norm_cdf(&self) -> Self {
1228 Self {
1229 inner: self.inner.inv_norm_cdf(),
1230 #[cfg(debug_assertions)]
1231 gen_id: self.gen_id,
1232 }
1233 }
1234
1235 #[inline]
1237 pub fn powf(&self, n: T) -> Self {
1238 Self {
1239 inner: self.inner.powf(n),
1240 #[cfg(debug_assertions)]
1241 gen_id: self.gen_id,
1242 }
1243 }
1244
1245 #[inline]
1247 pub fn powi(&self, n: i32) -> Self {
1248 Self {
1249 inner: self.inner.powi(n),
1250 #[cfg(debug_assertions)]
1251 gen_id: self.gen_id,
1252 }
1253 }
1254}
1255
1256impl<T: Passive> fmt::Display for NamedJet1Vec<T> {
1257 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1258 write!(f, "NamedJet1Vec({})", self.inner.real)
1259 }
1260}
1261
1262macro_rules! __named_dual_binop {
1273 ($trait:ident, $method:ident, $op:tt) => {
1274 impl<T: Passive> ::core::ops::$trait<NamedJet1Vec<T>> for NamedJet1Vec<T> {
1275 type Output = NamedJet1Vec<T>;
1276 #[inline]
1277 fn $method(self, rhs: NamedJet1Vec<T>) -> NamedJet1Vec<T> {
1278 #[cfg(debug_assertions)]
1279 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1280 NamedJet1Vec {
1281 inner: self.inner $op rhs.inner,
1282 #[cfg(debug_assertions)]
1283 gen_id: self.gen_id,
1284 }
1285 }
1286 }
1287 impl<T: Passive> ::core::ops::$trait<&NamedJet1Vec<T>> for &NamedJet1Vec<T> {
1288 type Output = NamedJet1Vec<T>;
1289 #[inline]
1290 fn $method(self, rhs: &NamedJet1Vec<T>) -> NamedJet1Vec<T> {
1291 #[cfg(debug_assertions)]
1292 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1293 NamedJet1Vec {
1294 inner: &self.inner $op &rhs.inner,
1295 #[cfg(debug_assertions)]
1296 gen_id: self.gen_id,
1297 }
1298 }
1299 }
1300 impl<T: Passive> ::core::ops::$trait<&NamedJet1Vec<T>> for NamedJet1Vec<T> {
1301 type Output = NamedJet1Vec<T>;
1302 #[inline]
1303 fn $method(self, rhs: &NamedJet1Vec<T>) -> NamedJet1Vec<T> {
1304 #[cfg(debug_assertions)]
1305 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1306 NamedJet1Vec {
1307 inner: self.inner $op &rhs.inner,
1308 #[cfg(debug_assertions)]
1309 gen_id: self.gen_id,
1310 }
1311 }
1312 }
1313 impl<T: Passive> ::core::ops::$trait<NamedJet1Vec<T>> for &NamedJet1Vec<T> {
1314 type Output = NamedJet1Vec<T>;
1315 #[inline]
1316 fn $method(self, rhs: NamedJet1Vec<T>) -> NamedJet1Vec<T> {
1317 #[cfg(debug_assertions)]
1318 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1319 NamedJet1Vec {
1320 inner: &self.inner $op rhs.inner,
1321 #[cfg(debug_assertions)]
1322 gen_id: self.gen_id,
1323 }
1324 }
1325 }
1326 impl<T: Passive> ::core::ops::$trait<T> for NamedJet1Vec<T> {
1327 type Output = NamedJet1Vec<T>;
1328 #[inline]
1329 fn $method(self, rhs: T) -> NamedJet1Vec<T> {
1330 NamedJet1Vec {
1331 inner: self.inner $op rhs,
1332 #[cfg(debug_assertions)]
1333 gen_id: self.gen_id,
1334 }
1335 }
1336 }
1337 impl<T: Passive> ::core::ops::$trait<T> for &NamedJet1Vec<T> {
1338 type Output = NamedJet1Vec<T>;
1339 #[inline]
1340 fn $method(self, rhs: T) -> NamedJet1Vec<T> {
1341 NamedJet1Vec {
1342 inner: &self.inner $op rhs,
1343 #[cfg(debug_assertions)]
1344 gen_id: self.gen_id,
1345 }
1346 }
1347 }
1348 };
1349}
1350
1351__named_dual_binop!(Add, add, +);
1352__named_dual_binop!(Sub, sub, -);
1353__named_dual_binop!(Mul, mul, *);
1354__named_dual_binop!(Div, div, /);
1355
1356impl<T: Passive> ::core::ops::AddAssign<NamedJet1Vec<T>> for NamedJet1Vec<T> {
1357 #[inline]
1358 fn add_assign(&mut self, rhs: NamedJet1Vec<T>) {
1359 #[cfg(debug_assertions)]
1360 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1361 self.inner += rhs.inner;
1362 }
1363}
1364impl<T: Passive> ::core::ops::AddAssign<&NamedJet1Vec<T>> for NamedJet1Vec<T> {
1365 #[inline]
1366 fn add_assign(&mut self, rhs: &NamedJet1Vec<T>) {
1367 #[cfg(debug_assertions)]
1368 crate::forward_tape::check_gen(self.gen_id, rhs.gen_id);
1369 self.inner += &rhs.inner;
1370 }
1371}
1372
1373impl<T: Passive> ::core::ops::Neg for NamedJet1Vec<T> {
1374 type Output = NamedJet1Vec<T>;
1375 #[inline]
1376 fn neg(self) -> NamedJet1Vec<T> {
1377 NamedJet1Vec {
1378 inner: -self.inner,
1379 #[cfg(debug_assertions)]
1380 gen_id: self.gen_id,
1381 }
1382 }
1383}
1384impl<T: Passive> ::core::ops::Neg for &NamedJet1Vec<T> {
1385 type Output = NamedJet1Vec<T>;
1386 #[inline]
1387 fn neg(self) -> NamedJet1Vec<T> {
1388 NamedJet1Vec {
1389 inner: -&self.inner,
1390 #[cfg(debug_assertions)]
1391 gen_id: self.gen_id,
1392 }
1393 }
1394}
1395
1396macro_rules! __named_dual_scalar_lhs {
1404 ($scalar:ty) => {
1405 impl ::core::ops::Add<NamedJet1Vec<$scalar>> for $scalar {
1406 type Output = NamedJet1Vec<$scalar>;
1407 #[inline]
1408 fn add(self, rhs: NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1409 NamedJet1Vec {
1410 inner: self + rhs.inner,
1411 #[cfg(debug_assertions)]
1412 gen_id: rhs.gen_id,
1413 }
1414 }
1415 }
1416 impl ::core::ops::Add<&NamedJet1Vec<$scalar>> for $scalar {
1417 type Output = NamedJet1Vec<$scalar>;
1418 #[inline]
1419 fn add(self, rhs: &NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1420 NamedJet1Vec {
1421 inner: self + &rhs.inner,
1422 #[cfg(debug_assertions)]
1423 gen_id: rhs.gen_id,
1424 }
1425 }
1426 }
1427 impl ::core::ops::Sub<NamedJet1Vec<$scalar>> for $scalar {
1428 type Output = NamedJet1Vec<$scalar>;
1429 #[inline]
1430 fn sub(self, rhs: NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1431 NamedJet1Vec {
1432 inner: self - rhs.inner,
1433 #[cfg(debug_assertions)]
1434 gen_id: rhs.gen_id,
1435 }
1436 }
1437 }
1438 impl ::core::ops::Sub<&NamedJet1Vec<$scalar>> for $scalar {
1439 type Output = NamedJet1Vec<$scalar>;
1440 #[inline]
1441 fn sub(self, rhs: &NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1442 NamedJet1Vec {
1443 inner: self - &rhs.inner,
1444 #[cfg(debug_assertions)]
1445 gen_id: rhs.gen_id,
1446 }
1447 }
1448 }
1449 impl ::core::ops::Mul<NamedJet1Vec<$scalar>> for $scalar {
1450 type Output = NamedJet1Vec<$scalar>;
1451 #[inline]
1452 fn mul(self, rhs: NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1453 NamedJet1Vec {
1454 inner: self * rhs.inner,
1455 #[cfg(debug_assertions)]
1456 gen_id: rhs.gen_id,
1457 }
1458 }
1459 }
1460 impl ::core::ops::Mul<&NamedJet1Vec<$scalar>> for $scalar {
1461 type Output = NamedJet1Vec<$scalar>;
1462 #[inline]
1463 fn mul(self, rhs: &NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1464 NamedJet1Vec {
1465 inner: self * &rhs.inner,
1466 #[cfg(debug_assertions)]
1467 gen_id: rhs.gen_id,
1468 }
1469 }
1470 }
1471 impl ::core::ops::Div<NamedJet1Vec<$scalar>> for $scalar {
1472 type Output = NamedJet1Vec<$scalar>;
1473 #[inline]
1474 fn div(self, rhs: NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1475 NamedJet1Vec {
1476 inner: self / rhs.inner,
1477 #[cfg(debug_assertions)]
1478 gen_id: rhs.gen_id,
1479 }
1480 }
1481 }
1482 impl ::core::ops::Div<&NamedJet1Vec<$scalar>> for $scalar {
1483 type Output = NamedJet1Vec<$scalar>;
1484 #[inline]
1485 fn div(self, rhs: &NamedJet1Vec<$scalar>) -> NamedJet1Vec<$scalar> {
1486 NamedJet1Vec {
1487 inner: self / &rhs.inner,
1488 #[cfg(debug_assertions)]
1489 gen_id: rhs.gen_id,
1490 }
1491 }
1492 }
1493 };
1494}
1495
1496__named_dual_scalar_lhs!(f64);
1497__named_dual_scalar_lhs!(f32);