georust / geographiclib-rs

A port of geographiclib in Rust.
MIT License
41 stars 9 forks source link

consider simplifying polyval #6

Closed stonylohr closed 2 years ago

stonylohr commented 4 years ago

A few changes to consider for geomath::polyval.

Consider removing the "s" argument, instead have the calling function slice the "p" argument any time "s" might be non-zero. I think this is more idiomatic rust, and allows simplifying polyval's loop. That in turn allows for dropping the shadowing of the "n" argument (see #3) without needing to even declare "n" as a mut argument. polyval's definition then becomes a more streamlined...

pub fn polyval(n: i64, p: &[f64], x: f64) -> f64 {
    let mut y = if n < 0 { 0.0 } else { p[0] };
    for i in 1..=n {
        y = y * x + p[i as usize];
    }
    y
}

Calls to polyval would need to be modified to pass a slice starting at index "s", rather than starting at 0 and passing "s".

Example 1 (called with fixed "s" of 0): Before: polyval(m, &COEFF, 0, sq(eps)) After: polyval(m, &COEFF, sq(eps))

Example 2 (called with variable "s" value): Before: geomath::polyval(m, &COEFF_A3, o as usize, _n) After: geomath::polyval(m, &COEFF_A3[o as usize..], _n) (Note that calling contexts like this would be further improved by declaring "o" as usize in the first place, to avoid the need to cast. This change should be straightforward in most cases.)

One note of caution on this change is that C++ also supports similar array offset mechanics, yet Charles Karney chose not to use them in the original C++. I suspect this is mostly a matter of style, but it could hint that there's a slight performance benefit to passing the "s" index, in which case that difference is likely (but not certain) to carry over into rust.

Depending on your change tolerance, I would further consider changing polyval's "n" argument from i64 to isize. The difference is academic for 64 bit platforms, but for values that are primarily used as an index into arrays, it's generally preferable to stick with the architecture's native size. This would involve a change at each calling context. While the changes should be fairly simple, the benefit is also modest.

stonylohr commented 4 years ago

I take back my cautionary note about the C++ version passing the index. Apparently my memory failed me on that point. It appears to use offset mechanics similar to those I proposed.

I'm also reminded by looking back at the C++ code that it might be worth adding a comment to consider using fma in polyval at some point in the future, as hardware support for it becomes common enough (apparently using it without hardware support results in a performance hit).

michaelkirk commented 4 years ago

Interesting, sounds good to me!

The original transcription was based on the python implementation -so I assume that explains the divergence from the C++ implementation.

Note: I edited your post to add some "```" in order to format the code blocks more readably.

michaelkirk commented 2 years ago

Fixed by #39