Skip to content

fix: make tanh and powf work on the autodiff scalars - #234

Open
DPS0340 wants to merge 1 commit into
kmolan:mainfrom
DPS0340:fix/228-tanh-powf-autodiff
Open

fix: make tanh and powf work on the autodiff scalars#234
DPS0340 wants to merge 1 commit into
kmolan:mainfrom
DPS0340:fix/228-tanh-powf-autodiff

Conversation

@DPS0340

@DPS0340 DPS0340 commented Jul 27, 2026

Copy link
Copy Markdown

What & why

Fixes #228. Both Numeric defaults broke on exactly the types the crate exists for — f32/f64 override them with libm and were fine, Dual/HyperDual/Jet took the defaults and were not.

Reproduced on 5c45541 before changing anything:

tanh(800):    value=NaN  deriv=NaN
powf(0, 2):   value=0    deriv=NaN   x*x gives (0, 0)
powf(-2, 2):  value=NaN  deriv=NaN   x*x gives (4, -4)

Worth noting powf(0, 2) is the nastier of the two: the value comes out right and only the derivative is NaN, so it fails silently.

tanh

sinh/cosh is inf/inf once |x| passes the exp overflow point (~88 in f32, ~709 in f64), where the true value has simply settled at ±1. Now evaluated from exp of a non-positive argument, which underflows to zero instead of overflowing.

before  tanh(800) = (NaN, NaN)
after   tanh(800) = (1, 0)

One detail that cost me a revision: I first wrote it as signum(x) · f(|x|), which is tidier but wrong at the origin — signum(0) is zero, so it returned tanh'(0) = 0 instead of 1. The branch is on the sign for that reason, and there is a test pinning the origin.

powf

exp(n · ln self) is NaN for every negative base and produces a NaN derivative at zero. Non-positive bases are now handled explicitly:

  • self == 0: 0^0 = 1, 0^1 = 0 with derivative 1, 0^n = 0 with derivative 0 for integral n > 1.
  • self < 0 with integral n: real-valued. Magnitude from |self|, sign from the parity of n, so the chain rule stays intact — (-2)^3 gives (-8, 12).
  • self < 0 with fractional n: no real value, stays NaN as libm::pow does.

I did not add a to_i32 to the Numeric trait to route through powi; the parity is recovered with n - 2·floor(n/2), which stays inside the existing trait surface.

Tests

Three, in tests/suite/scalar.rs:

  • tanh_saturates_on_the_autodiff_scalars — the issue's case, plus the ordinary range and the origin so the saturating form cannot cost accuracy where the old one worked.
  • powf_handles_non_positive_bases — asserts x.powf(2.0) agrees with x * x in both value and derivative, which is the property the issue asks for.
  • primitive_tanh_and_powf_still_match_stdguard rail. The f32/f64 impls must keep matching std, so the defaults cannot have been "fixed" by loosening what the primitives do.

Verified they bite rather than assuming it. Reverting each default fails only its own test:

tanh reverted  -> tanh_saturates...              FAILED
powf reverted  -> powf_handles_non_positive...   FAILED
both reverted  -> primitive_tanh_and_powf...     ok      (guard rail)

Checklist

  • cargo test (625 passed / 0 failed) + cargo clippy --all-targets clean locally
  • cargo fmt --check clean
  • CHANGELOG entry under ## [Unreleased]### Fixed
  • No unwrap/expect/panic on library paths

One thing I deliberately left alone

sinh and cosh return inf for large |x|, which is correct, but their derivatives come back NaN. That is not these defaults — it is Dual's product rule forming inf * 0 against the HALF constant:

sinh(800): value=inf  deriv=NaN     (true value is inf)

Different cause, different fix, and it would widen this PR past the issue. Happy to open a separate issue for it if you want it tracked.

Both Numeric defaults broke on exactly the types the crate exists for.
f32/f64 override them with libm and were fine; Dual, HyperDual and Jet
took the defaults and were not.

tanh was sinh/cosh, so once |x| passed the exp overflow point (about 88
in f32, 709 in f64) it computed inf/inf — NaN — where the true value has
simply settled at plus or minus one. It now evaluates from exp of a
non-positive argument, which underflows to zero instead of overflowing:

  before  tanh(800) = (NaN, NaN)
  after   tanh(800) = (1, 0)

The branch is on the sign rather than a multiply by signum, because
signum(0) is zero and would wipe out the derivative at the origin where
tanh'(0) = 1. Verified that case explicitly.

powf was exp(n · ln self), which is NaN for every negative base, and at
zero returns the right value with a NaN derivative — a silent failure
rather than an obvious one:

  before  powf(0, 2)  = (0,   NaN)     x*x gives (0, 0)
  before  powf(-2, 2) = (NaN, NaN)     x*x gives (4, -4)
  after   both agree with x*x

Non-positive bases are handled explicitly. An integral exponent on a
negative base is real-valued, so the magnitude comes from |self| and the
sign from the parity of n, which keeps the chain rule intact. A negative
base with a genuinely fractional exponent has no real value and stays
NaN, matching libm::pow.

Tests assert the issue's cases plus a guard rail: the f32/f64 impls must
still match std, so the defaults cannot have been 'fixed' by loosening
what the primitives do. Reverting either default fails only its own test
and leaves the guard rail passing.

