Mercurial > repos > non-riemannian-opt / file comparison
comparison: src/cylinder.rs
src/cylinder.rs
- changeset 41
- 56d609b70a3b
- parent 40
- 1d865db9d3e4
- child 46
- 90cc221eb52b
equal
deleted
inserted
replaced
| 40:1d865db9d3e4 | 41:56d609b70a3b |
|---|---|
| 499 let φ = self.side_cap_crossing(side, cap, z); | 499 let φ = self.side_cap_crossing(side, cap, z); |
| 500 let capedge = CapPoint{ angle : φ, r }; | 500 let capedge = CapPoint{ angle : φ, r }; |
| 501 let sideedge = SidePoint{ angle : φ, z }; | 501 let sideedge = SidePoint{ angle : φ, z }; |
| 502 let t1 = cap.log(&capedge); | 502 let t1 = cap.log(&capedge); |
| 503 let t2 = sideedge.log(r, side); | 503 let t2 = sideedge.log(r, side); |
| 504 // Either option should give the same result, but the first one avoids division. | 504 // Either option should give the same result. |
| 505 t1 + capedge.tangent_from_side(top, t2) | 505 // The first one avoids division, but seems numerically more unstable due to |
| 506 // let n = t1.norm(L2); | 506 // sines and cosines. |
| 507 // (t1/n)*(n + t2.norm(L2)) | 507 //t1 + capedge.tangent_from_side(top, t2) |
| 508 let n = t1.norm(L2); | |
| 509 if n < f64::EPSILON { | |
| 510 capedge.tangent_from_side(top, t2) | |
| 511 } else { | |
| 512 (t1/n)*(n + t2.norm(L2)) | |
| 513 } | |
| 508 } | 514 } |
| 509 } | 515 } |
| 510 | 516 |
| 511 fn side_cap_log( | 517 fn side_cap_log( |
| 512 &self, | 518 &self, |
| 526 let φ = self.side_cap_crossing(side, cap, z); | 532 let φ = self.side_cap_crossing(side, cap, z); |
| 527 let capedge = CapPoint{ angle : φ, r }; | 533 let capedge = CapPoint{ angle : φ, r }; |
| 528 let sideedge = SidePoint{ angle : φ, z }; | 534 let sideedge = SidePoint{ angle : φ, z }; |
| 529 let t1 = side.log(r, &sideedge); | 535 let t1 = side.log(r, &sideedge); |
| 530 let t2 = capedge.log(cap); | 536 let t2 = capedge.log(cap); |
| 531 // Either option should give the same result, but the first one avoids division. | 537 // Either option should give the same result. |
| 532 t1 + capedge.tangent_to_side(top, t2) | 538 // The first one avoids division, but seems numerically more unstable due to |
| 533 // let n = t1.norm(L2); | 539 // sines and cosines. |
| 534 // (t1/n)*(n + t2.norm(L2)) | 540 // t1 + capedge.tangent_to_side(top, t2) |
| 541 let n = t1.norm(L2); | |
| 542 if n < f64::EPSILON { | |
| 543 capedge.tangent_to_side(top, t2) | |
| 544 } else { | |
| 545 (t1/n)*(n + t2.norm(L2)) | |
| 546 } | |
| 535 } | 547 } |
| 536 } | 548 } |
| 537 | 549 |
| 538 fn side_cap_side_log( | 550 fn side_cap_side_log( |
| 539 &self, | 551 &self, |
| 561 let capedge2 = CapPoint{ angle : φ2, r }; | 573 let capedge2 = CapPoint{ angle : φ2, r }; |
| 562 let sideedge2 = SidePoint{ angle : φ2, z }; | 574 let sideedge2 = SidePoint{ angle : φ2, z }; |
| 563 let t1 = side1.log(r, &sideedge1); | 575 let t1 = side1.log(r, &sideedge1); |
| 564 let t2 = capedge1.log(&capedge2); | 576 let t2 = capedge1.log(&capedge2); |
| 565 let t3 = sideedge2.log(r, &side2); | 577 let t3 = sideedge2.log(r, &side2); |
| 566 // Any option should give the same result, but the first one avoids division. | 578 // Either option should give the same result. |
| 579 // The first one avoids division, but seems numerically more unstable due to | |
| 580 // sines and cosines. | |
| 581 // | |
| 567 // t1 + capedge1.tangent_to_side(top, t2 + capedge2.tangent_from_side(top, t3)) | 582 // t1 + capedge1.tangent_to_side(top, t2 + capedge2.tangent_from_side(top, t3)) |
| 568 // | |
| 569 // let n = t2.norm(L2); | |
| 570 // t1 + capedge1.tangent_to_side(top, (t2/n)*(n + t3.norm(L2))) | |
| 571 // | |
| 572 let n = t1.norm(L2); | 583 let n = t1.norm(L2); |
| 573 (t1/n)*(n + t2.norm(L2) + t3.norm(L2)) | 584 if n < f64::EPSILON { |
| 574 // | 585 let n2 = t2.norm(L2); |
| 575 // let n = t1.norm(L2); | 586 capedge1.tangent_to_side(top, |
| 576 // let t23 = t2 + capedge2.tangent_from_side(top, t3); | 587 if n2 < f64::EPSILON { |
| 577 // (t1/n)*(n + t23.norm(L2)) | 588 capedge2.tangent_from_side(top, t3) |
| 589 } else { | |
| 590 (t2/n2)*(n2 + t3.norm(L2)) | |
| 591 } | |
| 592 ) | |
| 593 } else { | |
| 594 (t1/n)*(n + t2.norm(L2) + t3.norm(L2)) | |
| 595 } | |
| 578 } | 596 } |
| 579 } | 597 } |
| 580 | 598 |
| 581 fn cap_side_cap_log( | 599 fn cap_side_cap_log( |
| 582 &self, | 600 &self, |
| 604 let sideedge2 = SidePoint{ angle : φ2, z : z2}; | 622 let sideedge2 = SidePoint{ angle : φ2, z : z2}; |
| 605 let capedge2 = CapPoint{ angle : φ2, r }; | 623 let capedge2 = CapPoint{ angle : φ2, r }; |
| 606 let t1 = cap1.log(&capedge1); | 624 let t1 = cap1.log(&capedge1); |
| 607 let t2 = sideedge1.log(r, &sideedge2); | 625 let t2 = sideedge1.log(r, &sideedge2); |
| 608 let t3 = capedge2.log(cap2); | 626 let t3 = capedge2.log(cap2); |
| 609 // Either option should give the same result, but the first one avoids division. | 627 // Either option should give the same result. |
| 610 t1 + capedge1.tangent_from_side(top1, t2 + capedge2.tangent_to_side(top2, t3)) | 628 // The first one avoids division, but seems numerically more unstable due to |
| 611 //let n = t1.norm(L2); | 629 // sines and cosines. |
| 612 //(t1/n)*(n + t2.norm(L2) + t3.norm(L2)) | 630 // t1 + capedge1.tangent_from_side(top1, t2 + capedge2.tangent_to_side(top2, t3)) |
| 631 let n = t1.norm(L2); | |
| 632 if n < f64::EPSILON { | |
| 633 let n2 = t2.norm(L2); | |
| 634 capedge1.tangent_from_side(top1, | |
| 635 if n2 < f64::EPSILON { | |
| 636 capedge2.tangent_to_side(top2, t3) | |
| 637 } else { | |
| 638 (t2/n2)*(n2 + t3.norm(L2)) | |
| 639 } | |
| 640 ) | |
| 641 } else { | |
| 642 (t1/n)*(n + t2.norm(L2) + t3.norm(L2)) | |
| 643 } | |
| 613 } | 644 } |
| 614 } | 645 } |
| 615 | 646 |
| 616 /// Calculates both the logarithmic map and distance to another point | 647 /// Calculates both the logarithmic map and distance to another point |
| 617 fn log_dist(&self, source : &Point, destination : &Point) -> (Tangent, f64) { | 648 fn log_dist(&self, source : &Point, destination : &Point) -> (Tangent, f64) { |
| 849 | 880 |
| 850 #[cfg(test)] | 881 #[cfg(test)] |
| 851 mod tests { | 882 mod tests { |
| 852 use super::*; | 883 use super::*; |
| 853 | 884 |
| 854 static TOL : f64 = 1e-8; | 885 static TOL : f64 = 1e-7; |
| 855 | 886 |
| 856 macro_rules! check_distance { | 887 macro_rules! check_distance { |
| 857 ($distance : expr, $expected : expr) => { | 888 ($distance : expr, $expected : expr) => { |
| 858 assert!( | 889 assert!( |
| 859 ($distance-$expected).abs() < TOL, | 890 ($distance-$expected).abs() < TOL, |