One thing this does not change: sinh and cosh return inf for large |x|,
which is correct, but their derivatives are NaN because Dual's product
rule forms inf * 0 against the HALF constant. That is a property of the
Dual arithmetic rather than these defaults, so it is left alone here.

cargo test 625 passed / 0 failed. clippy --all-targets clean.

Fixes kmolan#228.
@DPS0340
DPS0340 requested a review from kmolan as a code owner July 27, 2026 11:16

@github-actions github-actions Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the first PR! 🎉 You'll get a first review within 48 h. CI runs the 5-target matrix — a budget-gate failure isn't your fault, see ci/README.md.

Comment thread CHANGELOG.md
for every negative base and gives a NaN derivative at zero. `f32`/`f64` override both
with `libm` and were unaffected; `Dual`, `HyperDual` and `Jet` took the defaults. `tanh`
now saturates to ±1 with a vanishing derivative, and `powf` handles zero and negative
bases with integral exponents, keeping the derivatives correct. @DPS0340 (#228)

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

House style is the PR number, not the issue. This should be (#234) rather than (#228).

}
if n > Self::ONE && n == n.floor() {
return self * self * (n - Self::ONE);
}

@rtmongold rtmongold Jul 28, 2026

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This zero-base path claims to match powi for integral n, but
self * self * (n - ONE) only preserves first-order Dual derivatives.

On HyperDual / Jet it is wrong for n ≠ 2 (e.g. HyperDual 0^3 gives
f'' = 4 via powf, 0 via powi), and replaces a loud NaN with a
silent wrong higher derivative.

Please fix this for higher-order types. Prefer routing all integral
exponents through powi (zero and negative bases) — that removes both
hand-rolled special cases and the Dual-only formula in one go, even if it
needs a small to_i32 / to_f64 on Numeric. Don’t claim powi
equivalence for self * self * (n - ONE).

Comment on lines +398 to +435
fn powf_handles_non_positive_bases() {
for &x in &[0.0_f64, -2.0, -0.5, 2.0] {
let squared = Dual::variable(x).powf(Dual::constant(2.0));
let by_multiplication = {
let d = Dual::variable(x);
d * d
};
assert!(
f64::abs(squared.value - by_multiplication.value) < TOL,
"({x})^2 value: powf gave {}, x*x gave {}",
squared.value,
by_multiplication.value
);
assert!(
f64::abs(squared.deriv - by_multiplication.deriv) < TOL,
"({x})^2 derivative: powf gave {}, x*x gave {}",
squared.deriv,
by_multiplication.deriv
);
}

// Odd integer exponent on a negative base keeps the sign.
let cubed = Dual::variable(-2.0_f64).powf(Dual::constant(3.0));
assert!(f64::abs(cubed.value + 8.0) < TOL);
assert!(f64::abs(cubed.deriv - 12.0) < TOL); // 3·(−2)² = 12

// n = 1 and n = 0 at the origin, where the exp/ln form divides by -inf.
let linear = Dual::variable(0.0_f64).powf(Dual::constant(1.0));
assert!(f64::abs(linear.value) < TOL);
assert!(f64::abs(linear.deriv - 1.0) < TOL);
let unity = Dual::variable(0.0_f64).powf(Dual::constant(0.0));
assert!(f64::abs(unity.value - 1.0) < TOL);

// A negative base with a genuinely fractional exponent has no real
// value and must stay NaN rather than inventing one.
let undefined = Dual::variable(-2.0_f64).powf(Dual::constant(0.5));
assert!(undefined.value.is_nan());
}

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These only exercise Dual, and at zero only n ∈ {0, 1, 2}. That is exactly where the bad self * self * (n - 1) formula still looks correct.
Please add HyperDual (and ideally Jet) cases for at least:

  • 0^2 vs x * x (value + higher parts)
  • 0^3 / 0^4 vs powi — this is what catches the bug
  • (-2)^2 / (-2)^3 — negative-base path looks good on Dual; worth locking for HyperDual too

Changelog says Dual / HyperDual / Jet all took the broken defaults — tests should match that claim.

Copy link
Copy Markdown
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@rtmongold I was actually thinking of removing the CHANGELOG requirement from PRs. As we scale up, its just going to become a source of merge conflicts and pain. I anyway have to do it myself in bulk before bumping the release version. Thoughts?

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@kmolan Funny, I was thinking the same the other day — there’s already a lot of overlap on that file. One option that avoids the conflicts: a changelog.d/ directory where each PR adds its own file (e.g. 234.added.md) with the description inside. At release we assemble those into CHANGELOG.md.

I’ll open a new issue so we can compare approaches (fragments vs drafting from merged PRs at release time, etc.).

///
/// - `self == 0`: the value is `0` for `n > 0` and the derivative is `n · 0^(n−1)`,
/// which is `0` for `n > 1`, `1` for `n == 1`, and divergent below that. Computed
/// directly rather than through `ln 0 = -inf`.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same issue as the implementation: “derivative is 0 for n > 1” is true for Dual, but HyperDual/Jet still carry non-zero higher derivatives for 0^n (e.g. f''(0) for n = 2, third-order for n = 3). Worth tightening once the code path is fixed.

@rtmongold

Copy link
Copy Markdown
Collaborator

Will need to also rebase to resolve CHANGELOG.md conflicts, but can wait until we get everything else fixed.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Fix tanh and powf on the autodiff scalars

3 participants