From a4fb5ddd38eeb7cb8c59725307ca236ee6d2bf69 Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 17:40:56 +0000 Subject: [PATCH 1/9] Adopt the check_dist gate: commit Cargo.lock and the pinned toolchain and clippy config, and deny the lints it requires --- .clippy.toml | 3 - .gitignore | 1 - Cargo.lock | 697 ++++++++++++++++++++++++++++++++++++++++++++ Cargo.toml | 7 + clippy.toml | 11 + rust-toolchain.toml | 5 + 6 files changed, 720 insertions(+), 4 deletions(-) delete mode 100644 .clippy.toml create mode 100644 Cargo.lock create mode 100644 clippy.toml create mode 100644 rust-toolchain.toml diff --git a/.clippy.toml b/.clippy.toml deleted file mode 100644 index c4f4673..0000000 --- a/.clippy.toml +++ /dev/null @@ -1,3 +0,0 @@ -cognitive-complexity-threshold = 25 -too-many-arguments-threshold = 7 -type-complexity-threshold = 250 diff --git a/.gitignore b/.gitignore index 636c941..9f338bb 100644 --- a/.gitignore +++ b/.gitignore @@ -1,6 +1,5 @@ # Rust artifacts target/ -Cargo.lock # IDE .idea/ diff --git a/Cargo.lock b/Cargo.lock new file mode 100644 index 0000000..77f975d --- /dev/null +++ b/Cargo.lock @@ -0,0 +1,697 @@ +# This file is automatically @generated by Cargo. +# It is not intended for manual editing. +version = 4 + +[[package]] +name = "aho-corasick" +version = "1.1.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ddd31a130427c27518df266943a5308ed92d4b226cc639f5a8f1002816174301" +dependencies = [ + "memchr", +] + +[[package]] +name = "alloca" +version = "0.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e5a7d05ea6aea7e9e64d25b9156ba2fee3fdd659e34e41063cd2fc7cd020d7f4" +dependencies = [ + "cc", +] + +[[package]] +name = "anes" +version = "0.1.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "4b46cbb362ab8752921c97e041f5e366ee6297bd428a31275b9fcf1e380f7299" + +[[package]] +name = "anstyle" +version = "1.0.14" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "940b3a0ca603d1eade50a4846a2afffd5ef57a9feac2c0e2ec2e14f9ead76000" + +[[package]] +name = "autocfg" +version = "1.5.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f2032f911046de80f0a198e0901378627c33f59ea0ac00e363d481118bd70a53" + +[[package]] +name = "az" +version = "1.3.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "be5eb007b7cacc6c660343e96f650fedf4b5a77512399eb952ca6642cf8d13f7" + +[[package]] +name = "bumpalo" +version = "3.20.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "72f5acc6cb2ba439de613abc23857ec3d78374d8ed5ac84e9d11336e87da8649" + +[[package]] +name = "bytemuck" +version = "1.25.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "c8efb64bd706a16a1bdde310ae86b351e4d21550d98d056f22f8a7f7a2183fec" + +[[package]] +name = "cast" +version = "0.3.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "37b2a672a2cb129a2e41c10b1224bb368f9f37a2b16b612598138befd7b37eb5" + +[[package]] +name = "cc" +version = "1.2.66" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f5d6cac793997bd970000024b2934968efe83b382de4fdcf4fcb46b6ee4ad996" +dependencies = [ + "find-msvc-tools", + "shlex", +] + +[[package]] +name = "cfg-if" +version = "1.0.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9330f8b2ff13f34540b44e946ef35111825727b38d33286ef986142615121801" + +[[package]] +name = "ciborium" +version = "0.2.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "42e69ffd6f0917f5c029256a24d0161db17cea3997d185db0d35926308770f0e" +dependencies = [ + "ciborium-io", + "ciborium-ll", + "serde", +] + +[[package]] +name = "ciborium-io" +version = "0.2.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "05afea1e0a06c9be33d539b876f1ce3692f4afea2cb41f740e7743225ed1c757" + +[[package]] +name = "ciborium-ll" +version = "0.2.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "57663b653d948a338bfb3eeba9bb2fd5fcfaecb9e199e87e1eda4d9e8b240fd9" +dependencies = [ + "ciborium-io", + "half", +] + +[[package]] +name = "clap" +version = "4.6.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1ddb117e43bbf7dacf0a4190fef4d345b9bad68dfc649cb349e7d17d28428e51" +dependencies = [ + "clap_builder", +] + +[[package]] +name = "clap_builder" +version = "4.6.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "714a53001bf66416adb0e2ef5ac857140e7dc3a0c48fb28b2f10762fc4b5069f" +dependencies = [ + "anstyle", + "clap_lex", +] + +[[package]] +name = "clap_lex" +version = "1.1.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "c8d4a3bb8b1e0c1050499d1815f5ab16d04f0959b233085fb31653fbfc9d98f9" + +[[package]] +name = "criterion" +version = "0.8.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "950046b2aa2492f9a536f5f4f9a3de7b9e2476e575e05bd6c333371add4d98f3" +dependencies = [ + "alloca", + "anes", + "cast", + "ciborium", + "clap", + "criterion-plot", + "itertools", + "num-traits", + "oorandom", + "page_size", + "plotters", + "rayon", + "regex", + "serde", + "serde_json", + "tinytemplate", + "walkdir", +] + +[[package]] +name = "criterion-plot" +version = "0.8.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d8d80a2f4f5b554395e47b5d8305bc3d27813bacb73493eb1001e8f76dae29ea" +dependencies = [ + "cast", + "itertools", +] + +[[package]] +name = "crossbeam-deque" +version = "0.8.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9dd111b7b7f7d55b72c0a6ae361660ee5853c9af73f70c3c2ef6858b950e2e51" +dependencies = [ + "crossbeam-epoch", + "crossbeam-utils", +] + +[[package]] +name = "crossbeam-epoch" +version = "0.9.18" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5b82ac4a3c2ca9c3460964f020e1402edd5753411d7737aa39c3714ad1b5420e" +dependencies = [ + "crossbeam-utils", +] + +[[package]] +name = "crossbeam-utils" +version = "0.8.21" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d0a5c400df2834b80a4c3327b3aad3a4c4cd4de0629063962b03235697506a28" + +[[package]] +name = "crunchy" +version = "0.2.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "460fbee9c2c2f33933d720630a6a0bac33ba7053db5344fac858d4b8952d77d5" + +[[package]] +name = "either" +version = "1.16.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "91622ff5e7162018101f2fea40d6ebf4a78bbe5a49736a2020649edf9693679e" + +[[package]] +name = "find-msvc-tools" +version = "0.1.9" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5baebc0774151f905a1a2cc41989300b1e6fbb29aff0ceffa1064fdd3088d582" + +[[package]] +name = "fixed" +version = "1.31.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9af2cbf772fa6d1c11358f92ef554cb6b386201210bcf0e91fb7fba8a907fb40" +dependencies = [ + "az", + "bytemuck", + "half", + "typenum", +] + +[[package]] +name = "fixed_analytics" +version = "3.1.0" +dependencies = [ + "criterion", + "fixed", + "no-panic", +] + +[[package]] +name = "futures-core" +version = "0.3.32" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "7e3450815272ef58cec6d564423f6e755e25379b217b0bc688e295ba24df6b1d" + +[[package]] +name = "futures-task" +version = "0.3.32" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "037711b3d59c33004d3856fbdc83b99d4ff37a24768fa1be9ce3538a1cde4393" + +[[package]] +name = "futures-util" +version = "0.3.32" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "389ca41296e6190b48053de0321d02a77f32f8a5d2461dd38762c0593805c6d6" +dependencies = [ + "futures-core", + "futures-task", + "pin-project-lite", + "slab", +] + +[[package]] +name = "half" +version = "2.7.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "6ea2d84b969582b4b1864a92dc5d27cd2b77b622a8d79306834f1be5ba20d84b" +dependencies = [ + "cfg-if", + "crunchy", + "zerocopy", +] + +[[package]] +name = "itertools" +version = "0.13.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "413ee7dfc52ee1a4949ceeb7dbc8a33f2d6c088194d9f922fb8318faf1f01186" +dependencies = [ + "either", +] + +[[package]] +name = "itoa" +version = "1.0.18" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "8f42a60cbdf9a97f5d2305f08a87dc4e09308d1276d28c869c684d7777685682" + +[[package]] +name = "js-sys" +version = "0.3.103" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "53b44bfcdb3f8d5837a46dae1ca9660a837176eee74a28b229bc626816589102" +dependencies = [ + "cfg-if", + "futures-util", + "wasm-bindgen", +] + +[[package]] +name = "libc" +version = "0.2.186" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "68ab91017fe16c622486840e4c83c9a37afeff978bd239b5293d61ece587de66" + +[[package]] +name = "memchr" +version = "2.8.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "88904434abc2901f197fe8cc55f0445e7ded921dba5911dad2e2b39b48e663c4" + +[[package]] +name = "no-panic" +version = "0.1.37" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "fc80370b4544f28ffa317e3c3474ee3ecbfe269196c01ae657d9f837a7d944a1" +dependencies = [ + "proc-macro2", + "quote", + "syn 3.0.4", +] + +[[package]] +name = "num-traits" +version = "0.2.19" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "071dfc062690e90b734c0b2273ce72ad0ffa95f0c74596bc250dcfd960262841" +dependencies = [ + "autocfg", +] + +[[package]] +name = "once_cell" +version = "1.21.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9f7c3e4beb33f85d45ae3e3a1792185706c8e16d043238c593331cc7cd313b50" + +[[package]] +name = "oorandom" +version = "11.1.5" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d6790f58c7ff633d8771f42965289203411a5e5c68388703c06e14f24770b41e" + +[[package]] +name = "page_size" +version = "0.6.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "30d5b2194ed13191c1999ae0704b7839fb18384fa22e49b57eeaa97d79ce40da" +dependencies = [ + "libc", + "winapi", +] + +[[package]] +name = "pin-project-lite" +version = "0.2.17" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "a89322df9ebe1c1578d689c92318e070967d1042b512afbe49518723f4e6d5cd" + +[[package]] +name = "plotters" +version = "0.3.7" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5aeb6f403d7a4911efb1e33402027fc44f29b5bf6def3effcc22d7bb75f2b747" +dependencies = [ + "num-traits", + "plotters-backend", + "plotters-svg", + "wasm-bindgen", + "web-sys", +] + +[[package]] +name = "plotters-backend" +version = "0.3.7" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "df42e13c12958a16b3f7f4386b9ab1f3e7933914ecea48da7139435263a4172a" + +[[package]] +name = "plotters-svg" +version = "0.3.7" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "51bae2ac328883f7acdfea3d66a7c35751187f870bc81f94563733a154d7a670" +dependencies = [ + "plotters-backend", +] + +[[package]] +name = "proc-macro2" +version = "1.0.106" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "8fd00f0bb2e90d81d1044c2b32617f68fcb9fa3bb7640c23e9c748e53fb30934" +dependencies = [ + "unicode-ident", +] + +[[package]] +name = "quote" +version = "1.0.46" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "dfbc457d0c7a0759a614551b11a6409e5951f6c7537be1f1b7682b9ae9230368" +dependencies = [ + "proc-macro2", +] + +[[package]] +name = "rayon" +version = "1.12.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "fb39b166781f92d482534ef4b4b1b2568f42613b53e5b6c160e24cfbfa30926d" +dependencies = [ + "either", + "rayon-core", +] + +[[package]] +name = "rayon-core" +version = "1.13.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "22e18b0f0062d30d4230b2e85ff77fdfe4326feb054b9783a3460d8435c8ab91" +dependencies = [ + "crossbeam-deque", + "crossbeam-utils", +] + +[[package]] +name = "regex" +version = "1.12.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f1292b7759ae1cb9ec195452d1390a074f0cd8541ab7a5a8c31cd6db45d4a6ba" +dependencies = [ + "aho-corasick", + "memchr", + "regex-automata", + "regex-syntax", +] + +[[package]] +name = "regex-automata" +version = "0.4.14" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "6e1dd4122fc1595e8162618945476892eefca7b88c52820e74af6262213cae8f" +dependencies = [ + "aho-corasick", + "memchr", + "regex-syntax", +] + +[[package]] +name = "regex-syntax" +version = "0.8.11" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d6f6ff9a378485b298a5286656da665ba74413d36db0979633275d2e708145d4" + +[[package]] +name = "rustversion" +version = "1.0.22" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b39cdef0fa800fc44525c84ccb54a029961a8215f9619753635a9c0d2538d46d" + +[[package]] +name = "same-file" +version = "1.0.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "93fc1dc3aaa9bfed95e02e6eadabb4baf7e3078b0bd1b4d7b6b0b68378900502" +dependencies = [ + "winapi-util", +] + +[[package]] +name = "serde" +version = "1.0.228" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9a8e94ea7f378bd32cbbd37198a4a91436180c5bb472411e48b5ec2e2124ae9e" +dependencies = [ + "serde_core", + "serde_derive", +] + +[[package]] +name = "serde_core" +version = "1.0.228" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "41d385c7d4ca58e59fc732af25c3983b67ac852c1a25000afe1175de458b67ad" +dependencies = [ + "serde_derive", +] + +[[package]] +name = "serde_derive" +version = "1.0.228" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d540f220d3187173da220f885ab66608367b6574e925011a9353e4badda91d79" +dependencies = [ + "proc-macro2", + "quote", + "syn 2.0.118", +] + +[[package]] +name = "serde_json" +version = "1.0.150" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e8014e44b4736ed0538adeecded0fce2a272f22dc9578a7eb6b2d9993c74cfb9" +dependencies = [ + "itoa", + "memchr", + "serde", + "serde_core", + "zmij", +] + +[[package]] +name = "shlex" +version = "2.0.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f8fadd59c855ef2080decdef8ff161eb6661b86933c9d82e5ba29dc602a55aba" + +[[package]] +name = "slab" +version = "0.4.12" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "0c790de23124f9ab44544d7ac05d60440adc586479ce501c1d6d7da3cd8c9cf5" + +[[package]] +name = "syn" +version = "2.0.118" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1b9ae57f904213ebb649ce6895b8a66c66f0203b9319718f69a5612a065b1422" +dependencies = [ + "proc-macro2", + "quote", + "unicode-ident", +] + +[[package]] +name = "syn" +version = "3.0.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e6275cddf4610d1775e6d1fe9469b2e77d0f39fd98fb7450901b821e0c53649f" +dependencies = [ + "proc-macro2", + "quote", + "unicode-ident", +] + +[[package]] +name = "tinytemplate" +version = "1.2.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "be4d6b5f19ff7664e8c98d03e2139cb510db9b0a60b55f8e8709b689d939b6bc" +dependencies = [ + "serde", + "serde_json", +] + +[[package]] +name = "typenum" +version = "1.20.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b6f5e870be6c3b371b77fe0ee0bafb859fa4964b4404c27de1d380043c4dda20" + +[[package]] +name = "unicode-ident" +version = "1.0.24" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e6e4313cd5fcd3dad5cafa179702e2b244f760991f45397d14d4ebf38247da75" + +[[package]] +name = "walkdir" +version = "2.5.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "29790946404f91d9c5d06f9874efddea1dc06c5efe94541a7d6863108e3a5e4b" +dependencies = [ + "same-file", + "winapi-util", +] + +[[package]] +name = "wasm-bindgen" +version = "0.2.126" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "4b067c0c11094aef6b7a801c1e34a26affafdf3d051dba08456b868789aaf9a4" +dependencies = [ + "cfg-if", + "once_cell", + "rustversion", + "wasm-bindgen-macro", + "wasm-bindgen-shared", +] + +[[package]] +name = "wasm-bindgen-macro" +version = "0.2.126" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "167ce5e579f6bcf889c4f7175a8a5a585de84e8ff93976ce393efa5f2837aab1" +dependencies = [ + "quote", + "wasm-bindgen-macro-support", +] + +[[package]] +name = "wasm-bindgen-macro-support" +version = "0.2.126" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f3997c7839262f4ef12cf90b818d6340c18e80f263f1a94bf157d0ec4420380e" +dependencies = [ + "bumpalo", + "proc-macro2", + "quote", + "syn 2.0.118", + "wasm-bindgen-shared", +] + +[[package]] +name = "wasm-bindgen-shared" +version = "0.2.126" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "dc1b4cb0cc549fcf58d7dfc081778139b3d283a081644e833e84682ad71cea24" +dependencies = [ + "unicode-ident", +] + +[[package]] +name = "web-sys" +version = "0.3.103" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "8622dcb61c0bcc9fffa6938bed81210af2da9a7e4a1a834b2e37a59b6dfb6141" +dependencies = [ + "js-sys", + "wasm-bindgen", +] + +[[package]] +name = "winapi" +version = "0.3.9" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5c839a674fcd7a98952e593242ea400abe93992746761e38641405d28b00f419" +dependencies = [ + "winapi-i686-pc-windows-gnu", + "winapi-x86_64-pc-windows-gnu", +] + +[[package]] +name = "winapi-i686-pc-windows-gnu" +version = "0.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ac3b87c63620426dd9b991e5ce0329eff545bccbbb34f3be09ff6fb6ab51b7b6" + +[[package]] +name = "winapi-util" +version = "0.1.11" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "c2a7b1c03c876122aa43f3020e6c3c3ee5c05081c9a00739faf7503aeba10d22" +dependencies = [ + "windows-sys", +] + +[[package]] +name = "winapi-x86_64-pc-windows-gnu" +version = "0.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "712e227841d057c1ee1cd2fb22fa7e5a5461ae8e48fa2ca79ec42cfc1931183f" + +[[package]] +name = "windows-link" +version = "0.2.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f0805222e57f7521d6a62e36fa9163bc891acd422f971defe97d64e70d0a4fe5" + +[[package]] +name = "windows-sys" +version = "0.61.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ae137229bcbd6cdf0f7b80a31df61766145077ddf49416a728b02cb3921ff3fc" +dependencies = [ + "windows-link", +] + +[[package]] +name = "zerocopy" +version = "0.8.52" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ce1022995ff5ff5d841ad7d994facc23098cd40152f2c1d11cd607c6f530653f" +dependencies = [ + "zerocopy-derive", +] + +[[package]] +name = "zerocopy-derive" +version = "0.8.52" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1ae7f38b72ec2a254e2b87ef277cf2cd4fb97cbebf944faa6f33354da0867930" +dependencies = [ + "proc-macro2", + "quote", + "syn 2.0.118", +] + +[[package]] +name = "zmij" +version = "1.0.21" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b8848ee67ecc8aedbaf3e4122217aff892639231befc6a1b58d29fff4c2cabaa" diff --git a/Cargo.toml b/Cargo.toml index cb3c758..793984f 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -95,7 +95,14 @@ dbg_macro = "deny" panic = "deny" todo = "deny" exit = "deny" +allow_attributes = "deny" +unreachable = "deny" +let_underscore_must_use = "deny" +let_underscore_untyped = "deny" +unused_result_ok = "deny" # Allow specific pedantic lints that are too noisy module_name_repetitions = "allow" similar_names = "allow" +redundant_pub_crate = "allow" +multiple_crate_versions = "allow" diff --git a/clippy.toml b/clippy.toml new file mode 100644 index 0000000..992e80c --- /dev/null +++ b/clippy.toml @@ -0,0 +1,11 @@ +# Rendered for check.sh, which fails unless this file reads exactly so, and +# which bans every other clippy config in the repo. +# Banned calls stay banned in test code. Tests may print, as the harness +# captures what they print. +allow-unwrap-in-tests = false +allow-expect-in-tests = false +allow-panic-in-tests = false +allow-indexing-slicing-in-tests = false +allow-print-in-tests = true +# Lint public items too, even where the fix changes their API. +avoid-breaking-exported-api = false diff --git a/rust-toolchain.toml b/rust-toolchain.toml new file mode 100644 index 0000000..47d7ef4 --- /dev/null +++ b/rust-toolchain.toml @@ -0,0 +1,5 @@ +# Pinned for check.sh, which stops unless this file reads exactly so. +[toolchain] +channel = "1.99.1" +components = ["clippy", "rustfmt"] +profile = "minimal" From 6a6c081386f0d29ce86e0125609a35dff2b4e1c3 Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 17:40:56 +0000 Subject: [PATCH 2/9] Satisfy the lints without suppressions and cut comments to the density bound; long crate and module docs move to docs/ and are included, so docs.rs shows the same text --- README.md | 4 +- benches/benchmarks.rs | 2 - docs/chebyshev.md | 22 +++++ docs/crate.md | 51 ++++++++++ docs/exp.md | 40 ++++++++ docs/exponential.md | 32 ++++++ docs/hyperbolic.md | 17 ++++ docs/kernel.md | 36 +++++++ docs/sin_cos.md | 20 ++++ docs/sinh_cosh.md | 29 ++++++ docs/sqrt.md | 31 ++++++ src/bin/verify_no_panic.rs | 73 +++++++------- src/bounded.rs | 163 +++++++++--------------------- src/error.rs | 7 +- src/kernel/cordic.rs | 119 ++++------------------ src/kernel/mod.rs | 24 +---- src/lib.rs | 56 +---------- src/ops/algebraic.rs | 25 +++-- src/ops/circular.rs | 154 ++++++++++++----------------- src/ops/exponential.rs | 176 +++++++++++---------------------- src/ops/hyperbolic.rs | 193 +++++++++++++----------------------- src/ops/mod.rs | 13 +-- src/tables/chebyshev.rs | 92 ++++++----------- src/tables/circular.rs | 51 +++++----- src/tables/hyperbolic.rs | 60 ++++++----- src/tables/mod.rs | 13 +-- src/traits.rs | 197 +++++++++---------------------------- 27 files changed, 718 insertions(+), 982 deletions(-) create mode 100644 docs/chebyshev.md create mode 100644 docs/crate.md create mode 100644 docs/exp.md create mode 100644 docs/exponential.md create mode 100644 docs/hyperbolic.md create mode 100644 docs/kernel.md create mode 100644 docs/sin_cos.md create mode 100644 docs/sinh_cosh.md create mode 100644 docs/sqrt.md diff --git a/README.md b/README.md index da01a82..d71f1f5 100644 --- a/README.md +++ b/README.md @@ -63,9 +63,9 @@ The following total functions saturate, clamping to the representable range near | Function | I16F16 Threshold | I32F32 Threshold | Result | |----------|------------------|------------------|--------| | `exp` | x ≥ 10.4 | x ≥ 21.5 | `T::MAX` | -| `exp` | x ≤ -11.1 | x ≤ -22.2 | Zero | +| `exp` | x ≤ -11.8 | x ≤ -22.9 | Zero | | `pow2` | x ≥ 15.0 | x ≥ 31.0 | `T::MAX` | -| `pow2` | x ≤ -16.1 | x ≤ -32.1 | Zero | +| `pow2` | x ≤ -17.0 | x ≤ -33.0 | Zero | | `sinh` | \|x\| ≥ 11.1 | \|x\| ≥ 22.2 | `T::MAX` or `T::MIN` | | `cosh` | \|x\| ≥ 11.1 | \|x\| ≥ 22.2 | `T::MAX` | | `tan` | \|x - pole\| < 4e-5 | \|x - pole\| < 5e-10 | `T::MAX` or `T::MIN` | diff --git a/benches/benchmarks.rs b/benches/benchmarks.rs index 399fa3a..587c413 100644 --- a/benches/benchmarks.rs +++ b/benches/benchmarks.rs @@ -1,7 +1,5 @@ //! Benchmarks for CORDIC functions, on `I16F16` and on `I64F64`. -#![allow(missing_docs, reason = "benchmark code does not need documentation")] - use std::hint::black_box; use criterion::{Criterion, criterion_group, criterion_main}; diff --git a/docs/chebyshev.md b/docs/chebyshev.md new file mode 100644 index 0000000..0607706 --- /dev/null +++ b/docs/chebyshev.md @@ -0,0 +1,22 @@ +Minimax (Chebyshev) polynomial coefficients for sin/cos evaluation. + +The coefficients come from the Remez algorithm (`mpmath.chebyfit`). Each +table approximates a reduced remainder, factored so every coefficient fits in +I1F63 (magnitude below 1). + +Tables are stored highest degree first, `[cₙ, cₙ₋₁, …, c₀]`, and evaluated by +Horner's method: + +```text +P(x) = cₙ·xⁿ + ··· + c₁·x + c₀ = c₀ + x·(c₁ + x·(c₂ + ··· + x·cₙ)) +``` + +# Reconstruction + +| Function | Factored form | Variable | Domain | +|----------|---------------|----------|--------| +| sin | `x + x³·P(x²)` | u = x² | \[0, (π/4)²\] | +| cos | `1 + x²·Q(x²)` | u = x² | \[0, (π/4)²\] | + +The `_LO` tables suit I16F16-class precision (fewer than 24 fractional bits) +and the `_HI` tables I32F32-class and wider. diff --git a/docs/crate.md b/docs/crate.md new file mode 100644 index 0000000..3be52f7 --- /dev/null +++ b/docs/crate.md @@ -0,0 +1,51 @@ +Fixed-point mathematical functions: accurate, deterministic, and guaranteed not to panic. + +`no_std` compatible with no floating-point operations. + +# Quick Start + +```rust +use fixed::types::I16F16; +use fixed_analytics::{sin, cos, sqrt, ln}; + +let angle = I16F16::from_num(0.5); +let (s, c) = (sin(angle), cos(angle)); + +let root = sqrt(I16F16::from_num(2.0)).unwrap(); +assert!((root.to_num::() - 1.414).abs() < 0.001); + +let log = ln(I16F16::E).unwrap(); +assert!((log.to_num::() - 1.0).abs() < 0.01); +``` + +# Available Functions + +**Total functions** return `T` directly, saturating on overflow. +**Fallible functions** return [`Result`] on domain violations. + +| Category | Total | Fallible | +|--------------|-------|----------| +| Trigonometric | [`sin`], [`cos`], [`tan`], [`sin_cos`], [`atan`], [`atan2`] | [`asin`], [`acos`] | +| Hyperbolic | [`sinh`], [`cosh`], [`tanh`], [`sinh_cosh`], [`asinh`] | [`acosh`], [`atanh`], [`acoth`], [`coth`] | +| Exponential | [`exp`], [`pow2`] | [`ln`], [`log2`], [`log10`], [`pow`] | +| Algebraic | — | [`sqrt`] | + +Functions use polynomial evaluation, CORDIC, and Newton-Raphson techniques. +Complete absence of panic is verified at the linker level via the +[`no-panic`](https://github.com/dtolnay/no-panic) crate. + +# Accuracy + +| Type | Typical Accuracy | +|------|------------------| +| `I16F16` | ~4 decimal digits | +| `I32F32` | ~8 decimal digits | + +All functions are benchmarked against MPFR reference implementations. +Accuracy regressions are not permitted across releases. + +# Features + +- **`std`** (default): Enables `std::error::Error` impl on [`Error`] + +See the [`kernel`] module for algorithm details. diff --git a/docs/exp.md b/docs/exp.md new file mode 100644 index 0000000..c77c862 --- /dev/null +++ b/docs/exp.md @@ -0,0 +1,40 @@ +# How `exp` is computed + +This page records the derivations behind the comments in +`src/ops/exponential.rs`. + +## Argument reduction + +`exp(x) = 2^k · exp(r)` with `r ∈ [0, ln 2)`, so `exp(r) ∈ [1, 2)` carries a +full significand into the final shift. + +`k = ⌊x / ln 2⌋` is estimated to within one and then corrected once in either +direction. `to_i32` floors the quotient. When the integer bits do not +outnumber the fractional bits, the quotient comes from a multiply by +`1/ln 2 − 1` (as I1F63, since `1/ln 2` itself exceeds 1) plus `x`, which is +within `|x|·1.5·2^-frac + 2^-frac < 1` of the true value. Other types divide. + +`2^k · exp(r)` exceeds `MAX` once `k ≥ int_bits − 1` and rounds to zero once +`k ≤ −(frac_bits + 2)`. The early exits allow one of slack, and repeat the test +exactly after the correction. + +## Polynomial + +`exp(r) = 1 + r·(1 + r/2·(1 + r/3·(1 + … r/n)))`, evaluated innermost first +with integer division for each `r/K`. The omitted term `(ln 2)^(n+1)/(n+1)!` +is: + +| Degree | Omitted term | Below half an ulp up to | +|--------|--------------|-------------------------| +| 7 | 1.3e-6 | 18 fractional bits | +| 12 | 1.4e-12 | 38 fractional bits | +| 19 | 2.7e-22 | 70 fractional bits | + +The code uses degree 7 below 24 fractional bits, 12 from 24 and 19 from 40. + +## Final scaling + +A positive `k` shifts left after checking `exp(r)` against `MAX >> k`, since +the shift would wrap. A non-positive `k` shifts right with rounding to +nearest: truncation would bias results that are only a few ulps low by up to +a full ulp. diff --git a/docs/exponential.md b/docs/exponential.md new file mode 100644 index 0000000..f2e0e53 --- /dev/null +++ b/docs/exponential.md @@ -0,0 +1,32 @@ +Exponential and logarithmic functions. + +# Saturation + +[`exp`](self::exp), [`pow2`](self::pow2) and [`pow`](self::pow) saturate for +extreme inputs rather than returning an error. The result is `T::MAX` above +`ln(T::MAX)`, or `log2(T::MAX)` for `pow2`, and zero where the true value +rounds to nearest below half an ulp: + +| Function | Result | I16F16 | I32F32 | +|----------|--------|--------|--------| +| `exp` | `T::MAX` | x > ~10.4 | x > ~21.5 | +| `exp` | zero | x < ~-11.8 | x < ~-22.9 | +| `pow2` | `T::MAX` | x > ~15 | x > ~31 | +| `pow2` | zero | x < ~-17 | x < ~-33 | + +Saturation is silent and deterministic. To detect it, check the input range +before the call: + +``` +use fixed::types::I16F16; +use fixed_analytics::exp; + +let x = I16F16::from_num(5.0); +let max_safe = I16F16::from_num(10.0); + +if x < max_safe { + let result = exp(x); +} else { + // Handle potential saturation +} +``` diff --git a/docs/hyperbolic.md b/docs/hyperbolic.md new file mode 100644 index 0000000..e159845 --- /dev/null +++ b/docs/hyperbolic.md @@ -0,0 +1,17 @@ +Hyperbolic functions and their inverses. + +# Saturation + +[`sinh`](self::sinh) and [`cosh`](self::cosh) grow as `e^|x|/2` and saturate +for extreme inputs: `sinh` to `T::MAX` or `T::MIN` past `±asinh(T::MAX)`, and +`cosh`, which is always positive, to `T::MAX` past `acosh(T::MAX)`. +[`sinh_cosh`](self::sinh_cosh) saturates both together, so +`cosh²(x) − sinh²(x) = 1` does not hold for saturated outputs. + +| Type | Saturates for | +|------|---------------| +| I16F16 | \|x\| > ~11.1 | +| I32F32 | \|x\| > ~22.2 | + +Below that, both come from a Taylor polynomial on `|x| ≤ 1.118`, extended by +halving the argument and doubling back. diff --git a/docs/kernel.md b/docs/kernel.md new file mode 100644 index 0000000..6d23550 --- /dev/null +++ b/docs/kernel.md @@ -0,0 +1,36 @@ +CORDIC (Coordinate Rotation Digital Computer) kernels. + +Users should call the functions in [`crate::ops`] rather than these kernels. + +# Algorithm + +CORDIC rotates a vector through a fixed sequence of angles using only shifts +and adds. Each iteration performs one micro-rotation: + +```text +x' = x - σ·d·y·2^(-i) +y' = y + σ·x·2^(-i) +z' = z - σ·angle[i] +``` + +where + +- σ = ±1 is the direction of rotation, +- d = +1 for circular and −1 for hyperbolic mode, +- angle\[i\] = atan(2^-i) for circular and atanh(2^-i) for hyperbolic mode. + +**Vectoring mode** drives y toward zero while z accumulates the angle: + +| Mode | Vectoring (y → 0) | +|------|-------------------| +| Circular | atan | +| Hyperbolic | atanh, ln | + +Hyperbolic mode uses the `atanh(2^-i)` table and repeats iterations 4, 13, +40, … so that the sequence converges. + +# Termination + +Both kernels stop after about a third of the fractional bits and close the +residual angle with one division: for `|w| = |y/x| < 2^-⌈frac/3⌉`, `atan(w)` +and `atanh(w)` differ from `w` by under `|w|³/3 < 2^-frac/3`. diff --git a/docs/sin_cos.md b/docs/sin_cos.md new file mode 100644 index 0000000..fc7c9c9 --- /dev/null +++ b/docs/sin_cos.md @@ -0,0 +1,20 @@ +# How `sin_cos` reduces its argument + +This page records the derivation behind the comments in +`src/ops/circular.rs`. + +An angle outside `[-π, π]` is reduced by `n = round(angle / 2π)` to +`angle − n·2π`, at any magnitude and without iteration limits. + +When the integer bits do not outnumber the fractional bits, the quotient +comes from a multiply by `1/(2π)` as I1F63. That is within +`|angle|·1.5·2^-frac + 2^-frac < 1` of the true quotient, so `n` is off by at +most one, and the clamp to `[-π, π]` that follows absorbs it bit for bit. +Other types divide by 2π. + +`n·2π` is exact and the difference is small, so wrapping arithmetic recovers +it even when `n·2π` exceeds `MAX`, for angles within π of `MAX`. + +The reduced angle then folds to `[-π/2, π/2]`, negating both results. Past +`π/4` the polynomials take `π/2 − |x|` and swap roles, sin from the cosine +polynomial and cos from the sine one, which avoids cancellation near π/2. diff --git a/docs/sinh_cosh.md b/docs/sinh_cosh.md new file mode 100644 index 0000000..382b8fc --- /dev/null +++ b/docs/sinh_cosh.md @@ -0,0 +1,29 @@ +# How `sinh_cosh` is computed + +This page records the derivation behind the comments in +`src/ops/hyperbolic.rs`. + +The argument is halved until `|x| ≤ 1.1182`, the sum of the atanh table, and +both functions come from factored Taylor series in `u = x²`: + +```text +sinh(x) = x·(1 + u/6·(1 + u/20·(1 + u/42·…))) +cosh(x) = 1 + u/2·(1 + u/12·(1 + u/30·…)) +``` + +The divisors are `(2k)(2k+1)` for sinh and `(2k−1)(2k)` for cosh. Integer +division computes each `u/K` in one step without pre-rounding a reciprocal, +and the factored form keeps each multiplicand `u/K` below 1 for `K ≥ 2`. + +At `|x| ≤ 1.1182` the omitted relative term is: + +| Degrees (sinh / cosh) | Omitted term | Below half an ulp up to | +|-----------------------|--------------|-------------------------| +| 9 / 10 | 7.7e-8 / 8.0e-9 | 22 fractional bits | +| 13 / 14 | 3.7e-12 / 2.9e-13 | 37 fractional bits | +| 21 / 22 | 4.5e-22 / 2.4e-23 | 70 fractional bits | + +The code uses the first row below 24 fractional bits, the second from 24 and +the third from 40. The results are rebuilt by doubling, +`sinh(2x) = 2·sinh(x)·cosh(x)` and `cosh(2x) = cosh²(x) + sinh²(x)`, once +per halving. diff --git a/docs/sqrt.md b/docs/sqrt.md new file mode 100644 index 0000000..3c9b6b4 --- /dev/null +++ b/docs/sqrt.md @@ -0,0 +1,31 @@ +# Correctly rounded square root for 128-bit types + +`CordicNumber::sqrt_round` returns the raw bits `⌊√N + ½⌋` with `N = X·2^f`, +where `X` is the raw value and `f` the number of fractional bits. + +## Types up to 64 bits + +When `X << (f + 2)` fits in a wider machine integer, the root comes from one +integer square root: `⌊√N + ½⌋ = (⌊2√N⌋ + 1) >> 1`, and `⌊2√N⌋ = isqrt(4N)`. +The root has at most `(total_bits + f) / 2 + 1` bits, so it fits the raw type. + +## 128-bit types + +`N` can need up to 254 bits, so `sqrt_round_u128` in `src/traits.rs` works in +two steps. + +1. Shift `x` left by the largest `s ≤ leading_zeros(x)` with `s ≡ f (mod 2)`, + giving `q = isqrt(x·2^s)` with 64 significant bits. +2. If `s ≥ f`, `q >> ((s − f)/2)` is exactly `⌊√N⌋`, and the rounding bit + comes from one more bit of the root, or from the remainder `N − q²` when + `s = f`. +3. Otherwise `seed = (q + 1) << m`, with `m = (f − s)/2`, lies above `√N` by + less than `2^(m+1)`. One Newton step from above, + `(seed + ⌊N / seed⌋) / 2`, lands on `⌊√N⌋` or `⌊√N⌋ + 1`, since its error + is below `2^(m − 64)` and `m ≤ 63`. +4. A 256-bit product `newton²` against `N` settles which, and the remainder + `rem = N − root²`, in `[0, 2·root]`, rounds up exactly when `rem > root`. + +The division `⌊x·2^f / seed⌋` is the fixed-point division of the raw values. +It is called only with `seed ≥ √N`, so the quotient is below `2^127` and +cannot overflow; the saturating fallback only keeps the function total. diff --git a/src/bin/verify_no_panic.rs b/src/bin/verify_no_panic.rs index 47d38c9..9c2daaa 100644 --- a/src/bin/verify_no_panic.rs +++ b/src/bin/verify_no_panic.rs @@ -1,10 +1,6 @@ -//! Binary that instantiates every public function with concrete types. -//! -//! This exists solely to trigger monomorphization so that `no_panic`'s -//! linker-level check can verify that no panic paths survive optimization. -//! It is only compiled under the `verify-no-panic` feature, for three -//! layouts that compile to different code: `I16F16`, `I24F8` (more integer -//! than fractional bits) and `I64F64` (128-bit). +//! Instantiates every public function so `no_panic`'s link check can prove +//! it panic-free, for three layouts that compile differently: `I16F16`, +//! `I24F8` (more integer than fractional bits) and `I64F64` (128-bit). #[cfg(not(feature = "verify-no-panic"))] compile_error!("this binary should only be built with --features verify-no-panic"); @@ -19,47 +15,48 @@ use fixed_analytics::{ pow, pow2, sin, sin_cos, sinh, sinh_cosh, sqrt, tan, tanh, }; +/// Hands `v` to `black_box`, so the optimizer must compute it in full. +fn keep(v: V) { + std::hint::black_box(v); +} + fn exercise(x: T, y: T, two: T) { - // Total functions (return T) - let _ = std::hint::black_box(sin(x)); - let _ = std::hint::black_box(cos(x)); - let _ = std::hint::black_box(tan(x)); - let _ = std::hint::black_box(sin_cos(x)); - let _ = std::hint::black_box(atan(x)); - let _ = std::hint::black_box(atan2(y, x)); - let _ = std::hint::black_box(exp(x)); - let _ = std::hint::black_box(pow2(x)); - let _ = std::hint::black_box(sinh(x)); - let _ = std::hint::black_box(cosh(x)); - let _ = std::hint::black_box(tanh(x)); - let _ = std::hint::black_box(sinh_cosh(x)); - let _ = std::hint::black_box(asinh(x)); - let _ = std::hint::black_box(asinh(two)); + keep(sin(x)); + keep(cos(x)); + keep(tan(x)); + keep(sin_cos(x)); + keep(atan(x)); + keep(atan2(y, x)); + keep(exp(x)); + keep(pow2(x)); + keep(sinh(x)); + keep(cosh(x)); + keep(tanh(x)); + keep(sinh_cosh(x)); + keep(asinh(x)); + keep(asinh(two)); - // Fallible functions (return Result) - let _ = std::hint::black_box(asin(x)); - let _ = std::hint::black_box(acos(x)); - let _ = std::hint::black_box(sqrt(x)); - let _ = std::hint::black_box(ln(x)); - let _ = std::hint::black_box(log2(x)); - let _ = std::hint::black_box(log10(x)); - let _ = std::hint::black_box(acosh(two)); - let _ = std::hint::black_box(atanh(x)); - let _ = std::hint::black_box(coth(x)); - let _ = std::hint::black_box(acoth(two)); - let _ = std::hint::black_box(pow(two, x)); + keep(asin(x)); + keep(acos(x)); + keep(sqrt(x)); + keep(ln(x)); + keep(log2(x)); + keep(log10(x)); + keep(acosh(two)); + keep(atanh(x)); + keep(coth(x)); + keep(acoth(two)); + keep(pow(two, x)); - // Type-safe wrapper functions if let Some(nn) = NonNegative::new(x) { - let _ = std::hint::black_box(sqrt_nonneg(nn)); + keep(sqrt_nonneg(nn)); } if let Some(ou) = OpenUnitInterval::new(x) { - let _ = std::hint::black_box(atanh_open(ou)); + keep(atanh_open(ou)); } } fn main() { - // Use black_box to prevent the optimizer from eliminating calls entirely. exercise( std::hint::black_box(I16F16::from_num(0.5)), std::hint::black_box(I16F16::from_num(0.25)), diff --git a/src/bounded.rs b/src/bounded.rs index 7da2ecd..17cfa50 100644 --- a/src/bounded.rs +++ b/src/bounded.rs @@ -1,53 +1,24 @@ -//! Bounded value types that encode mathematical invariants at the type level. +//! Value types whose range is known from how they were built. //! -//! These types provide compile-time guarantees about value ranges, eliminating -//! the need for runtime checks in internal computations where the domain is -//! already validated. -//! -//! # Types -//! -//! - [`NonNegative`]: Values >= 0 (for sqrt inputs) -//! - [`UnitInterval`]: Values in [-1, 1] (for asin/acos inputs) -//! - [`OpenUnitInterval`]: Values in (-1, 1) (for atanh inputs) -//! -//! # Design Philosophy -//! -//! Rather than using `unsafe` or `expect`, these types encode the mathematical -//! relationships between operations. For example: -//! -//! - `1 + x^2` is always >= 1, so `NonNegative::one_plus_square(x)` is infallible -//! - If `|x| <= 1`, then `1 - x^2` is in [0, 1], so `NonNegative::one_minus_square(unit_x)` is infallible -//! - `x / sqrt(1 + x^2)` is always in (-1, 1), so `OpenUnitInterval::x_div_sqrt_one_plus_x_sq(x)` is infallible +//! Each holds a value checked by `new`, or the result of an identity that +//! cannot leave the range, such as [`NonNegative::one_plus_square`] for +//! `1 + x²`. Internal code then needs no domain checks, `unsafe` or `expect`. use crate::traits::CordicNumber; -/// A value guaranteed to be non-negative (>= 0). -/// -/// This type enables infallible sqrt operations by encoding the non-negativity -/// constraint at the type level. -/// -/// # Construction -/// -/// - [`NonNegative::new`]: Checked construction, returns `Option` -/// - [`NonNegative::one_plus_square`]: From `1 + x^2`, always valid -/// - [`NonNegative::one_minus_square`]: From `1 - x^2` where `|x| <= 1`, always valid -/// - [`NonNegative::square_minus_one`]: From `x^2 - 1` where `|x| >= 1`, always valid +/// A value `≥ 0`, the input to an infallible square root. #[derive(Clone, Copy, Debug)] pub struct NonNegative(T); impl NonNegative { - /// Creates a new `NonNegative` value if the input is >= 0. - /// - /// Returns `None` if the value is negative. + /// Returns `Some` if `value ≥ 0`. #[inline] #[must_use] pub fn new(value: T) -> Option { (value >= T::zero()).then_some(Self(value)) } - /// Constructs from `1 + x^2`, which is always >= 1. - /// - /// This is mathematically infallible: for any real `x`, `1 + x^2 >= 1`. + /// Builds `1 + x²`, which is at least 1. #[inline] #[must_use] pub fn one_plus_square(x: T) -> Self { @@ -55,10 +26,7 @@ impl NonNegative { Self(T::one().saturating_add(x_sq)) } - /// Constructs from `1 - x^2` where `|x| <= 1`. - /// - /// This is mathematically infallible: if `|x| <= 1`, then `x^2 <= 1`, - /// so `1 - x^2 >= 0`. + /// Builds `1 − x²`, which is non-negative for `|x| ≤ 1`. #[inline] #[must_use] pub fn one_minus_square(x: UnitInterval) -> Self { @@ -66,10 +34,7 @@ impl NonNegative { Self(T::one().saturating_sub(x_sq)) } - /// Constructs from `x^2 - 1` where `|x| >= 1`. - /// - /// This is mathematically infallible: if `|x| >= 1`, then `x^2 >= 1`, - /// so `x^2 - 1 >= 0`. + /// Builds `x² − 1`, which is non-negative for `x ≥ 1`. #[inline] #[must_use] pub fn square_minus_one(x: AtLeastOne) -> Self { @@ -77,7 +42,7 @@ impl NonNegative { Self(x_sq.saturating_sub(T::one())) } - /// Unwraps the inner value. + /// Returns the inner value. #[inline] #[must_use] pub const fn get(self) -> T { @@ -85,17 +50,12 @@ impl NonNegative { } } -/// A value guaranteed to be in the closed interval [-1, 1]. -/// -/// This type is used for inputs to functions like asin and acos that -/// require their argument to be in this range. +/// A value in the closed interval `[-1, 1]`, as asin and acos take. #[derive(Clone, Copy, Debug)] pub struct UnitInterval(T); impl UnitInterval { - /// Creates a new `UnitInterval` value if the input is in [-1, 1]. - /// - /// Returns `None` if the value is outside the interval. + /// Returns `Some` if `value` is in `[-1, 1]`. #[inline] #[must_use] pub fn new(value: T) -> Option { @@ -103,14 +63,14 @@ impl UnitInterval { (value >= -one && value <= one).then_some(Self(value)) } - /// Constructs `1 / x` where `x >= 1`, which is always in (0, 1]. + /// Builds `1 / x`, which is in `(0, 1]` for `x ≥ 1`. #[inline] #[must_use] pub fn from_reciprocal(x: AtLeastOne) -> Self { Self(T::one().div(x.0)) } - /// Unwraps the inner value. + /// Returns the inner value. #[inline] #[must_use] pub const fn get(self) -> T { @@ -118,16 +78,12 @@ impl UnitInterval { } } -/// A value guaranteed to be in the open interval (-1, 1). -/// -/// This type is used for inputs to atanh, which requires strict inequality. +/// A value in the open interval `(-1, 1)`, as atanh takes. #[derive(Clone, Copy, Debug)] pub struct OpenUnitInterval(T); impl OpenUnitInterval { - /// Creates a new `OpenUnitInterval` value if the input is in (-1, 1). - /// - /// Returns `None` if the value is outside the interval or on the boundary. + /// Returns `Some` if `value` is in `(-1, 1)`. #[inline] #[must_use] pub fn new(value: T) -> Option { @@ -135,34 +91,23 @@ impl OpenUnitInterval { (value > -one && value < one).then_some(Self(value)) } - /// Constructs from `x / sqrt(1 + x^2)`, which is always in (-1, 1). - /// - /// This is mathematically infallible: for any real `x`, - /// `|x / sqrt(1 + x^2)| < 1` because `sqrt(1 + x^2) > |x|`. - /// - /// Note: Requires the sqrt to be computed first. + /// Builds `x / √(1 + x²)` from the computed root; it is below 1 in magnitude as + /// `√(1 + x²) > |x|`. #[inline] #[must_use] pub fn from_div_by_sqrt_one_plus_square(x: T, sqrt_one_plus_x_sq: T) -> Self { Self(x.div(sqrt_one_plus_x_sq)) } - /// Constructs from `sqrt(x^2 - 1) / x` where `|x| >= 1`. - /// - /// This is mathematically infallible: for `|x| >= 1`, - /// `|sqrt(x^2 - 1) / x| < 1` because `sqrt(x^2 - 1) < |x|` for `|x| > 1`. - /// - /// Note: Requires the sqrt to be computed first. + /// Builds `√(x² − 1) / x` from the computed root; it is below 1 in magnitude as + /// `√(x² − 1) < x`. #[inline] #[must_use] pub fn from_sqrt_square_minus_one_div(sqrt_x_sq_minus_one: T, x: AtLeastOne) -> Self { Self(sqrt_x_sq_minus_one.div(x.0)) } - /// Constructs from `(x - 1) / (x + 1)` where `x` is in [0.5, 2]. - /// - /// This is mathematically infallible: for `x` in [0.5, 2], - /// `(x - 1) / (x + 1)` is in `(-1/3, 1/3)` which is a subset of `(-1, 1)`. + /// Builds `(x − 1) / (x + 1)`, which lies in `[-1/3, 1/3]` for `x` in `[0.5, 2]`. #[inline] #[must_use] pub(crate) fn from_normalized_ln_arg(x: NormalizedLnArg) -> Self { @@ -171,7 +116,7 @@ impl OpenUnitInterval { Self(x_minus_1.div(x_plus_1)) } - /// Unwraps the inner value. + /// Returns the inner value. #[inline] #[must_use] pub const fn get(self) -> T { @@ -179,23 +124,19 @@ impl OpenUnitInterval { } } -/// A value guaranteed to be >= 1 (or <= -1 for the absolute value). -/// -/// This type is used for inputs to acosh which requires x >= 1. +/// A value `≥ 1`, as acosh takes. #[derive(Clone, Copy, Debug)] pub struct AtLeastOne(T); impl AtLeastOne { - /// Creates a new `AtLeastOne` value if the input is >= 1. - /// - /// Returns `None` if the value is less than 1. + /// Returns `Some` if `value ≥ 1`. #[inline] #[must_use] pub fn new(value: T) -> Option { (value >= T::one()).then_some(Self(value)) } - /// Unwraps the inner value. + /// Returns the inner value. #[inline] #[must_use] pub const fn get(self) -> T { @@ -203,18 +144,12 @@ impl AtLeastOne { } } -/// A value guaranteed to be in [0.5, 2], used for ln argument normalization. -/// -/// After normalizing the input for ln computation, the value is always -/// in this range, which guarantees that `(x-1)/(x+1)` is in `(-1/3, 1/3)`. +/// The ln argument after normalization, in `[0.5, 2]`. #[derive(Clone, Copy, Debug)] pub(crate) struct NormalizedLnArg(T); impl NormalizedLnArg { - /// Creates a new `NormalizedLnArg` from the normalization loop result. - /// - /// The ln function's normalization loop guarantees the result is in [0.5, 2]. - /// This constructor trusts that invariant. + /// Wraps the normalized value unchecked; `ln_positive` guarantees the range. #[inline] #[must_use] pub(crate) const fn from_normalized(value: T) -> Self { @@ -223,7 +158,6 @@ impl NormalizedLnArg { } #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod tests { use super::*; use fixed::types::I16F16; @@ -243,18 +177,20 @@ mod tests { } #[test] - fn non_negative_one_minus_square() { - let unit = UnitInterval::new(I16F16::from_num(0.5)).unwrap(); + fn non_negative_one_minus_square() -> Result<(), &'static str> { + let unit = UnitInterval::new(I16F16::from_num(0.5)).ok_or("out of range")?; let nn = NonNegative::one_minus_square(unit); let val: f32 = nn.get().to_num(); assert!((val - 0.75).abs() < 0.01); + Ok(()) } #[test] - fn non_negative_square_minus_one() { - let at_least = AtLeastOne::new(I16F16::from_num(2)).unwrap(); + fn non_negative_square_minus_one() -> Result<(), &'static str> { + let at_least = AtLeastOne::new(I16F16::from_num(2)).ok_or("out of range")?; let nn = NonNegative::square_minus_one(at_least); assert_eq!(nn.get(), I16F16::from_num(3)); + Ok(()) } #[test] @@ -267,18 +203,20 @@ mod tests { } #[test] - fn unit_interval_from_reciprocal() { - let at_least = AtLeastOne::new(I16F16::from_num(4)).unwrap(); + fn unit_interval_from_reciprocal() -> Result<(), &'static str> { + let at_least = AtLeastOne::new(I16F16::from_num(4)).ok_or("out of range")?; assert_eq!( UnitInterval::from_reciprocal(at_least).get(), I16F16::from_num(0.25) ); + Ok(()) } #[test] - fn unit_interval_get() { - let unit = UnitInterval::new(I16F16::from_num(0.5)).unwrap(); + fn unit_interval_get() -> Result<(), &'static str> { + let unit = UnitInterval::new(I16F16::from_num(0.5)).ok_or("out of range")?; assert_eq!(unit.get(), I16F16::from_num(0.5)); + Ok(()) } #[test] @@ -290,9 +228,10 @@ mod tests { } #[test] - fn open_unit_interval_get() { - let open = OpenUnitInterval::new(I16F16::from_num(0.5)).unwrap(); + fn open_unit_interval_get() -> Result<(), &'static str> { + let open = OpenUnitInterval::new(I16F16::from_num(0.5)).ok_or("out of range")?; assert_eq!(open.get(), I16F16::from_num(0.5)); + Ok(()) } #[test] @@ -305,21 +244,14 @@ mod tests { } #[test] - fn open_unit_interval_from_sqrt_div() { - let at_least = AtLeastOne::new(I16F16::from_num(2)).unwrap(); - #[allow( - clippy::approx_constant, - reason = "testing with known approximation of sqrt(3)" - )] + fn open_unit_interval_from_sqrt_div() -> Result<(), &'static str> { + let at_least = AtLeastOne::new(I16F16::from_num(2)).ok_or("out of range")?; let sqrt_x_sq_minus_one = I16F16::from_num(1.732); let open = OpenUnitInterval::from_sqrt_square_minus_one_div(sqrt_x_sq_minus_one, at_least); let val: f32 = open.get().to_num(); - #[allow( - clippy::approx_constant, - reason = "testing with known approximation of sqrt(3)/2" - )] let expected = 0.866_f32; assert!((val - expected).abs() < 0.01); + Ok(()) } #[test] @@ -330,9 +262,10 @@ mod tests { } #[test] - fn at_least_one_get() { - let at_least = AtLeastOne::new(I16F16::from_num(2)).unwrap(); + fn at_least_one_get() -> Result<(), &'static str> { + let at_least = AtLeastOne::new(I16F16::from_num(2)).ok_or("out of range")?; assert_eq!(at_least.get(), I16F16::from_num(2)); + Ok(()) } #[test] diff --git a/src/error.rs b/src/error.rs index 93c07f3..78adc0f 100644 --- a/src/error.rs +++ b/src/error.rs @@ -5,10 +5,7 @@ use core::fmt; /// Errors that can occur during CORDIC computations. #[derive(Debug, Clone, Copy, PartialEq, Eq)] pub enum Error { - /// Input value is outside the valid domain for the function. - /// - /// For example, `asin(2.0)` would produce this error since - /// arcsine is only defined for inputs in [-1, 1]. + /// The input is outside the function's domain, as for `asin(2.0)`. DomainError { /// Name of the function that encountered the error. function: &'static str, @@ -18,7 +15,7 @@ pub enum Error { } impl Error { - /// Create a domain error for the given function. + /// A domain error for the given function. #[must_use] pub const fn domain(function: &'static str, expected: &'static str) -> Self { Self::DomainError { function, expected } diff --git a/src/kernel/cordic.rs b/src/kernel/cordic.rs index f193c6e..afeb413 100644 --- a/src/kernel/cordic.rs +++ b/src/kernel/cordic.rs @@ -1,70 +1,34 @@ -//! Core CORDIC iteration implementations. -//! -//! CORDIC vectoring mode drives y toward zero while accumulating angles: -//! -//! | Mode | Vectoring (y → 0) | -//! |------|-------------------| -//! | Circular | atan | -//! | Hyperbolic | atanh, ln | -//! -//! # Algorithm -//! -//! Each iteration performs a micro-rotation: -//! ```text -//! x' = x - σ * d * y * 2^(-i) -//! y' = y + σ * x * 2^(-i) -//! z' = z - σ * angle[i] -//! ``` -//! -//! Where: -//! - σ = ±1 (direction of rotation) -//! - d = +1 for circular, -1 for hyperbolic, 0 for linear -//! - angle[i] = atan(2^-i) for circular, atanh(2^-i) for hyperbolic -//! -//! Both kernels stop after about a third of the fractional bits and close -//! the residual angle with one division: for `|w| = |y/x| < 2^-⌈frac/3⌉`, -//! `atan(w)` and `atanh(w)` differ from `w` by under `|w|³/3 < 2^-frac/3`. +//! CORDIC vectoring kernels; the [kernel module docs](super) cover the algorithm. use crate::tables::hyperbolic::needs_repeat; use crate::tables::{ATAN_TABLE, ATANH_TABLE}; use crate::traits::CordicNumber; -/// Index of the last shift stage before the closing division: after the -/// stage with shift `2^-k` the residual angle, and so `|y/x|`, is below -/// `2^-k`. At most 42, so every table lookup is in bounds. +/// Last shift stage before the closing division, at most 42. After stage `k` +/// the residual `|y/x|` is below `2^-k`. const fn last_stage(frac_bits: u32) -> u32 { frac_bits.div_ceil(3) } -/// Table lookup for CORDIC iteration. The modulo is a no-op that keeps the -/// access provably in bounds for the no-panic check. +/// Table lookup; the modulo keeps the index in bounds, so the fallback is dead. #[inline] -const fn table_lookup(table: &[i64; 64], index: u32) -> i64 { - #[allow(clippy::indexing_slicing, reason = "index reduced modulo table length")] - table[index as usize % table.len()] +fn table_lookup(table: &[i64; 64], index: u32) -> i64 { + match table.get(index as usize % table.len()) { + Some(&v) => v, + None => 0, + } } -/// Converts an I1F63 table constant to `T`, rounding to nearest. -/// -/// `CordicNumber::from_i1f63` truncates, which biases every converted -/// angle low by up to one ulp. For circular CORDIC the accumulated z is -/// dominated by this table quantization, so rounding to nearest halves -/// the worst case and removes the systematic bias. -/// -/// Hyperbolic CORDIC deliberately keeps the truncating conversion: its -/// datapath shifts carry their own bias, which the truncated table -/// angles partially cancel (measured, not designed — see the accuracy -/// baseline history). +/// Converts an I1F63 table angle to `T`, rounding to nearest: truncation would +/// bias circular z low. Hyperbolic mode keeps truncation, which partly cancels +/// its shift bias (measured in the accuracy baseline). fn from_i1f63_rounded(bits: i64) -> T { let truncated = T::from_i1f63(bits); let frac = T::frac_bits(); if frac >= 63 { - // Conversion is exact; nothing to round. return truncated; } - // The conversion discards the low (63 - frac) bits. If the highest - // discarded bit is set, the true value is at least half an output - // ulp above the truncated result, so round up by one ulp. + // Round up when the highest discarded bit is set. let shift = 63 - frac; if (bits >> (shift - 1)) & 1 == 1 { let ulp = T::one() >> frac; @@ -74,27 +38,10 @@ fn from_i1f63_rounded(bits: i64) -> T { } } -/// Performs circular CORDIC in vectoring mode. -/// -/// Given an initial vector (x, y) with x > 0, rotates it toward the x axis -/// through `⌈frac/3⌉ + 1` shift stages, then closes the residual angle -/// with one division. -/// -/// # Arguments -/// -/// * `x` - Initial x coordinate (must be positive) -/// * `y` - Initial y coordinate -/// * `z` - Initial angle accumulator (usually 0) -/// -/// # Returns -/// -/// Tuple of (x, y, z): z ≈ z₀ + atan(y₀/x₀); x ≈ K·sqrt(x₀² + y₀²) with K -/// the gain of the stages run; y is the residual the shift stages left, -/// whose angle is already in z. +/// Circular CORDIC in vectoring mode for `x > 0`; `(x, y, 0)` gives atan(y/x). /// -/// # Note -/// -/// For computing atan(y/x), pass (1, y/x, 0) or (x, y, 0). +/// Returns z ≈ z₀ + atan(y₀/x₀), x ≈ K·√(x₀² + y₀²) with K the stages' gain, +/// and the residual y, whose angle is already in z. #[must_use] pub fn circular_vectoring(mut x: T, mut y: T, mut z: T) -> (T, T, T) { let zero = T::zero(); @@ -103,13 +50,11 @@ pub fn circular_vectoring(mut x: T, mut y: T, mut z: T) -> (T, let angle = from_i1f63_rounded::(table_lookup(&ATAN_TABLE, i)); if y < zero { - // y is negative, rotate counter-clockwise to bring y toward zero let x_new = x.saturating_sub(y >> i); y = y.saturating_add(x >> i); x = x_new; z -= angle; } else { - // y is positive or zero, rotate clockwise let x_new = x.saturating_add(y >> i); y = y.saturating_sub(x >> i); x = x_new; @@ -123,28 +68,10 @@ pub fn circular_vectoring(mut x: T, mut y: T, mut z: T) -> (T, (x, y, z) } -/// Performs hyperbolic CORDIC in vectoring mode. -/// -/// Drives y toward zero through shift stages 1..=⌈frac/3⌉ (repeating 4, -/// 13, 40, … for convergence), then closes the residual angle with one -/// division. -/// -/// # Arguments -/// -/// * `x` - Initial x value (should satisfy |x| > |y|) -/// * `y` - Initial y value -/// * `z` - Initial angle accumulator -/// -/// # Returns -/// -/// Tuple of (x, y, z): z ≈ z₀ + atanh(y₀/x₀); x ≈ `K_h`·sqrt(x₀² - y₀²) -/// with `K_h` the gain of the stages run; y is the residual the shift -/// stages left, whose angle is already in z. -/// -/// # Note +/// Hyperbolic CORDIC in vectoring mode for `|x| > |y|`; `(1, v, 0)` gives atanh(v). /// -/// For computing atanh(v), pass (1, v, 0) where |v| < 1. -/// For computing ln(x), use the identity: ln(x) = 2 * atanh((x-1)/(x+1)) +/// Returns z ≈ z₀ + atanh(y₀/x₀), x ≈ `K_h`·√(x₀² − y₀²) with `K_h` the stages' +/// gain, and the residual y, whose angle is already in z. #[must_use] pub fn hyperbolic_vectoring(mut x: T, mut y: T, mut z: T) -> (T, T, T) { let zero = T::zero(); @@ -155,19 +82,13 @@ pub fn hyperbolic_vectoring(mut x: T, mut y: T, mut z: T) -> (T let passes = if needs_repeat(i) { 2 } else { 1 }; for _ in 0..passes { - // Hyperbolic pseudo-rotation equations: - // x' = x + σ*y*2^(-i) - // y' = y + σ*x*2^(-i) - // z' = z + σ*angle (accumulating for vectoring) - // where σ = -sign(y) to drive y toward zero + // σ = −sign(y) drives y toward zero. if y < zero { - // y is negative: σ = +1 let x_new = x.saturating_add(y >> i); y = y.saturating_add(x >> i); x = x_new; z -= angle; } else { - // y is positive or zero: σ = -1 let x_new = x.saturating_sub(y >> i); y = y.saturating_sub(x >> i); x = x_new; diff --git a/src/kernel/mod.rs b/src/kernel/mod.rs index 9fc237c..697ef48 100644 --- a/src/kernel/mod.rs +++ b/src/kernel/mod.rs @@ -1,26 +1,4 @@ -//! CORDIC (Coordinate Rotation Digital Computer) kernels. -//! -//! # Algorithm -//! -//! Iteratively rotates vectors using only shifts and adds: -//! -//! ```text -//! x' = x - σ·y·2^(-i) -//! y' = y + σ·x·2^(-i) -//! z' = z - σ·atan(2^-i) -//! ``` -//! -//! **Vectoring mode** (y→0): computes atan/atanh from coordinates. -//! -//! Hyperbolic mode uses `atanh(2^-i)` tables and requires iteration repeats -//! at indices 4, 13, 40, ... for convergence. -//! -//! | Mode | Vectoring (y → 0) | -//! |------|-------------------| -//! | Circular | atan | -//! | Hyperbolic | atanh, ln | -//! -//! Users should call functions in [`crate::ops`] rather than kernels directly. +#![doc = include_str!("../../docs/kernel.md")] mod cordic; diff --git a/src/lib.rs b/src/lib.rs index f667fd6..86e9cf4 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -1,55 +1,4 @@ -//! Fixed-point mathematical functions: accurate, deterministic, and guaranteed not to panic. -//! -//! `no_std` compatible with no floating-point operations. -//! -//! # Quick Start -//! -//! ```rust -//! use fixed::types::I16F16; -//! use fixed_analytics::{sin, cos, sqrt, ln}; -//! -//! let angle = I16F16::from_num(0.5); -//! let (s, c) = (sin(angle), cos(angle)); -//! -//! let root = sqrt(I16F16::from_num(2.0)).unwrap(); -//! assert!((root.to_num::() - 1.414).abs() < 0.001); -//! -//! let log = ln(I16F16::E).unwrap(); -//! assert!((log.to_num::() - 1.0).abs() < 0.01); -//! ``` -//! -//! # Available Functions -//! -//! **Total functions** return `T` directly, saturating on overflow. -//! **Fallible functions** return [`Result`] on domain violations. -//! -//! | Category | Total | Fallible | -//! |--------------|-------|----------| -//! | Trigonometric | [`sin`], [`cos`], [`tan`], [`sin_cos`], [`atan`], [`atan2`] | [`asin`], [`acos`] | -//! | Hyperbolic | [`sinh`], [`cosh`], [`tanh`], [`sinh_cosh`], [`asinh`] | [`acosh`], [`atanh`], [`acoth`], [`coth`] | -//! | Exponential | [`exp`], [`pow2`] | [`ln`], [`log2`], [`log10`], [`pow`] | -//! | Algebraic | — | [`sqrt`] | -//! -//! Functions use polynomial evaluation, CORDIC, and Newton-Raphson techniques. -//! Complete absence of panic is verified at the linker level via the -//! [`no-panic`](https://github.com/dtolnay/no-panic) crate. -//! -//! # Accuracy -//! -//! | Type | Typical Accuracy | -//! |------|------------------| -//! | `I16F16` | ~4 decimal digits | -//! | `I32F32` | ~8 decimal digits | -//! -//! All functions are benchmarked against MPFR reference implementations. -//! Accuracy regressions are not permitted across releases. -//! -//! # Features -//! -//! - **`std`** (default): Enables `std::error::Error` impl on [`Error`] -//! -//! See the [`kernel`] module for algorithm details. - +#![doc = include_str!("../docs/crate.md")] #![no_std] #![cfg_attr(docsrs, feature(doc_cfg))] @@ -60,14 +9,11 @@ pub mod ops; pub mod tables; pub mod traits; -// Re-export the fixed crate for convenience pub use fixed; -// Re-export main types pub use error::{Error, Result}; pub use traits::CordicNumber; -// Re-export all mathematical functions at crate root for convenience pub use ops::algebraic::sqrt; pub use ops::circular::{acos, asin, atan, atan2, cos, sin, sin_cos, tan}; pub use ops::exponential::{exp, ln, log2, log10, pow, pow2}; diff --git a/src/ops/algebraic.rs b/src/ops/algebraic.rs index 5c17890..b13e4aa 100644 --- a/src/ops/algebraic.rs +++ b/src/ops/algebraic.rs @@ -5,31 +5,28 @@ use crate::error::{Error, Result}; use crate::traits::CordicNumber; /// Square root, rounded to the nearest representable value. Domain: `x ≥ 0`. -/// /// # Errors /// Returns `DomainError` if `x < 0`. #[must_use = "returns the square root result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sqrt(x: T) -> Result { NonNegative::new(x) .map(sqrt_nonneg) .ok_or_else(|| Error::domain("sqrt", "non-negative value")) } -/// Infallible square root for non-negative values, rounded to the nearest -/// representable value. -/// -/// This function takes a [`NonNegative`] wrapper, guaranteeing at the type -/// level that the input is valid. No domain check is performed at runtime. -/// -/// Use this when the non-negativity of the input is already established -/// through mathematical invariants (e.g., `1 + x²`, `1 - x²` for `|x| ≤ 1`). +/// Square root of a [`NonNegative`] value, with no runtime domain check. /// -/// The result is correctly rounded at every magnitude (raw bits -/// `⌊√(X·2^f) + ½⌋`), computed from an integer square root of the raw -/// representation; see [`CordicNumber::sqrt_round`]. +/// The result is correctly rounded at every magnitude: raw bits +/// `⌊√(X·2^f) + ½⌋`, as [`CordicNumber::sqrt_round`] computes them. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sqrt_nonneg(x: NonNegative) -> T { x.get().sqrt_round() } diff --git a/src/ops/circular.rs b/src/ops/circular.rs index 153939c..679ce79 100644 --- a/src/ops/circular.rs +++ b/src/ops/circular.rs @@ -10,37 +10,34 @@ use crate::traits::CordicNumber; /// `1/(2π)` as I1F63, for the angle reduction multiply in [`sin_cos`]. const FRAC_1_2PI_I1F63: i64 = 0x145F_306D_C9C8_82A5; -/// Sine and cosine. More efficient than separate calls. Accepts any angle. +/// Sine and cosine of any angle, cheaper than separate calls. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sin_cos(angle: T) -> (T, T) { let pi = T::pi(); let frac_pi_2 = T::frac_pi_2(); let two_pi = pi + pi; - // Reduce angle to [-π, π] using direct quotient computation. - // This handles arbitrarily large angles without iteration limits. + // Reduce to [-π, π] by n = round(angle / 2π) at any magnitude. The + // reciprocal multiply leaves n off by at most one, which the clamp below + // absorbs bit for bit (docs/sin_cos.md); other types divide. let reduced = if angle > pi || angle < -pi { - // Compute n = round(angle / 2π), then reduced = angle - n * 2π. - // The reciprocal multiply is within |angle|·1.5·2^-frac + 2^-frac < 1 - // of the true quotient when the integer bits do not outnumber the - // fractional bits, so n is off by at most one and the clamp below - // absorbs it bit for bit; other types keep the division. let quotient = if T::total_bits() <= 2 * T::frac_bits() { angle.saturating_mul(T::from_i1f63(FRAC_1_2PI_I1F63)) } else { angle.div(two_pi) }; let n = quotient.round(); - // n·2π is exact and the difference small, so wrapping arithmetic - // recovers it even when n·2π exceeds MAX (angles within π of MAX). + // n·2π is exact and the difference small, so wrapping recovers it even + // when n·2π exceeds MAX. angle.wrapping_sub(n.wrapping_mul(two_pi)) } else { angle }; - // Clamp to [-π, π] to handle any residual from saturation. - // This is a safety net; mathematically unnecessary for valid inputs. let reduced = if reduced > pi { reduced.saturating_sub(two_pi) } else if reduced < -pi { @@ -49,7 +46,6 @@ pub fn sin_cos(angle: T) -> (T, T) { reduced }; - // Further reduce to [-π/2, π/2] and track sign let (reduced, negate) = if reduced > frac_pi_2 { (reduced - pi, true) } else if reduced < -frac_pi_2 { @@ -58,10 +54,8 @@ pub fn sin_cos(angle: T) -> (T, T) { (reduced, false) }; - // Polynomial evaluation via factored Horner form. - // To avoid catastrophic cancellation near π/2, reduce to [0, π/4]: - // For |x| ∈ [0, π/4]: sin(x) = sin_poly(x), cos(x) = cos_poly(x) - // For |x| ∈ (π/4, π/2]: sin(x) = cos_poly(π/2-|x|), cos(x) = sin_poly(π/2-|x|) + // Past π/4 the polynomials take π/2 − |x| and swap roles, which avoids + // cancellation near π/2. let one = T::one(); let frac_pi_4 = T::frac_pi_4(); let abs_reduced = reduced.abs(); @@ -72,28 +66,20 @@ pub fn sin_cos(angle: T) -> (T, T) { }; let u = poly_arg.saturating_mul(poly_arg); - // Evaluate sin and cos polynomials over [0, π/4] using minimax - // (Chebyshev) coefficients. Uses multiply-by-constant instead of - // division, avoiding cumulative rounding error from per-step divides. - // - // sin(x) = x + x³·P(x²) where P = minimax poly of (sin(x)-x)/x³ - // cos(x) = 1 + x²·Q(x²) where Q = minimax poly of (cos(x)-1)/x² + // Minimax polynomials, multiplying by constants rather than dividing at + // each step: degree 15 and 14 from 24 fractional bits, else 9 and 8. let (sp_val, cp_val) = if T::frac_bits() >= 24 { - // High precision: degree 15 sin, degree 14 cos let sp = horner(&SIN_P_HI, u); let sin_approx = poly_arg.saturating_add(poly_arg.saturating_mul(u).saturating_mul(sp)); let cp = horner(&COS_Q_HI, u); (sin_approx, one.saturating_add(u.saturating_mul(cp))) } else { - // Low precision: degree 9 sin, degree 8 cos let sp = horner(&SIN_P_LO, u); let sin_approx = poly_arg.saturating_add(poly_arg.saturating_mul(u).saturating_mul(sp)); let cp = horner(&COS_Q_LO, u); (sin_approx, one.saturating_add(u.saturating_mul(cp))) }; - // Map back: if we swapped, sin(x) = cos_poly, cos(x) = sin_poly - // Also restore sign of sin for negative angles. let (sin_unsigned, cos_val) = if swapped { (cp_val, sp_val) } else { @@ -112,44 +98,32 @@ pub fn sin_cos(angle: T) -> (T, T) { } } -/// Sine. Accepts any angle (reduced internally). +/// Sine of any angle. #[inline] #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sin(angle: T) -> T { sin_cos(angle).0 } -/// Cosine. Accepts any angle (reduced internally). +/// Cosine of any angle. #[inline] #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn cos(angle: T) -> T { sin_cos(angle).1 } -/// Tangent. Returns `sin(angle) / cos(angle)`. -/// -/// # Overflow Behavior -/// -/// Tangent has poles at ±π/2, ±3π/2, etc. where it approaches ±∞. -/// Since these poles occur at irrational values that cannot be exactly -/// represented in fixed-point, this function will never divide by -/// exactly zero. However, near poles the result may: +/// Tangent, as `sin(angle) / cos(angle)`. /// -/// - Saturate to `T::MAX` or `T::MIN` for very small denominators -/// - Produce very large finite values that may overflow in subsequent operations -/// -/// The threshold for potential overflow is approximately: -/// - `|angle - π/2| < 2^(-frac_bits/2)` for the nearest pole -/// -/// For I16F16, this means angles within ~0.004 radians of π/2 may overflow. -/// For I32F32, within ~0.00003 radians. -/// -/// If you need to detect near-pole conditions, check `cos(angle).abs()` -/// against a threshold before calling `tan`. -/// -/// # Example +/// It saturates to `T::MAX` or `T::MIN` within about 4e-5 of a pole for +/// I16F16 and 5e-10 for I32F32. Check the cosine to catch near-pole angles: /// /// ``` /// use fixed::types::I16F16; @@ -157,39 +131,43 @@ pub fn cos(angle: T) -> T { /// /// let angle = I16F16::from_num(1.5); // Close to π/2 ≈ 1.571 /// -/// // Safe: check cosine magnitude first /// let c = cos(angle); /// if c.abs() > I16F16::from_num(0.01) { /// let t = tan(angle); -/// // Use t safely /// } else { /// // Handle near-pole case /// } /// ``` #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn tan(angle: T) -> T { let (s, c) = sin_cos(angle); s.div(c) } -/// Arcsine. Domain: `[-1, 1]`. Returns angle in `[-π/2, π/2]`. -/// +/// Arcsine. Domain: `[-1, 1]`. Returns an angle in `[-π/2, π/2]`. /// # Errors /// Returns `DomainError` if `|x| > 1`. -#[allow( - clippy::inline_always, - reason = "inlined into acos so the no-panic link check sees only nounwind callees" -)] -#[inline(always)] #[must_use = "returns the arcsine result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn asin(x: T) -> Result { + asin_unguarded(x) +} + +/// [`asin`] without the no-panic guard, which [`acos`] calls so the link +/// check sees no guarded callee inside it. +#[inline] +fn asin_unguarded(x: T) -> Result { let Some(unit_x) = UnitInterval::new(x) else { return Err(Error::domain("asin", "value in range [-1, 1]")); }; - // Special cases if x == T::one() { return Ok(T::frac_pi_2()); } @@ -200,13 +178,11 @@ pub fn asin(x: T) -> Result { return Ok(T::zero()); } - // Use the identity: asin(x) = atan(x / sqrt(1 - x²)) - // NonNegative::one_minus_square gives 1 - x², which is ≥ 0 since |x| ≤ 1 + // asin(x) = atan(x / √(1 − x²)). let sqrt_term = sqrt_nonneg(NonNegative::one_minus_square(unit_x)); - // Handle case where sqrt_term is very small (x close to ±1) + // A root below 2^-15 puts x all but at ±1, so the result is ±π/2. if sqrt_term < T::from_i1f63(0x0001_0000_0000_0000) { - // Very close to ±1, return ±π/2 return if x.is_positive() { Ok(T::frac_pi_2()) } else { @@ -217,31 +193,33 @@ pub fn asin(x: T) -> Result { Ok(atan(x.div(sqrt_term))) } -/// Arccosine. Domain: `[-1, 1]`. Returns angle in `[0, π]`. -/// +/// Arccosine. Domain: `[-1, 1]`. Returns an angle in `[0, π]`. /// # Errors /// Returns `DomainError` if `|x| > 1`. #[must_use = "returns the arccosine result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn acos(x: T) -> Result { - // acos(x) = π/2 - asin(x) - asin(x).map(|a| T::frac_pi_2().saturating_sub(a)) + asin_unguarded(x).map(|a| T::frac_pi_2().saturating_sub(a)) } -/// Arctangent. Accepts any value. Returns angle in `(-π/2, π/2)`. +/// Arctangent of any value. Returns an angle in `(-π/2, π/2)`. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn atan(x: T) -> T { let zero = T::zero(); let one = T::one(); - // Special cases if x == zero { return zero; } - // For |x| > 1, use atan(x) = sign(x) * π/2 - atan(1/x) - // This keeps the argument in the convergent range + // Past 1, atan(x) = sign(x)·π/2 − atan(1/x) keeps the CORDIC convergent. let abs_x = x.abs(); if abs_x > one { let recip = one.div(x); @@ -253,25 +231,26 @@ pub fn atan(x: T) -> T { -T::frac_pi_2() - atan_recip } } else { - // |x| <= 1, use CORDIC directly circular_vectoring(one, x, zero).2 } } -/// Four-quadrant arctangent. Returns angle in `[-π, π]`. Returns 0 for (0, 0). +/// Four-quadrant arctangent. Returns an angle in `[-π, π]`, and 0 for (0, 0). #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn atan2(y: T, x: T) -> T { let zero = T::zero(); let pi = T::pi(); let frac_pi_2 = T::frac_pi_2(); - // Handle special cases if x == zero { return if y.is_negative() { -frac_pi_2 } else if y == zero { - zero // Undefined, but return 0 + zero } else { frac_pi_2 }; @@ -281,19 +260,12 @@ pub fn atan2(y: T, x: T) -> T { return if x.is_negative() { pi } else { zero }; } - // Compute atan(|y|/|x|) using CORDIC vectoring mode - // Using absolute values ensures the base angle is always positive let (_, _, base_angle) = circular_vectoring(x.abs(), y.abs(), zero); - // Adjust for quadrant based on signs of original x and y match (x.is_negative(), y.is_negative()) { - // Q1: x positive, y positive -> angle is base_angle (false, false) => base_angle, - // Q4: x positive, y negative -> angle is -base_angle (false, true) => -base_angle, - // Q2: x negative, y positive -> angle is π - base_angle (true, false) => pi - base_angle, - // Q3: x negative, y negative -> angle is -(π - base_angle) = base_angle - π (true, true) => base_angle - pi, } } diff --git a/src/ops/exponential.rs b/src/ops/exponential.rs index 77bfc01..326722f 100644 --- a/src/ops/exponential.rs +++ b/src/ops/exponential.rs @@ -1,4 +1,4 @@ -//! Exponential and logarithmic functions. +#![doc = include_str!("../../docs/exponential.md")] use crate::bounded::{NormalizedLnArg, OpenUnitInterval}; use crate::error::{Error, Result}; @@ -8,18 +8,14 @@ use crate::traits::CordicNumber; /// `1/ln 2 − 1` as I1F63 (`1/ln 2` itself exceeds 1). const FRAC_1_LN2_MINUS_1_I1F63: i64 = 0x38AA_3B29_5C17_F0BC; -/// Right-shifts a non-negative value by `n`, rounding to nearest. -/// -/// A plain `>>` truncates, which biases results low by up to a full ulp. -/// That matters for the final scaling of `exp`, where the true result may -/// be only a few ulps: rounding halves the worst case and removes the bias. +/// Right-shifts a non-negative value by `n`, rounding to nearest. A plain `>>` +/// would bias small `exp` results, only a few ulps, low by up to an ulp. fn shr_rounded(x: T, n: u32) -> T { if n == 0 { return x; } let ulp = T::one() >> T::frac_bits(); let shifted = x >> n; - // Reconstruct the discarded low bits: rem ∈ [0, 2^n) ulps. let rem = x.saturating_sub(shifted << n); let half = ulp << (n - 1); if rem >= half { @@ -29,60 +25,27 @@ fn shr_rounded(x: T, n: u32) -> T { } } -/// Exponential function (e^x). -/// -/// # Saturation Behavior -/// -/// This function saturates for extreme inputs rather than returning an error: -/// -/// | Condition | Result | Example (I16F16) | -/// |-----------|--------|------------------| -/// | x > `ln(T::MAX)` | `T::MAX` | x > ~10.4 → 32767.99 | -/// | x < `ln(T::MIN_POSITIVE / 2)` | `T::ZERO` | x < ~-11.8 → 0 | -/// -/// The zero threshold is where e^x rounds to nearest below half an ulp. -/// The exact thresholds depend on the type's range: -/// - **I16F16:** Saturates to MAX for x > ~10.4, to zero for x < ~-11.8 -/// - **I32F32:** Saturates to MAX for x > ~21.5, to zero for x < ~-22.9 -/// -/// Saturation is silent and deterministic. If you need to detect overflow, -/// check the input range before calling: -/// -/// ``` -/// use fixed::types::I16F16; -/// use fixed_analytics::exp; -/// -/// let x = I16F16::from_num(5.0); -/// let max_safe = I16F16::from_num(10.0); -/// -/// if x < max_safe { -/// let result = exp(x); // Safe -/// } else { -/// // Handle potential saturation -/// } -/// ``` +/// Exponential function e^x. Saturates to `T::MAX` above `ln(T::MAX)` and to +/// zero far below 0; the [module docs](self) list the thresholds. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn exp(x: T) -> T { let zero = T::zero(); let one = T::one(); let ln2 = T::ln_2(); - // Handle special case if x == zero { return one; } - // Argument reduction: exp(x) = 2^k * exp(r), where r ∈ [0, ln2), so - // exp(r) ∈ [1, 2) carries a full significand into the scaling shift. - // k = ⌊x / ln2⌋ is estimated to within one (`to_i32` floors; the - // reciprocal multiply is within |x|·1.5·2^-frac + 2^-frac < 1 when the - // integer bits do not outnumber the fractional bits, and other types - // divide) and then corrected once in either direction. - #[allow(clippy::cast_possible_wrap, reason = "bit counts bounded by type size")] + // exp(x) = 2^k·exp(r) with r in [0, ln 2). The estimate of k = ⌊x / ln 2⌋ + // is within one, so one correction below fixes it (docs/exp.md). let (int_bits, frac_bits) = ( - (T::total_bits() - T::frac_bits()) as i32, - T::frac_bits() as i32, + (T::total_bits() - T::frac_bits()).cast_signed(), + T::frac_bits().cast_signed(), ); let quotient = if T::total_bits() <= 2 * T::frac_bits() { x.saturating_add(x.saturating_mul(T::from_i1f63(FRAC_1_LN2_MINUS_1_I1F63))) @@ -91,8 +54,8 @@ pub fn exp(x: T) -> T { }; let mut scale = quotient.to_i32(); - // 2^k·exp(r) exceeds MAX once k ≥ int_bits − 1 and rounds to zero once - // k ≤ −(frac_bits + 2); allow one of slack until k is corrected. + // The result exceeds MAX once k ≥ int_bits − 1 and rounds to zero once + // k ≤ −(frac_bits + 2); these exits allow one of slack before correction. if scale >= int_bits { return T::max_value(); } @@ -115,10 +78,7 @@ pub fn exp(x: T) -> T { return zero; } - // Factored Taylor: exp(r) = 1 + r*(1 + r/2*(1 + r/3*(1 + ... r/n))). - // The omitted term (ln2)^(n+1)/(n+1)! is 1.3e-6 at degree 7, 1.4e-12 at - // degree 12 and 2.7e-22 at degree 19: below half an ulp up to 18, 38 - // and 70 fractional bits respectively. + // Factored Taylor series, of degree 7, 12 or 19 by precision (docs/exp.md). let mut p = one; if T::frac_bits() >= 40 { p = one.saturating_add(r.div_int(19).saturating_mul(p)); @@ -144,16 +104,9 @@ pub fn exp(x: T) -> T { p = one.saturating_add(r.div_int(2).saturating_mul(p)); let exp_r = one.saturating_add(r.saturating_mul(p)); - // Scale by 2^scale using bit shifts. - // scale is already bounded to (-(frac_bits + 2), int_bits - 1) by the exits above. - #[allow( - clippy::cast_sign_loss, - reason = "sign of scale checked before each cast" - )] if scale > 0 { - let shift = scale as u32; - // Detect overflow before shifting: if exp_r > MAX >> shift, - // the left shift would wrap, so saturate to MAX instead. + let shift = scale.unsigned_abs(); + // A left shift past the headroom would wrap. let headroom = T::max_value() >> shift; if exp_r > headroom { T::max_value() @@ -161,24 +114,22 @@ pub fn exp(x: T) -> T { exp_r << shift } } else { - // scale ≤ 0: rounded right shift (shr_rounded(x, 0) is the identity) - shr_rounded(exp_r, (-scale) as u32) + shr_rounded(exp_r, scale.unsigned_abs()) } } -/// Power function `base^exponent`, computed as `exp(exponent · ln|base|)`. -/// -/// Domain: `base > 0`; `base = 0` with `exponent ≥ 0`; or `base < 0` with -/// an integer exponent, where the sign follows the exponent's parity. -/// The logarithm's rounding error is multiplied by the exponent, so the -/// relative error grows with `|exponent|`. Saturates like [`exp`]. +/// Power `base^exponent`, as `exp(exponent · ln|base|)`, saturating like [`exp`]. /// +/// A negative base takes its sign from the exponent's parity. The relative +/// error grows with `|exponent|`, which scales the logarithm's rounding error. /// # Errors -/// Returns `DomainError` for a negative base with a non-integer exponent -/// (no real result, as for [`sqrt`](crate::sqrt)) or a zero base with a -/// negative exponent. +/// Returns `DomainError` for a zero base with a negative exponent, or a +/// negative base with a non-integer exponent, which has no real result. #[must_use = "returns the power result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn pow(base: T, exponent: T) -> Result { let zero = T::zero(); if exponent == zero { @@ -198,7 +149,7 @@ pub fn pow(base: T, exponent: T) -> Result { if !base.is_negative() { return Ok(magnitude); } - // Shifting the raw bits isolates the integer part exactly at any magnitude. + // The shift isolates the integer part exactly at any magnitude. let integer_part = exponent >> T::frac_bits(); if integer_part << T::frac_bits() != exponent { return Err(Error::domain("pow", "integer exponent for a negative base")); @@ -208,11 +159,13 @@ pub fn pow(base: T, exponent: T) -> Result { } /// Natural logarithm. Domain: `x > 0`. -/// /// # Errors /// Returns `DomainError` if `x ≤ 0`. #[must_use = "returns the natural logarithm result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn ln(x: T) -> Result { if x <= T::zero() { return Err(Error::domain("ln", "positive value")); @@ -230,55 +183,43 @@ pub(crate) fn ln_positive(x: T) -> T { return zero; } - // For x far from 1, use argument reduction: - // ln(x) = ln(x * 2^(-k)) + k * ln(2) - // where k is chosen so that x * 2^(-k) is close to 1 - let ln2 = T::ln_2(); - // Reduce to [0.5, 2] from the leading bit e = ⌊log₂ x⌋: x > 2 shifts - // right until ≤ 2 (stopping early when a shift lands exactly on 2), - // x < 0.5 doubles until ≥ 0.5. + // ln(x) = ln(x·2^-k) + k·ln 2, with x·2^-k in [0.5, 2] chosen from the + // leading bit e = ⌊log₂ x⌋; a value landing exactly on 2 stops there. let e = x.checked_int_log2().unwrap_or(0); - #[allow( - clippy::cast_sign_loss, - reason = "shift counts are non-negative by construction" - )] let (normalized, k) = if e >= 1 { - let candidate = x >> (e - 1) as u32; + let candidate = x >> (e - 1).unsigned_abs(); if candidate == two { (candidate, e - 1) } else { - (x >> e as u32, e) + (x >> e.unsigned_abs(), e) } } else if e <= -2 { - (x << (-1 - e) as u32, e + 1) + (x << (-1 - e).unsigned_abs(), e + 1) } else { (x, 0) }; let k_ln2 = ln2.mul_int(k); - // Now compute ln(normalized) where 0.5 <= normalized <= 2 - // Using ln(x) = 2 * atanh((x-1)/(x+1)) - // NormalizedLnArg encodes that normalized ∈ [0.5, 2] + // ln(x) = 2·atanh((x − 1)/(x + 1)), whose argument is within 1/3 of 0. let norm = NormalizedLnArg::from_normalized(normalized); - - // OpenUnitInterval::from_normalized_ln_arg computes (x-1)/(x+1), - // which is in (-1/3, 1/3) ⊂ (-1, 1) for x ∈ [0.5, 2] let arg = OpenUnitInterval::from_normalized_ln_arg(norm); let atanh_val = atanh_open(arg); - let ln_normalized = atanh_val.saturating_add(atanh_val); // 2 * atanh + let ln_normalized = atanh_val.saturating_add(atanh_val); ln_normalized.saturating_add(k_ln2) } /// Base-2 logarithm. Domain: `x > 0`. -/// /// # Errors /// Returns `DomainError` if `x ≤ 0`. #[must_use = "returns the base-2 logarithm result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn log2(x: T) -> Result { let ln_x = ln(x)?; let ln_2 = T::ln_2(); @@ -286,34 +227,27 @@ pub fn log2(x: T) -> Result { } /// Base-10 logarithm. Domain: `x > 0`. -/// /// # Errors /// Returns `DomainError` if `x ≤ 0`. #[must_use = "returns the base-10 logarithm result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn log10(x: T) -> Result { let ln_x = ln(x)?; let ln_10 = T::ln_10(); Ok(ln_x.div(ln_10)) } -/// Power of 2 (2^x). Computed as exp(x × ln(2)). -/// -/// # Saturation Behavior -/// -/// Saturates for extreme inputs: -/// -/// | Condition | Result | Example (I16F16) | -/// |-----------|--------|------------------| -/// | x > `log2(T::MAX)` | `T::MAX` | x > ~15 → 32767.99 | -/// | x < `log2(T::MIN_POSITIVE / 2)` | `T::ZERO` | x < ~-17 → 0 | -/// -/// The zero threshold is where 2^x rounds to nearest below half an ulp. -/// The exact thresholds: -/// - **I16F16:** Saturates for x > ~15 or x < ~-17 -/// - **I32F32:** Saturates for x > ~31 or x < ~-33 +/// Power of two 2^x, as `exp(x · ln 2)`. Saturates to `T::MAX` above +/// `log2(T::MAX)` and to zero far below 0; the [module docs](self) list the +/// thresholds. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn pow2(x: T) -> T { let ln_2 = T::ln_2(); exp(x.saturating_mul(ln_2)) diff --git a/src/ops/hyperbolic.rs b/src/ops/hyperbolic.rs index 54df348..b4a7176 100644 --- a/src/ops/hyperbolic.rs +++ b/src/ops/hyperbolic.rs @@ -1,4 +1,4 @@ -//! Hyperbolic functions via hyperbolic CORDIC. +#![doc = include_str!("../../docs/hyperbolic.md")] use crate::bounded::{AtLeastOne, NonNegative, OpenUnitInterval, UnitInterval}; use crate::error::{Error, Result}; @@ -7,68 +7,38 @@ use crate::ops::algebraic::sqrt_nonneg; use crate::ops::exponential::ln_positive; use crate::traits::CordicNumber; -/// Hyperbolic CORDIC converges for |x| < sum of atanh table ≈ 1.1182. -/// Stored as fractional part (0.1182) since I1F63 can't hold 1.x. -/// Full limit = 1 + this value. +/// The fractional part of 1.1182, the sum of the atanh table, which bounds the +/// polynomial's argument; I1F63 cannot hold 1.1182 itself. const HYPERBOLIC_CONVERGENCE_LIMIT_FRAC_I1F63: i64 = 0x0F22_3D70_A3D7_0A3D; -/// Argument reduction threshold for atanh. -/// -/// For |x| > this threshold, atanh uses the identity: -/// ```text -/// atanh(x) = atanh(0.5) + atanh((x - 0.5) / (1 - 0.5x)) -/// ``` -/// to reduce the argument into CORDIC's optimal convergence range. -/// -/// Value 0.75 keeps reduced arguments within the convergent region. +/// 0.75, above which atanh moves its argument into the CORDIC's convergent +/// range with `atanh(x) = atanh(0.5) + atanh((x − 0.5) / (1 − 0.5x))`. const ATANH_REDUCTION_THRESHOLD_I1F63: i64 = 0x6000_0000_0000_0000; -/// Hyperbolic sine and cosine. More efficient than separate calls. -/// -/// # Saturation Behavior -/// -/// Both outputs can saturate for large |x|: -/// - sinh saturates to `T::MAX` or `T::MIN` -/// - cosh saturates to `T::MAX` -/// -/// See [`sinh`] and [`cosh`] for threshold details. -/// -/// When saturation occurs, both values saturate together (they grow -/// at the same rate), so the relationship cosh²(x) - sinh²(x) = 1 -/// will not hold for saturated outputs. +/// Hyperbolic sine and cosine, cheaper than separate calls. Both saturate +/// together for large `|x|`, as the [module docs](self) describe. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sinh_cosh(x: T) -> (T, T) { let one = T::one(); - // Compute limit as 1 + fractional_part (~1.1182) let limit = one.saturating_add(T::from_i1f63(HYPERBOLIC_CONVERGENCE_LIMIT_FRAC_I1F63)); - // Iterative argument reduction for large values. - // Count how many halvings are needed, then rebuild on the way back. + // Halve until within the limit, then double back; total_bits halvings + // reach zero, so the bound never cuts the loop short. let mut reduced = x; let mut depth: u32 = 0; - // Max depth bounded by bit width: each halving shifts right by 1, - // so after total_bits iterations the value is zero and below the limit. while reduced.abs() > limit && depth < T::total_bits() { reduced = reduced >> 1; depth += 1; } - // Factored Taylor evaluation. After argument reduction, |reduced| ≤ 1.1182. - // - // sinh(x) = x * (1 + u/6 * (1 + u/20 * (1 + u/42 * ...))) - // cosh(x) = 1 + u/2 * (1 + u/12 * (1 + u/30 * ...)) - // where u = x². - // - // Integer division (u/K) is used because it computes the quotient in one - // step without pre-rounding a reciprocal, and the factored form keeps each - // step's multiplicand u/K below 1 for K ≥ 2. + // Factored Taylor series in u = x², with divisors (2k)(2k+1) for sinh and + // (2k−1)(2k) for cosh; docs/sinh_cosh.md gives the degrees and errors. let u = reduced.saturating_mul(reduced); - // Divisors are (2k)(2k+1) for sinh and (2k-1)(2k) for cosh. At |x| ≤ - // 1.1182 the omitted relative term is 7.7e-8 / 8.0e-9 at degrees 9 / 10, - // 3.7e-12 / 2.9e-13 at 13 / 14 and 4.5e-22 / 2.4e-23 at 21 / 22: below - // half an ulp up to 22, 37 and 70 fractional bits respectively. let mut sp = one; let mut cp = one; if T::frac_bits() >= 40 { @@ -98,8 +68,7 @@ pub fn sinh_cosh(x: T) -> (T, T) { cp = one.saturating_add(u.div_int(2).saturating_mul(cp)); let (mut sh, mut ch) = (reduced.saturating_mul(sp), cp); - // Reconstruct via doubling: sinh(2x) = 2·sinh(x)·cosh(x), - // cosh(2x) = cosh²(x) + sinh²(x) + // sinh(2x) = 2·sinh(x)·cosh(x) and cosh(2x) = cosh²(x) + sinh²(x). for _ in 0..depth { let new_sh = sh.saturating_mul(ch).saturating_mul(T::two()); let new_ch = ch.saturating_mul(ch).saturating_add(sh.saturating_mul(sh)); @@ -110,68 +79,49 @@ pub fn sinh_cosh(x: T) -> (T, T) { (sh, ch) } -/// Hyperbolic sine. -/// -/// # Saturation Behavior -/// -/// sinh(x) grows exponentially: sinh(x) ≈ e^x/2 for large |x|. -/// This function saturates for extreme inputs: -/// -/// | Condition | Result | Example (I16F16) | -/// |-----------|--------|------------------| -/// | x > `asinh(T::MAX)` | `T::MAX` | x > ~11.1 → 32767.99 | -/// | x < `-asinh(T::MAX)` | `T::MIN` | x < ~-11.1 → -32768.0 | -/// -/// The exact thresholds: -/// - **I16F16:** Saturates for |x| > ~11.1 -/// - **I32F32:** Saturates for |x| > ~22.2 -/// -/// Within the non-saturating range, sinh is computed via polynomial -/// evaluation with argument reduction (halving/doubling) for |x| > 1.118. +/// Hyperbolic sine. Saturates to `T::MAX` or `T::MIN` past `±asinh(T::MAX)`; +/// the [module docs](self) list the thresholds. #[inline] #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn sinh(x: T) -> T { sinh_cosh(x).0 } -/// Hyperbolic cosine. Always ≥ 1. -/// -/// # Saturation Behavior -/// -/// cosh(x) grows exponentially: cosh(x) ≈ e^|x|/2 for large |x|. -/// This function saturates for extreme inputs: -/// -/// | Condition | Result | Example (I16F16) | -/// |-----------|--------|------------------| -/// | \|x\| > `acosh(T::MAX)` | `T::MAX` | \|x\| > ~11.1 → 32767.99 | -/// -/// Unlike sinh, cosh is always positive, so it only saturates to `T::MAX`. -/// -/// The exact thresholds: -/// - **I16F16:** Saturates for |x| > ~11.1 -/// - **I32F32:** Saturates for |x| > ~22.2 +/// Hyperbolic cosine, always at least 1. Saturates to `T::MAX` past +/// `acosh(T::MAX)`; the [module docs](self) list the thresholds. #[inline] #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn cosh(x: T) -> T { sinh_cosh(x).1 } /// Hyperbolic tangent. Result in `(-1, 1)`. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn tanh(x: T) -> T { let (s, c) = sinh_cosh(x); s.div(c) } /// Hyperbolic cotangent. Domain: `x ≠ 0`. -/// /// # Errors /// Returns `DomainError` if `x = 0`. #[must_use = "returns the hyperbolic cotangent result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn coth(x: T) -> Result { if x == T::zero() { return Err(Error::domain("coth", "non-zero value")); @@ -182,7 +132,10 @@ pub fn coth(x: T) -> Result { /// Inverse hyperbolic sine. Accepts any value. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn asinh(x: T) -> T { let zero = T::zero(); let one = T::one(); @@ -193,21 +146,18 @@ pub fn asinh(x: T) -> T { let abs_x = x.abs(); if abs_x <= one { - // asinh(x) = atanh(x / sqrt(1 + x²)); the argument is at most 1/√2, - // inside the CORDIC's direct range, and 1 + x² is always ≥ 1. + // asinh(x) = atanh(x / √(1 + x²)), whose argument is at most 1/√2. let sqrt_term = sqrt_nonneg(NonNegative::one_plus_square(x)); let arg = OpenUnitInterval::from_div_by_sqrt_one_plus_square(x, sqrt_term); return atanh_open(arg); } - // Beyond 1 that argument approaches 1, where each step of atanh's - // argument reduction triples the error, so use - // asinh(x) = sign(x)·ln(|x|·(1 + sqrt(1 + 1/x²))), - // with 1/x² formed as (1/|x|)², which cannot overflow. + // Past 1 that argument nears 1, where each atanh reduction step triples + // the error, so use sign(x)·ln(|x|·(1 + √(1 + (1/|x|)²))). let recip = one.div(abs_x); let factor = one.saturating_add(sqrt_nonneg(NonNegative::one_plus_square(recip))); let magnitude = if abs_x > T::max_value() >> 2 { - // |x|·factor would overflow (factor < 2.42): split the logarithm. + // |x|·factor would overflow, as factor < 2.42. ln_positive(abs_x).saturating_add(ln_positive(factor)) } else { ln_positive(abs_x.saturating_mul(factor)) @@ -220,11 +170,13 @@ pub fn asinh(x: T) -> T { } /// Inverse hyperbolic cosine. Domain: `x ≥ 1`. -/// /// # Errors /// Returns `DomainError` if `x < 1`. #[must_use = "returns the inverse hyperbolic cosine result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn acosh(x: T) -> Result { let at_least_one = AtLeastOne::new(x).ok_or_else(|| Error::domain("acosh", "value >= 1"))?; let zero = T::zero(); @@ -235,21 +187,18 @@ pub fn acosh(x: T) -> Result { } if x < one.saturating_add(T::half()) { - // acosh(x) = atanh(sqrt(x² - 1) / x); below 1.5 the argument is - // under 0.75, inside the CORDIC's direct range, and x² - 1 ≥ 0. + // acosh(x) = atanh(√(x² − 1) / x), whose argument is under 0.75 here. let sqrt_term = sqrt_nonneg(NonNegative::square_minus_one(at_least_one)); let arg = OpenUnitInterval::from_sqrt_square_minus_one_div(sqrt_term, at_least_one); return Ok(atanh_open(arg)); } - // Beyond 1.5 that argument approaches 1, where each step of atanh's - // argument reduction triples the error, so use - // acosh(x) = ln(x·(1 + sqrt(1 - 1/x²))), - // with 1/x² formed as (1/x)², which cannot overflow. + // Past 1.5 that argument nears 1, where each atanh reduction step triples + // the error, so use ln(x·(1 + √(1 − (1/x)²))). let recip = UnitInterval::from_reciprocal(at_least_one); let factor = one.saturating_add(sqrt_nonneg(NonNegative::one_minus_square(recip))); let result = if x > T::max_value() >> 1 { - // x·factor would overflow (factor < 2): split the logarithm. + // x·factor would overflow, as factor < 2. ln_positive(x).saturating_add(ln_positive(factor)) } else { ln_positive(x.saturating_mul(factor)) @@ -258,31 +207,31 @@ pub fn acosh(x: T) -> Result { } /// Inverse hyperbolic tangent. Domain: `(-1, 1)`. -/// /// # Errors /// Returns `DomainError` if `|x| ≥ 1`. #[must_use = "returns the inverse hyperbolic tangent result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn atanh(x: T) -> Result { OpenUnitInterval::new(x) .map(atanh_open) .ok_or_else(|| Error::domain("atanh", "value in range (-1, 1)")) } -/// Infallible inverse hyperbolic tangent for values in (-1, 1). -/// -/// This function takes an [`OpenUnitInterval`] wrapper, guaranteeing at the -/// type level that the input is valid. No domain check is performed at runtime. -/// -/// Use this when the input is known to be in (-1, 1) through mathematical -/// invariants (e.g., `x / sqrt(1 + x²)`). +/// Inverse hyperbolic tangent of an [`OpenUnitInterval`] value, with no +/// runtime domain check. #[must_use] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn atanh_open(x: OpenUnitInterval) -> T { atanh_core(x.get()) } -/// Core atanh implementation. Caller must ensure |x| < 1. +/// atanh for `|x| < 1`, which the caller must ensure. fn atanh_core(x: T) -> T { let zero = T::zero(); let one = T::one(); @@ -293,25 +242,21 @@ fn atanh_core(x: T) -> T { let threshold = T::from_i1f63(ATANH_REDUCTION_THRESHOLD_I1F63); - // Fast path: no argument reduction needed, use CORDIC directly. if x.abs() <= threshold { let (_, _, z) = hyperbolic_vectoring(one, x, zero); return z; } - // Argument reduction needed. Work with |x| and track sign. let half = T::half(); let atanh_half = T::from_i1f63(crate::tables::hyperbolic::ATANH_HALF); let sign = if x.is_negative() { -one } else { one }; let mut abs_x = x.abs(); - // Iterative argument reduction: each step reduces |x| and accumulates - // atanh(0.5) into the result. Max iterations bounded by frac_bits since - // each reduction roughly halves the distance from the threshold. + // Each step takes atanh(0.5) off and roughly halves the distance to the + // threshold, so frac_bits steps suffice. let mut accumulated = zero; let mut i: u32 = 0; while abs_x > threshold && i < T::frac_bits() { - // atanh(x) = atanh(0.5) + atanh((x - 0.5) / (1 - 0.5*x)) let numerator = abs_x.saturating_sub(half); let denominator = one.saturating_sub(half.saturating_mul(abs_x)); abs_x = numerator.div(denominator); @@ -319,18 +264,19 @@ fn atanh_core(x: T) -> T { i += 1; } - // Direct CORDIC computation on the reduced argument let (_, _, z) = hyperbolic_vectoring(one, abs_x, zero); sign.saturating_mul(accumulated.saturating_add(z)) } /// Inverse hyperbolic cotangent. Domain: `|x| > 1`. -/// /// # Errors /// Returns `DomainError` if `|x| ≤ 1`. #[must_use = "returns the inverse hyperbolic cotangent result which should be handled"] -#[cfg_attr(feature = "verify-no-panic", no_panic::no_panic)] +#[cfg_attr( + all(feature = "verify-no-panic", not(debug_assertions)), + no_panic::no_panic +)] pub fn acoth(x: T) -> Result { let one = T::one(); @@ -338,7 +284,6 @@ pub fn acoth(x: T) -> Result { return Err(Error::domain("acoth", "|value| > 1")); } - // acoth(x) = atanh(1/x) let recip = one.div(x); Ok(atanh_core(recip)) } diff --git a/src/ops/mod.rs b/src/ops/mod.rs index 9d7739b..f13260e 100644 --- a/src/ops/mod.rs +++ b/src/ops/mod.rs @@ -1,21 +1,10 @@ -//! High-level mathematical operations built on CORDIC kernels. -//! -//! This module provides user-friendly functions for computing various -//! mathematical operations using the CORDIC algorithm. -//! -//! # Modules -//! -//! - [`circular`]: Trigonometric functions (sin, cos, tan, asin, acos, atan, atan2) -//! - [`hyperbolic`]: Hyperbolic functions (sinh, cosh, tanh, asinh, acosh, atanh, acoth) -//! - [`exponential`]: Exponential and logarithmic functions (exp, ln, log2, log10, pow, pow2) -//! - [`algebraic`]: Algebraic functions (sqrt) +//! The mathematical functions, each also re-exported at the crate root. pub mod algebraic; pub mod circular; pub mod exponential; pub mod hyperbolic; -// Re-export all public functions pub use algebraic::sqrt; pub use circular::{acos, asin, atan, atan2, cos, sin, sin_cos, tan}; pub use exponential::{exp, ln, log2, log10, pow, pow2}; diff --git a/src/tables/chebyshev.rs b/src/tables/chebyshev.rs index 8424c30..b39f665 100644 --- a/src/tables/chebyshev.rs +++ b/src/tables/chebyshev.rs @@ -1,38 +1,14 @@ -//! Minimax (Chebyshev) polynomial coefficients for sin/cos evaluation. -//! -//! Generated via the Remez algorithm (`mpmath.chebyfit`). Each coefficient -//! approximates a reduced-form remainder, -//! factored so all values fit in I1F63 (magnitude < 1). -//! -//! Stored highest-degree first: `[cₙ, cₙ₋₁, …, c₀]`. -//! Evaluated via Horner's method: -//! -//! ```text -//! P(x) = cₙ·xⁿ + ··· + c₁·x + c₀ = c₀ + x·(c₁ + x·(c₂ + ··· + x·cₙ)) -//! ``` -//! -//! # Reconstruction -//! -//! | Function | Factored form | Variable | Domain | -//! |----------|---------------|----------|--------| -//! | sin | `x + x³·P(x²)` | u = x² | \[0, (π/4)²\] | -//! | cos | `1 + x²·Q(x²)` | u = x² | \[0, (π/4)²\] | +#![doc = include_str!("../../docs/chebyshev.md")] use crate::traits::CordicNumber; -/// Evaluate polynomial via Horner's method using precomputed I1F63 coefficients. +/// Evaluates I1F63 coefficients, highest degree first, by Horner's method. /// -/// Stored highest-degree first: `[cₙ, cₙ₋₁, …, c₀]`. -/// Evaluates P(x) = cₙ·xⁿ + cₙ₋₁·xⁿ⁻¹ + ··· + c₀ via Horner's method. -/// -/// Uses wrapping arithmetic, so the caller must keep the evaluation in -/// range: for the tables in this module and `0 ≤ x ≤ (π/4)²` every -/// intermediate stays below 1 in magnitude, so wrapping and saturating -/// arithmetic are bit-identical and the overflow checks would be pure cost. +/// The arithmetic wraps. For these tables and `0 ≤ x ≤ (π/4)²` every +/// intermediate stays below 1, where that matches saturating arithmetic. #[inline] pub fn horner(coeffs: &[i64; N], x: T) -> T { let mut iter = coeffs.iter(); - // First element is the highest-degree coefficient (N ≥ 3 for all tables). let mut result = T::from_i1f63(*iter.next().unwrap_or(&0)); for &coeff in iter { result = T::from_i1f63(coeff).wrapping_add(x.wrapping_mul(result)); @@ -40,52 +16,44 @@ pub fn horner(coeffs: &[i64; N], x: T) -> T { result } -// ============================================================================= -// Sin: P(u) = (sin(x) - x) / x³, u = x², domain [0, (π/4)²] -// ============================================================================= - -/// Minimax coeffs for (sin(x)-x)/x³ in u=x² on [0,(π/4)²]. Low precision (I16F16-class). +/// (sin(x) − x)/x³ in u = x², low precision. #[rustfmt::skip] pub const SIN_P_LO: [i64; 4] = [ - 0x0000_16DB_E083_8A07, // +2.725e-06 - -0x0006_804E_9E4E_E633, // -1.984e-04 - 0x0111_110D_EF2E_1C96, // +8.333e-03 - -0x1555_5555_45E0_ABF9, // -1.667e-01 + 0x0000_16DB_E083_8A07, + -0x0006_804E_9E4E_E633, + 0x0111_110D_EF2E_1C96, + -0x1555_5555_45E0_ABF9, ]; -/// Minimax coeffs for (sin(x)-x)/x³ in u=x² on [0,(π/4)²]. High precision (I32F32-class). +/// (sin(x) − x)/x³ in u = x², high precision. #[rustfmt::skip] pub const SIN_P_HI: [i64; 7] = [ - -0x0000_0000_006A_C5F5, // -7.587e-13 - 0x0000_0000_5848_5FC3, // +1.606e-10 - -0x0000_0035_CC8A_8259, // -2.505e-08 - 0x0000_171D_E3A5_4609, // +2.756e-06 - -0x0006_8068_0680_664E, // -1.984e-04 - 0x0111_1111_1111_1100, // +8.333e-03 - -0x1555_5555_5555_5555, // -1.667e-01 + -0x0000_0000_006A_C5F5, + 0x0000_0000_5848_5FC3, + -0x0000_0035_CC8A_8259, + 0x0000_171D_E3A5_4609, + -0x0006_8068_0680_664E, + 0x0111_1111_1111_1100, + -0x1555_5555_5555_5555, ]; -// ============================================================================= -// Cos: Q(u) = (cos(x) - 1) / x², u = x², domain [0, (π/4)²] -// ============================================================================= - -/// Minimax coeffs for (cos(x)-1)/x² in u=x² on [0,(π/4)²]. Low precision. +/// (cos(x) − 1)/x² in u = x², low precision. #[rustfmt::skip] pub const COS_Q_LO: [i64; 4] = [ - 0x0000_CD37_95D8_0CFA, // +2.446e-05 - -0x002D_81C1_0FEA_FDD5, // -1.389e-03 - 0x0555_5532_ED10_3C5F, // +4.167e-02 - -0x3FFF_FFFF_563B_4B19, // -5.000e-01 + 0x0000_CD37_95D8_0CFA, + -0x002D_81C1_0FEA_FDD5, + 0x0555_5532_ED10_3C5F, + -0x3FFF_FFFF_563B_4B19, ]; -/// Minimax coeffs for (cos(x)-1)/x² in u=x² on [0,(π/4)²]. High precision. +/// (cos(x) − 1)/x² in u = x², high precision. #[rustfmt::skip] pub const COS_Q_HI: [i64; 7] = [ - -0x0000_0000_063F_E751, // -1.137e-11 - 0x0000_0004_7BA9_FC96, // +2.088e-09 - -0x0000_024F_C9F1_C946, // -2.756e-07 - 0x0000_D00D_00CE_F099, // +2.480e-05 - -0x002D_82D8_2D82_BAF2, // -1.389e-03 - 0x0555_5555_5555_5435, // +4.167e-02 - -0x3FFF_FFFF_FFFF_FFFE, // -5.000e-01 + -0x0000_0000_063F_E751, + 0x0000_0004_7BA9_FC96, + -0x0000_024F_C9F1_C946, + 0x0000_D00D_00CE_F099, + -0x002D_82D8_2D82_BAF2, + 0x0555_5555_5555_5435, + -0x3FFF_FFFF_FFFF_FFFE, ]; diff --git a/src/tables/circular.rs b/src/tables/circular.rs index 9e5b651..0b4a2d3 100644 --- a/src/tables/circular.rs +++ b/src/tables/circular.rs @@ -1,33 +1,30 @@ -//! Arctangent table for circular CORDIC. Values: `atan(2^-i)` as I1F63. -//! -//! Gain factor K = ∏ sqrt(1 + 2^(-2i)) ≈ 1.6468. We store 1/K ≈ 0.6073. -//! For i ≥ 21, atan(x) ≈ x to I1F63 precision. +//! Arctangent table for circular CORDIC. -/// `atan(2^-i)` as I1F63. Index 0 = atan(1) = π/4. +/// `atan(2^-i)` as I1F63 at index i, starting at atan(1) = π/4. From 2^-21 +/// on, each entry equals its argument to I1F63 precision. #[rustfmt::skip] pub const ATAN_TABLE: [i64; 64] = [ - 0x6487_ED51_10B4_611A, // atan(2^0) = π/4 - 0x3B58_CE0A_C376_9ED1, // atan(2^-1) - 0x1F5B_75F9_2C80_DD63, // atan(2^-2) - 0x0FEA_DD4D_5617_B6E3, // atan(2^-3) - 0x07FD_56ED_CB3F_7A72, // atan(2^-4) - 0x03FF_AAB7_752E_C495, // atan(2^-5) - 0x01FF_F555_BBB7_29AB, // atan(2^-6) - 0x00FF_FEAA_ADDD_D4B9, // atan(2^-7) - 0x007F_FFD5_556E_EEDD, // atan(2^-8) - 0x003F_FFFA_AAAB_7777, // atan(2^-9) - 0x001F_FFFF_5555_5BBC, // atan(2^-10) - 0x000F_FFFF_EAAA_AADE, // atan(2^-11) - 0x0007_FFFF_FD55_5557, // atan(2^-12) - 0x0003_FFFF_FFAA_AAAB, // atan(2^-13) - 0x0001_FFFF_FFF5_5555, // atan(2^-14) - 0x0000_FFFF_FFFE_AAAB, // atan(2^-15) - 0x0000_7FFF_FFFF_D555, // atan(2^-16) - 0x0000_3FFF_FFFF_FAAB, // atan(2^-17) - 0x0000_1FFF_FFFF_FF55, // atan(2^-18) - 0x0000_0FFF_FFFF_FFEB, // atan(2^-19) - 0x0000_07FF_FFFF_FFFD, // atan(2^-20) - // For i >= 21: atan(2^-i) ≈ 2^-i to I1F63 precision + 0x6487_ED51_10B4_611A, + 0x3B58_CE0A_C376_9ED1, + 0x1F5B_75F9_2C80_DD63, + 0x0FEA_DD4D_5617_B6E3, + 0x07FD_56ED_CB3F_7A72, + 0x03FF_AAB7_752E_C495, + 0x01FF_F555_BBB7_29AB, + 0x00FF_FEAA_ADDD_D4B9, + 0x007F_FFD5_556E_EEDD, + 0x003F_FFFA_AAAB_7777, + 0x001F_FFFF_5555_5BBC, + 0x000F_FFFF_EAAA_AADE, + 0x0007_FFFF_FD55_5557, + 0x0003_FFFF_FFAA_AAAB, + 0x0001_FFFF_FFF5_5555, + 0x0000_FFFF_FFFE_AAAB, + 0x0000_7FFF_FFFF_D555, + 0x0000_3FFF_FFFF_FAAB, + 0x0000_1FFF_FFFF_FF55, + 0x0000_0FFF_FFFF_FFEB, + 0x0000_07FF_FFFF_FFFD, 0x0000_0400_0000_0000, 0x0000_0200_0000_0000, 0x0000_0100_0000_0000, diff --git a/src/tables/hyperbolic.rs b/src/tables/hyperbolic.rs index 0a93074..6f838cc 100644 --- a/src/tables/hyperbolic.rs +++ b/src/tables/hyperbolic.rs @@ -1,32 +1,29 @@ -//! Hyperbolic arctangent table. Values: `atanh(2^-i)` as I1F63. -//! -//! Starts at i=1 (atanh(1) undefined). Requires repeats at 4, 13, 40, ... -//! For i ≥ 21, atanh(x) ≈ x to I1F63 precision. +//! Hyperbolic arctangent table for hyperbolic CORDIC. -/// `atanh(2^-i)` as I1F63. Index 0 = atanh(2^-1) = atanh(0.5). +/// `atanh(2^-(k+1))` as I1F63 at index k (atanh(1) is infinite). From 2^-21 +/// on, each entry equals its argument to I1F63 precision. #[rustfmt::skip] pub const ATANH_TABLE: [i64; 64] = [ - 0x464F_A9EA_B40C_2A5E, // atanh(2^-1) = 0.549 - 0x20B1_5DF5_0228_A34E, // atanh(2^-2) - 0x1015_891C_9EAE_F76A, // atanh(2^-3) - 0x0802_AC45_69BA_D66E, // atanh(2^-4) - repeat this - 0x0400_5562_246B_B893, // atanh(2^-5) - 0x0200_0AAB_1115_A393, // atanh(2^-6) - 0x0100_0155_5888_91AD, // atanh(2^-7) - 0x0080_002A_AAC4_4457, // atanh(2^-8) - 0x0040_0005_5556_2222, // atanh(2^-9) - 0x0020_0000_AAAA_B111, // atanh(2^-10) - 0x0010_0000_1555_5589, // atanh(2^-11) - 0x0008_0000_02AA_AAAC, // atanh(2^-12) - 0x0004_0000_0055_5555, // atanh(2^-13) - repeat this - 0x0002_0000_000A_AAAB, // atanh(2^-14) - 0x0001_0000_0001_5555, // atanh(2^-15) - 0x0000_8000_0000_2AAB, // atanh(2^-16) - 0x0000_4000_0000_0555, // atanh(2^-17) - 0x0000_2000_0000_00AB, // atanh(2^-18) - 0x0000_1000_0000_0015, // atanh(2^-19) - 0x0000_0800_0000_0003, // atanh(2^-20) - // For i >= 21: atanh(2^-i) ≈ 2^-i to I1F63 precision + 0x464F_A9EA_B40C_2A5E, + 0x20B1_5DF5_0228_A34E, + 0x1015_891C_9EAE_F76A, + 0x0802_AC45_69BA_D66E, + 0x0400_5562_246B_B893, + 0x0200_0AAB_1115_A393, + 0x0100_0155_5888_91AD, + 0x0080_002A_AAC4_4457, + 0x0040_0005_5556_2222, + 0x0020_0000_AAAA_B111, + 0x0010_0000_1555_5589, + 0x0008_0000_02AA_AAAC, + 0x0004_0000_0055_5555, + 0x0002_0000_000A_AAAB, + 0x0001_0000_0001_5555, + 0x0000_8000_0000_2AAB, + 0x0000_4000_0000_0555, + 0x0000_2000_0000_00AB, + 0x0000_1000_0000_0015, + 0x0000_0800_0000_0003, 0x0000_0400_0000_0000, 0x0000_0200_0000_0000, 0x0000_0100_0000_0000, @@ -46,7 +43,7 @@ pub const ATANH_TABLE: [i64; 64] = [ 0x0000_0000_0400_0000, 0x0000_0000_0200_0000, 0x0000_0000_0100_0000, - 0x0000_0000_0080_0000, // atanh(2^-40) - repeat this + 0x0000_0000_0080_0000, 0x0000_0000_0040_0000, 0x0000_0000_0020_0000, 0x0000_0000_0010_0000, @@ -70,17 +67,16 @@ pub const ATANH_TABLE: [i64; 64] = [ 0x0000_0000_0000_0004, 0x0000_0000_0000_0002, 0x0000_0000_0000_0001, - 0x0000_0000_0000_0001, // rounds to 1 LSB + 0x0000_0000_0000_0001, ]; -/// Returns true if iteration `i` must be repeated for hyperbolic CORDIC convergence. -/// -/// The repeat sequence is 4, 13, 40, 121, 364, ... (each term is 3×previous + 1). +/// Whether hyperbolic CORDIC must repeat iteration `index` to converge: 4, 13, +/// 40, 121, 364, each term three times the last plus one. #[inline] #[must_use] pub const fn needs_repeat(index: u32) -> bool { matches!(index, 4 | 13 | 40 | 121 | 364) } -/// atanh(0.5) ≈ 0.549 (I1F63). Used for argument reduction. +/// atanh(0.5) ≈ 0.549 as I1F63, for argument reduction. pub const ATANH_HALF: i64 = 0x464F_A9EA_B40C_2A5E; diff --git a/src/tables/mod.rs b/src/tables/mod.rs index dc1806f..e0301d9 100644 --- a/src/tables/mod.rs +++ b/src/tables/mod.rs @@ -1,14 +1,5 @@ -//! Precomputed lookup tables for CORDIC algorithms and polynomial evaluation. -//! -//! Tables are stored as `i64` values representing signed I1F63 fixed-point -//! numbers (1 sign bit, 63 fractional bits), which are then converted to -//! the target type at runtime. -//! -//! # Table Contents -//! -//! - [`ATAN_TABLE`]: `atan(2^-i)` values for circular CORDIC mode -//! - [`ATANH_TABLE`]: `atanh(2^-i)` values for hyperbolic CORDIC mode -//! - [`chebyshev`]: Minimax polynomial coefficients for sin/cos evaluation +//! Lookup tables for CORDIC and polynomial evaluation, as `i64` holding signed +//! I1F63 values (1 sign bit, 63 fractional bits) that convert to the target type. pub mod chebyshev; pub mod circular; diff --git a/src/traits.rs b/src/traits.rs index 8a465c3..a31e668 100644 --- a/src/traits.rs +++ b/src/traits.rs @@ -5,24 +5,11 @@ use fixed::traits::{Fixed, FixedSigned}; use fixed::types::extra::{IsLessOrEqual, LeEqU128, True, Unsigned}; use fixed::{FixedI8, FixedI16, FixedI32, FixedI64, FixedI128}; -/// A number type that can be used with CORDIC-based algorithms. +/// The operations and constants the algorithms need from a fixed-point type. /// -/// This trait abstracts over fixed-point number types, providing the -/// operations and constants necessary for CORDIC computations. -/// -/// # Implementors -/// -/// This trait is implemented generically for all signed fixed-point types -/// from the `fixed` crate that have sufficient fractional bits to represent -/// the required constants (π, e, etc.): -/// -/// - [`FixedI8`](fixed::FixedI8) where Fract ≤ 5 (for π to fit) -/// - [`FixedI16`](fixed::FixedI16) where Fract ≤ 13 -/// - [`FixedI32`](fixed::FixedI32) where Fract ≤ 29 -/// - [`FixedI64`](fixed::FixedI64) where Fract ≤ 61 -/// - [`FixedI128`](fixed::FixedI128) where Fract ≤ 125 -/// -/// Common type aliases like `I16F16`, `I32F32`, `I8F24`, `I24F8` all work. +/// Every signed `fixed` type with three or more integer bits implements it, +/// from [`FixedI8`] to [`FixedI128`] with up +/// to 5, 13, 29, 61 or 125 fractional bits; `I16F16` and `I24F8` qualify. pub trait CordicNumber: Copy + PartialEq @@ -68,12 +55,11 @@ pub trait CordicNumber: /// Absolute value. #[must_use] fn abs(self) -> Self; - /// Fractional bits. Determines CORDIC iteration count. + /// Fractional bits, which set the CORDIC iteration count. fn frac_bits() -> u32; /// Total bits. fn total_bits() -> u32; - /// Converts from a raw I1F63 representation (1 sign bit, 63 fractional bits). - /// For constants in (-1, 1). + /// Converts a constant in (-1, 1) from raw I1F63 bits, truncating. fn from_i1f63(bits: i64) -> Self; /// Returns true if negative. fn is_negative(self) -> bool; @@ -90,9 +76,8 @@ pub trait CordicNumber: /// Saturating subtraction. #[must_use] fn saturating_sub(self, rhs: Self) -> Self; - /// Wrapping multiplication, for use only where overflow is provably - /// impossible: then it matches [`saturating_mul`](Self::saturating_mul) - /// bit for bit without the overflow check. + /// Wrapping multiplication, only for where overflow is impossible and it + /// matches [`saturating_mul`](Self::saturating_mul) without the check. #[must_use] fn wrapping_mul(self, rhs: Self) -> Self; /// Wrapping addition. Same contract as [`wrapping_mul`](Self::wrapping_mul). @@ -104,20 +89,18 @@ pub trait CordicNumber: /// Saturating multiplication by an integer. Exact unless it saturates. #[must_use] fn mul_int(self, rhs: i32) -> Self; - /// Division by a positive integer, truncated toward zero. Matches - /// [`div`](Self::div) bit for bit for non-negative `self` but divides the - /// raw representation directly, which is several times cheaper on - /// 128-bit types. A zero divisor saturates like [`div`](Self::div). + /// Division by a positive integer, truncated toward zero. For non-negative + /// `self` it matches [`div`](Self::div) bit for bit, several times cheaper + /// on 128-bit types; a zero divisor saturates the same way. #[must_use] fn div_int(self, divisor: u32) -> Self; /// Division. #[must_use] fn div(self, rhs: Self) -> Self; - /// Integer part of the base-2 logarithm, `⌊log₂(self)⌋`, or `None` if - /// `self ≤ 0`. + /// `⌊log₂(self)⌋`, or `None` if `self ≤ 0`. fn checked_int_log2(self) -> Option; - /// Square root, rounded to the nearest representable value. `self` must - /// be non-negative; the `fixed` types return the root of the magnitude. + /// Square root, rounded to nearest, of a non-negative `self`; the `fixed` + /// types return the root of the magnitude. #[must_use] fn sqrt_round(self) -> Self; /// Convert from numeric type. @@ -126,52 +109,35 @@ pub trait CordicNumber: fn max_value() -> Self; /// Minimum value. fn min_value() -> Self; - /// Round to nearest integer (half away from zero), saturating at the - /// type's bounds. + /// Round to the nearest integer, half away from zero, saturating. #[must_use] fn round(self) -> Self; - /// Convert to i32, rounding toward −∞ and saturating if the value does - /// not fit. + /// Convert to i32, rounding toward −∞ and saturating. #[must_use] fn to_i32(self) -> i32; } -// ============================================================================= -// Square root helpers -// ============================================================================= +// The square root derivations are in docs/sqrt.md. -/// Round-to-nearest square root when `X << (FRAC_NBITS + 2)` fits in -/// `$wide`: the raw result is `⌊√N + ½⌋ = (⌊2√N⌋ + 1) >> 1` with -/// `N = X·2^f`, and `⌊2√N⌋ = isqrt(4N)`. +/// Rounded square root when `X << (FRAC_NBITS + 2)` fits in `$wide`. macro_rules! sqrt_round_widened { ($self:expr, $bits:ty, $wide:ty) => {{ let n4 = <$wide>::from($self.to_bits().unsigned_abs()) << (Self::FRAC_NBITS + 2); - #[allow( - clippy::cast_possible_truncation, - clippy::cast_possible_wrap, - reason = "the root has at most (total_bits + frac_bits) / 2 + 1 bits" - )] - Self::from_bits(((n4.isqrt() + 1) >> 1) as $bits) + Self::from_bits(<$bits>::try_from((n4.isqrt() + 1) >> 1).unwrap_or(<$bits>::MAX)) }}; } -/// Round-to-nearest square root for 128-bit raw representations, where -/// `X·2^f` does not fit in a machine integer. See [`sqrt_round_u128`]. +/// Rounded square root for 128-bit types, through [`sqrt_round_u128`]. macro_rules! sqrt_round_i128 { ($self:expr, $bits:ty, $wide:ty) => {{ let root = sqrt_round_u128($self.to_bits().unsigned_abs(), Self::FRAC_NBITS, |a, b| { - // ⌊a·2^f / b⌋ for b ≥ √N is below 2^127, so this cannot - // overflow; the fallback only keeps the function total. - #[allow(clippy::cast_possible_wrap, reason = "operands below 2^127")] - let quotient = - Fixed::checked_div(Self::from_bits(a as i128), Self::from_bits(b as i128)); - #[allow(clippy::cast_sign_loss, reason = "quotient of positive operands")] - { - quotient.map_or(i128::MAX, Self::to_bits) as u128 - } + let quotient = Fixed::checked_div( + Self::from_bits(a.cast_signed()), + Self::from_bits(b.cast_signed()), + ); + quotient.map_or(i128::MAX, Self::to_bits).cast_unsigned() }); - #[allow(clippy::cast_possible_wrap, reason = "root below 2^127")] - Self::from_bits(root as $bits) + Self::from_bits(root.cast_signed()) }}; } @@ -208,17 +174,8 @@ const fn wide_sub_lo(a: U256, b: U256) -> u128 { a.1.wrapping_sub(b.1) } -/// Round-to-nearest square root of `x / 2^frac` (`x < 2^127`), in raw units: -/// `⌊√N + ½⌋` with `N = x·2^frac`, up to 254 bits. -/// -/// `x` is shifted left by the largest `s ≤ leading_zeros(x)` with -/// `s ≡ frac (mod 2)`, giving `q = isqrt(x·2^s)` with 64 significant bits. -/// If `s ≥ frac`, `q >> ((s − frac)/2)` is exactly `⌊√N⌋`. Otherwise -/// `seed = (q + 1) << m` lies above `√N` by less than `2^(m+1)`, one Newton -/// step from above lands on `⌊√N⌋` or `⌊√N⌋ + 1` (its error is below -/// `2^(m − 64)`, `m ≤ 63`), and a 256-bit remainder settles floor and -/// rounding. `div_scaled(a, b)` must return `⌊a·2^frac / b⌋`; it is only -/// called with `b ≥ √N`. +/// `⌊√N + ½⌋` with `N = x·2^frac` for `x < 2^127`, by one Newton step and a +/// 256-bit remainder. `div_scaled(a, b)` must return `⌊a·2^frac / b⌋`. fn sqrt_round_u128(x: u128, frac: u32, div_scaled: impl FnOnce(u128, u128) -> u128) -> u128 { if x == 0 { return 0; @@ -230,54 +187,42 @@ fn sqrt_round_u128(x: u128, frac: u32, div_scaled: impl FnOnce(u128, u128) -> u1 if shift >= frac { let excess = (shift - frac) / 2; if excess >= 1 { - // ⌊2√N⌋ = root_hi >> (excess − 1); round half up. return ((root_hi >> (excess - 1)) + 1) >> 1; } - // root_hi = ⌊√N⌋, and N fits in 128 bits since shift = frac ≤ lz. + // N fits in 128 bits since shift = frac ≤ lz. let root_sq = root_hi * root_hi; let rem = (x << frac) - root_sq; return if rem > root_hi { root_hi + 1 } else { root_hi }; } let missing = (frac - shift) / 2; - #[allow(clippy::cast_sign_loss, reason = "i128::MAX is positive")] - let seed = ((root_hi + 1) << missing).min(i128::MAX as u128); + let seed = ((root_hi + 1) << missing).min(i128::MAX.cast_unsigned()); let newton = u128::midpoint(seed, div_scaled(x, seed)); let n = wide_shl(x, frac); let newton_sq = wide_mul(newton, newton); let (root, rem) = if wide_gt(newton_sq, n) { - // newton = ⌊√N⌋ + 1: N − (newton − 1)² = (2·newton − 1) − (newton² − N). + // N − (newton − 1)² = (2·newton − 1) − (newton² − N). (newton - 1, (2 * newton - 1) - wide_sub_lo(newton_sq, n)) } else { (newton, wide_sub_lo(n, newton_sq)) }; - // rem = N − root² ∈ [0, 2·root]; round up exactly when rem > root. if rem > root { root + 1 } else { root } } -// ============================================================================= -// Generic implementations using macros -// ============================================================================= - -/// Macro to implement `CordicNumber` for `FixedI*` types generically. -/// -/// The bounds ensure: -/// - `Fract` fits within the type (e.g., ≤ 8 for `FixedI8`) -/// - `Fract` allows π to be represented (needs ~2 integer bits) -/// - `Fract` allows π/2 to be represented (needs ~1 integer bit) -/// - `Fract` allows π/4 and ln(2) to be represented (needs ~1 integer bit) +/// Implements `CordicNumber` for a `FixedI*` family, bounding `Fract` so the +/// type holds its constants: π needs two integer bits, π/2, π/4 and ln 2 one. macro_rules! impl_cordic_generic { ( $fixed_type:ident, $bits_type:ty, - $wide_type:ty, // Unsigned type holding bits << (frac + 2) for sqrt - $sqrt_impl:ident, // sqrt_round_widened or sqrt_round_i128 + $wide_type:ty, + $sqrt_impl:ident, $total_bits:expr, - $max_frac:ty, // Maximum fractional bits for the type - $pi_frac:ty, // Max frac bits where PI fits (total - 2) - $frac_pi_2:ty, // Max frac bits where FRAC_PI_2 fits (total - 1) - $frac_pi_4:ty // Max frac bits where FRAC_PI_4 and LN_2 fit + $max_frac:ty, + $pi_frac:ty, + $frac_pi_2:ty, + $frac_pi_4:ty ) => { impl CordicNumber for $fixed_type where @@ -339,35 +284,15 @@ macro_rules! impl_cordic_generic { } #[inline] - // Casts are safe: frac_bits ≤ 128, shift amounts bounded by type size - #[allow( - clippy::cast_possible_wrap, - clippy::cast_lossless, - reason = "frac_bits bounded by type size" - )] fn from_i1f63(bits: i64) -> Self { - // Convert from I1F63 representation to our type. - // I1F63 has 63 fractional bits. - // FRAC_NBITS is at most 128, which fits in i32. - let our_frac = Self::FRAC_NBITS as i32; - let shift = 63 - our_frac; - + let shift = 63 - Self::FRAC_NBITS.cast_signed(); + // The trait bounds leave integer bits to spare, so a value in + // (-1, 1) always fits and the fallback is dead. if shift >= 0 { - // We have fewer frac bits than I1F63, shift right - #[allow( - clippy::cast_possible_truncation, - reason = "intentional truncation to target type" - )] - Self::from_bits((bits >> shift) as $bits_type) + Self::from_bits(<$bits_type>::try_from(bits >> shift).unwrap_or(0)) } else { - // We have more frac bits than I1F63, shift left - // Must cast first to avoid losing sign bit - #[allow( - clippy::cast_possible_truncation, - reason = "intentional truncation to target type" - )] - let wide = bits as $bits_type; - Self::from_bits(wide << (-shift)) + // Only 128-bit types get here; widen before the left shift. + Self::from_bits(<$bits_type>::try_from(bits).unwrap_or(0) << (-shift)) } } @@ -410,8 +335,7 @@ macro_rules! impl_cordic_generic { fn mul_int(self, rhs: i32) -> Self { match <$bits_type>::try_from(rhs) { Ok(k) => Fixed::saturating_mul_int(self, k), - // |rhs| exceeds the raw type's range (8- and 16-bit types - // only), so the exact product saturates unless self is 0. + // |rhs| is past the raw range, so only zero escapes saturation. Err(_) => { if self == Self::ZERO { Self::ZERO @@ -429,7 +353,6 @@ macro_rules! impl_cordic_generic { match <$bits_type>::try_from(divisor) { Ok(k) => match Fixed::checked_div_int(self, k) { Some(v) => v, - // Division by zero: saturate based on sign. None => { if self.is_negative() { Self::MIN @@ -438,8 +361,7 @@ macro_rules! impl_cordic_generic { } } }, - // The divisor exceeds the raw type's range, so it exceeds - // |self| in raw units and the truncated quotient is zero. + // A divisor past the raw range exceeds |self|: the quotient is 0. Err(_) => Self::ZERO, } } @@ -448,7 +370,6 @@ macro_rules! impl_cordic_generic { fn div(self, rhs: Self) -> Self { match Fixed::checked_div(self, rhs) { Some(v) => v, - // Division by zero or overflow: saturate based on sign agreement. None => { if self.is_negative() != rhs.is_negative() { Self::MIN @@ -497,23 +418,14 @@ macro_rules! impl_cordic_generic { }; } -// Import the specific U* types we need for bounds use fixed::types::extra::{ U5, U6, U7, U8, U13, U14, U15, U16, U29, U30, U31, U32, U61, U62, U63, U64, U125, U126, U127, U128, }; -// FixedI8: 8 total bits -// - Max Fract: U8 (8 fractional bits = I0F8) -// - For PI (~3.14), need 2 integer bits, so Fract ≤ 6 (I2F6) -// - For FRAC_PI_2, FRAC_PI_4, LN_2, need 1 integer bit, so Fract ≤ 7 (I1F7) -// Being conservative: require Fract ≤ 5 so we have headroom +// Each type's π bound is one bit stricter than π needs, for headroom. impl_cordic_generic!(FixedI8, i8, u16, sqrt_round_widened, 8, U8, U5, U6, U7); -// FixedI16: 16 total bits -// - For PI, need Fract ≤ 14 (I2F14) -// - For FRAC_PI_2, FRAC_PI_4, LN_2, need Fract ≤ 15 (I1F15) -// - Conservative: Fract ≤ 13 impl_cordic_generic!( FixedI16, i16, @@ -526,10 +438,6 @@ impl_cordic_generic!( U15 ); -// FixedI32: 32 total bits -// - For PI, need Fract ≤ 30 -// - For FRAC_PI_2, FRAC_PI_4, LN_2, need Fract ≤ 31 -// - Conservative: Fract ≤ 29 impl_cordic_generic!( FixedI32, i32, @@ -542,10 +450,6 @@ impl_cordic_generic!( U31 ); -// FixedI64: 64 total bits -// - For PI, need Fract ≤ 62 -// - For FRAC_PI_2, FRAC_PI_4, LN_2, need Fract ≤ 63 -// - Conservative: Fract ≤ 61 impl_cordic_generic!( FixedI64, i64, @@ -558,11 +462,6 @@ impl_cordic_generic!( U63 ); -// FixedI128: 128 total bits -// - For PI, need Fract ≤ 126 -// - For FRAC_PI_2, FRAC_PI_4, LN_2, need Fract ≤ 127 -// - Conservative: Fract ≤ 125 -// No wider machine integer exists, so sqrt takes the Newton path. impl_cordic_generic!( FixedI128, i128, From 142a8e4e7511858a2d9916ccd9f94e5b07d1bef9 Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 17:40:56 +0000 Subject: [PATCH 3/9] Tests return Result instead of unwrapping, satisfy the lints without suppressions, and keep only comments that explain why --- tests/main.rs | 2 +- tests/unit/error.rs | 2 - tests/unit/kernel/cordic.rs | 7 +- tests/unit/kernel/mod.rs | 2 - tests/unit/mod.rs | 2 - tests/unit/ops/algebraic.rs | 160 ++++++++++----------- tests/unit/ops/circular.rs | 163 +++++++-------------- tests/unit/ops/exponential.rs | 243 +++++++++++--------------------- tests/unit/ops/hyperbolic.rs | 131 +++++------------ tests/unit/ops/mod.rs | 2 - tests/unit/smoke.rs | 113 ++++++++------- tests/unit/support.rs | 13 +- tests/unit/tables/chebyshev.rs | 34 +---- tests/unit/tables/circular.rs | 47 +++--- tests/unit/tables/hyperbolic.rs | 57 +++----- tests/unit/tables/mod.rs | 2 - tests/unit/traits.rs | 24 +--- tests/unit/verification.rs | 176 ++++++----------------- 18 files changed, 393 insertions(+), 787 deletions(-) diff --git a/tests/main.rs b/tests/main.rs index 6e57913..624e7f5 100644 --- a/tests/main.rs +++ b/tests/main.rs @@ -1,3 +1,3 @@ -//! Test harness for unit tests. +//! Unit tests for `fixed_analytics`. mod unit; diff --git a/tests/unit/error.rs b/tests/unit/error.rs index 8010ed4..0c01216 100644 --- a/tests/unit/error.rs +++ b/tests/unit/error.rs @@ -1,5 +1,3 @@ -//! Tests for error types - #[cfg(test)] mod tests { use fixed_analytics::Error; diff --git a/tests/unit/kernel/cordic.rs b/tests/unit/kernel/cordic.rs index bb0d56d..3140e89 100644 --- a/tests/unit/kernel/cordic.rs +++ b/tests/unit/kernel/cordic.rs @@ -1,5 +1,3 @@ -//! Tests for CORDIC core algorithms - #[cfg(test)] mod tests { use fixed::types::{I16F16, I64F64}; @@ -7,7 +5,6 @@ mod tests { #[test] fn circular_vectoring_atan_one() { - // vectoring mode with x=1, y=1 should give z ≈ π/4 let (_, _, z) = circular_vectoring(I16F16::ONE, I16F16::ONE, I16F16::ZERO); let z_f32: f32 = z.to_num(); let expected = core::f32::consts::FRAC_PI_4; @@ -16,8 +13,7 @@ mod tests { #[test] fn circular_vectoring_high_precision_type() { - // I64F64 has 64 fractional bits, more than the I1F63 tables: - // exercises the exact (no rounding needed) conversion path. + // I64F64 has more fractional bits than the I1F63 tables: the exact conversion path. let (_, _, z) = circular_vectoring(I64F64::ONE, I64F64::ONE, I64F64::ZERO); let z_f64: f64 = z.to_num(); let expected = core::f64::consts::FRAC_PI_4; @@ -28,7 +24,6 @@ mod tests { } } -/// Sweep both early-terminating kernels against `f64` at three widths. #[cfg(test)] mod early_termination { use crate::unit::support::Lcg; diff --git a/tests/unit/kernel/mod.rs b/tests/unit/kernel/mod.rs index c720457..dd387ca 100644 --- a/tests/unit/kernel/mod.rs +++ b/tests/unit/kernel/mod.rs @@ -1,3 +1 @@ -//! Tests for CORDIC kernel operations - mod cordic; diff --git a/tests/unit/mod.rs b/tests/unit/mod.rs index f360249..723743d 100644 --- a/tests/unit/mod.rs +++ b/tests/unit/mod.rs @@ -1,5 +1,3 @@ -//! Unit tests mirroring the src directory structure - mod error; mod kernel; mod ops; diff --git a/tests/unit/ops/algebraic.rs b/tests/unit/ops/algebraic.rs index fa032b5..d99a47f 100644 --- a/tests/unit/ops/algebraic.rs +++ b/tests/unit/ops/algebraic.rs @@ -1,8 +1,6 @@ -//! Tests for algebraic functions (sqrt) - #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod tests { + use crate::unit::support::TestResult; use fixed::types::I16F16; use fixed_analytics::sqrt; @@ -13,27 +11,29 @@ mod tests { } #[test] - fn sqrt_perfect_squares() { - assert!(approx_eq(sqrt(I16F16::from_num(0.0)).unwrap(), 0.0)); - assert!(approx_eq(sqrt(I16F16::from_num(1.0)).unwrap(), 1.0)); - assert!(approx_eq(sqrt(I16F16::from_num(4.0)).unwrap(), 2.0)); - assert!(approx_eq(sqrt(I16F16::from_num(9.0)).unwrap(), 3.0)); - assert!(approx_eq(sqrt(I16F16::from_num(16.0)).unwrap(), 4.0)); - assert!(approx_eq(sqrt(I16F16::from_num(25.0)).unwrap(), 5.0)); + fn sqrt_perfect_squares() -> TestResult { + assert!(approx_eq(sqrt(I16F16::from_num(0.0))?, 0.0)); + assert!(approx_eq(sqrt(I16F16::from_num(1.0))?, 1.0)); + assert!(approx_eq(sqrt(I16F16::from_num(4.0))?, 2.0)); + assert!(approx_eq(sqrt(I16F16::from_num(9.0))?, 3.0)); + assert!(approx_eq(sqrt(I16F16::from_num(16.0))?, 4.0)); + assert!(approx_eq(sqrt(I16F16::from_num(25.0))?, 5.0)); + Ok(()) } #[test] - fn sqrt_common_values() { + fn sqrt_common_values() -> TestResult { assert!(approx_eq( - sqrt(I16F16::from_num(2.0)).unwrap(), + sqrt(I16F16::from_num(2.0))?, core::f32::consts::SQRT_2 )); - assert!(approx_eq(sqrt(I16F16::from_num(3.0)).unwrap(), 1.7321)); + assert!(approx_eq(sqrt(I16F16::from_num(3.0))?, 1.7321)); assert!(approx_eq( - sqrt(I16F16::from_num(0.5)).unwrap(), + sqrt(I16F16::from_num(0.5))?, core::f32::consts::FRAC_1_SQRT_2 )); - assert!(approx_eq(sqrt(I16F16::from_num(0.25)).unwrap(), 0.5)); + assert!(approx_eq(sqrt(I16F16::from_num(0.25))?, 0.5)); + Ok(()) } #[test] @@ -43,10 +43,10 @@ mod tests { } #[test] - fn sqrt_squared_gives_original() { + fn sqrt_squared_gives_original() -> TestResult { for i in 1..20 { let x = I16F16::from_num(i) * I16F16::from_num(0.5); - let root = sqrt(x).unwrap(); + let root = sqrt(x)?; let squared: f32 = (root * root).to_num(); let original: f32 = x.to_num(); assert!( @@ -54,25 +54,18 @@ mod tests { "sqrt({original})² = {squared}, expected {original}" ); } + Ok(()) } } -/// Exactness: every result is the representable value nearest the true -/// root, checked against the `fixed` crate's rounded-down square root plus -/// a 256-bit remainder test for the rounding direction. +/// Each root must be the representable value nearest the true root. #[cfg(test)] -#[allow( - clippy::unwrap_used, - clippy::cast_sign_loss, - clippy::cast_possible_wrap, - clippy::cast_possible_truncation, - reason = "test code" -)] mod exact { use crate::unit::support::Lcg; + use crate::unit::support::TestResult; use fixed::traits::Fixed; use fixed::types::{ - I3F125, I4F4, I4F60, I8F8, I8F24, I16F16, I24F8, I32F32, I48F16, I64F64, I128F0, + I3F125, I4F4, I4F60, I8F8, I8F24, I16F16, I24F8, I32F32, I48F16, I64F64, I128F0, U128F0, }; use fixed_analytics::sqrt; @@ -96,28 +89,24 @@ mod exact { } } - /// `a - b` for `a ≥ b`. fn wide_sub(a: U256, b: U256) -> U256 { let (lo, borrow) = a.1.overflowing_sub(b.1); (a.0 - b.0 - u128::from(borrow), lo) } - /// Raw bits of the correctly rounded root of `x`. fn expected_bits(x: T) -> i128 where T::Bits: Into, { let floor: i128 = x.sqrt().to_bits().into(); - let n = wide_shl(x.to_bits().into() as u128, T::FRAC_NBITS); - let f = floor as u128; + let xb: i128 = x.to_bits().into(); + let n = wide_shl(xb.cast_unsigned(), T::FRAC_NBITS); + let f = floor.cast_unsigned(); let rem = wide_sub(n, wide_mul(f, f)); assert_eq!(rem.0, 0, "remainder exceeds 2·root + 1"); if rem.1 > f { floor + 1 } else { floor } } - /// Deterministic sample set spanning the whole non-negative range of `T`: - /// special values, powers of two and neighbours, near-squares, and - /// uniformly random raw bits at several magnitudes. fn samples(random: usize) -> Vec where T::Bits: TryFrom + Into, @@ -130,8 +119,7 @@ mod exact { raw.push((1i128 << e) + d); } } - // X = k² ± j lands √(X·2^f) near an integer, exercising both - // rounding directions and the 128-bit Newton correction step. + // X = k² ± j puts the root near an integer: both rounding directions and the Newton fix. for p in 1..64u32 { let k = 1i128 << p; if let Some(k2) = k.checked_mul(k) { @@ -142,7 +130,7 @@ mod exact { } let mut rng = Lcg(0xABCD_EF01 ^ u64::from(bits)); for _ in 0..random { - let r = (u128::from(rng.next_u64()) | (u128::from(rng.next_u64()) << 64)) as i128; + let r = (u128::from(rng.next_u64()) | (u128::from(rng.next_u64()) << 64)).cast_signed(); let full = (r >> (128 - bits)).abs().min(max); raw.extend([ full, @@ -162,111 +150,107 @@ mod exact { .collect() } - fn check(random: usize) + fn check(random: usize) -> TestResult where T: Fixed + fixed_analytics::CordicNumber + core::fmt::Display, T::Bits: TryFrom + Into, { for x in samples::(random) { - let got: i128 = sqrt(x).unwrap().to_bits().into(); + let got: i128 = sqrt(x)?.to_bits().into(); let want = expected_bits(x); assert_eq!(got, want, "sqrt({x}) raw bits: got {got}, want {want}"); } + Ok(()) } #[test] - fn i4f4_and_i8f8_are_correctly_rounded() { - check::(500); - check::(2000); + fn i4f4_and_i8f8_are_correctly_rounded() -> TestResult { + check::(500)?; + check::(2000)?; + Ok(()) } #[test] - fn i16f16_i8f24_i24f8_are_correctly_rounded() { - check::(5000); - check::(2000); - check::(2000); + fn i16f16_i8f24_i24f8_are_correctly_rounded() -> TestResult { + check::(5000)?; + check::(2000)?; + check::(2000)?; + Ok(()) } #[test] - fn i32f32_i4f60_i48f16_are_correctly_rounded() { - check::(5000); - check::(2000); - check::(2000); + fn i32f32_i4f60_i48f16_are_correctly_rounded() -> TestResult { + check::(5000)?; + check::(2000)?; + check::(2000)?; + Ok(()) } #[test] - fn i64f64_i3f125_i128f0_are_correctly_rounded() { + fn i64f64_i3f125_i128f0_are_correctly_rounded() -> TestResult { // 128-bit types take the Newton path above 2^(128 - 2·frac). - check::(5000); - check::(3000); - check::(1000); + check::(5000)?; + check::(3000)?; + check::(1000)?; + Ok(()) } /// Below 0.25 on I64F64 the exact root fits in `u128` arithmetic. #[test] - fn i64f64_sub_unit_matches_integer_root() { + fn i64f64_sub_unit_matches_integer_root() -> TestResult { let mut rng = Lcg(0x5EED); let mut raw: Vec = vec![1, 2, 3, u64::MAX.into()]; - // Log-spaced from one ulp (2^-64 ≈ 5e-20) up to 0.25. for i in 0..4000 { let exponent = 62.0 * f64::from(i) / 4000.0; - raw.push((exponent).exp2().max(1.0) as u128); + raw.push(U128F0::saturating_from_num(exponent.exp2().max(1.0).trunc()).to_bits()); raw.push(rng.next_u64().into()); raw.push(u128::from(rng.next_u64()) >> (rng.next_u64() % 62)); } for x_raw in raw { let x_raw = x_raw.min((1u128 << 62) - 1); - let x = I64F64::from_bits(x_raw as i128); - // round(√(X·2^64)) = (⌊2·√(X·2^64)⌋ + 1) >> 1 = (isqrt(X·2^66) + 1) >> 1 - let want = (((x_raw << 66).isqrt() + 1) >> 1) as i128; - let got = sqrt(x).unwrap().to_bits(); + let x = I64F64::from_bits(x_raw.cast_signed()); + // round(√(X·2^64)) = (isqrt(X·2^66) + 1) >> 1 + let want = (((x_raw << 66).isqrt() + 1) >> 1).cast_signed(); + let got = sqrt(x)?.to_bits(); assert_eq!(got, want, "sqrt({x}) raw bits: got {got}, want {want}"); } + Ok(()) } - /// Inputs just below a perfect square make the Newton step land one - /// above the floor. + /// Just below a perfect square, the Newton step lands one above the floor. #[test] - fn i64f64_newton_correction_cases() { + fn i64f64_newton_correction_cases() -> TestResult { for p in 33..63u32 { for j in 1..=3i128 { let x = I64F64::from_bits((1i128 << (2 * p - 64)) - j); - let got: i128 = sqrt(x).unwrap().to_bits(); + let got: i128 = sqrt(x)?.to_bits(); let want = expected_bits(x); assert_eq!(got, want, "sqrt({x}) raw bits: got {got}, want {want}"); } } + Ok(()) } #[test] - fn special_values() { - assert_eq!(sqrt(I64F64::ZERO).unwrap(), I64F64::ZERO); - assert_eq!(sqrt(I64F64::ONE).unwrap(), I64F64::ONE); - assert_eq!(sqrt(I64F64::from_num(4)).unwrap(), I64F64::from_num(2)); - assert_eq!(sqrt(I64F64::from_num(0.25)).unwrap(), I64F64::from_num(0.5)); - assert_eq!( - sqrt(I64F64::from_bits(1)).unwrap(), - I64F64::from_bits(1 << 32) - ); - assert_eq!( - sqrt(I64F64::MAX).unwrap().to_bits(), - expected_bits(I64F64::MAX) - ); - assert_eq!( - sqrt(I3F125::MAX).unwrap().to_bits(), - expected_bits(I3F125::MAX) - ); + fn special_values() -> TestResult { + assert_eq!(sqrt(I64F64::ZERO)?, I64F64::ZERO); + assert_eq!(sqrt(I64F64::ONE)?, I64F64::ONE); + assert_eq!(sqrt(I64F64::from_num(4))?, I64F64::from_num(2)); + assert_eq!(sqrt(I64F64::from_num(0.25))?, I64F64::from_num(0.5)); + assert_eq!(sqrt(I64F64::from_bits(1))?, I64F64::from_bits(1 << 32)); + assert_eq!(sqrt(I64F64::MAX)?.to_bits(), expected_bits(I64F64::MAX)); + assert_eq!(sqrt(I3F125::MAX)?.to_bits(), expected_bits(I3F125::MAX)); assert_eq!( - i128::from(sqrt(I16F16::MAX).unwrap().to_bits()), + i128::from(sqrt(I16F16::MAX)?.to_bits()), expected_bits(I16F16::MAX) ); + Ok(()) } #[test] - fn asin_acos_still_saturate_near_one() { + fn asin_acos_still_saturate_near_one() -> TestResult { let x = I64F64::ONE - I64F64::from_bits(1); - assert!( - (fixed_analytics::asin(x).unwrap() - I64F64::FRAC_PI_2).abs() < I64F64::from_num(1e-9) - ); + assert!((fixed_analytics::asin(x)? - I64F64::FRAC_PI_2).abs() < I64F64::from_num(1e-9)); + Ok(()) } } diff --git a/tests/unit/ops/circular.rs b/tests/unit/ops/circular.rs index 4338348..27a240e 100644 --- a/tests/unit/ops/circular.rs +++ b/tests/unit/ops/circular.rs @@ -1,13 +1,6 @@ -//! Tests for circular trigonometric functions - #[cfg(test)] -#[allow( - clippy::unwrap_used, - clippy::cast_precision_loss, - clippy::cast_possible_truncation, - reason = "test code uses unwrap and f32/f64 casts for conciseness" -)] mod tests { + use crate::unit::support::TestResult; use fixed::types::{I16F16, I32F32}; use fixed_analytics::{acos, asin, atan, atan2, cos, sin, sin_cos, tan}; @@ -17,7 +10,6 @@ mod tests { (a.to_num::() - b).abs() < TOLERANCE } - // Tests sin(0) = 0, sin(π/2) = 1, sin(-π/2) = -1, sin(π) = 0 #[test] fn sin_special_values() { assert!(approx_eq(sin(I16F16::ZERO), 0.0)); @@ -39,14 +31,13 @@ mod tests { assert!(approx_eq(tan(I16F16::FRAC_PI_4), 1.0)); } - // Tests Pythagorean identity: sin²(x) + cos²(x) = 1 #[test] fn sin_cos_pythagorean_identity() { - for i in -20..=20 { + for i in -20i16..=20 { let angle = I16F16::from_num(i) * I16F16::from_num(0.1); let (s, c) = sin_cos(angle); let sum_sq: f32 = (s * s + c * c).to_num(); - let angle_val = i as f32 * 0.1; + let angle_val = f32::from(i) * 0.1; assert!( (sum_sq - 1.0).abs() < 0.02, "sin²({angle_val}) + cos²({angle_val}) = {sum_sq}, expected ~1.0" @@ -61,33 +52,26 @@ mod tests { assert!(approx_eq(atan(-I16F16::ONE), -core::f32::consts::FRAC_PI_4)); } - // Tests atan2 quadrant correctness #[test] fn atan2_quadrants() { let one = I16F16::ONE; let neg_one = -I16F16::ONE; - // Q1: both positive let q1 = atan2(one, one); assert!(q1 > I16F16::ZERO && q1 < I16F16::FRAC_PI_2); - // Q2: y positive, x negative let q2 = atan2(one, neg_one); assert!(q2 > I16F16::FRAC_PI_2); - // Q3: both negative let q3 = atan2(neg_one, neg_one); assert!(q3 < -I16F16::FRAC_PI_2); - // Q4: y negative, x positive let q4 = atan2(neg_one, one); assert!(q4 < I16F16::ZERO && q4 > -I16F16::FRAC_PI_2); } #[test] fn atan2_precise_quadrant_values() { - // Test precise values for each quadrant - // atan2(1, 1) should be π/4 ≈ 0.785 let q1: f32 = atan2(I16F16::ONE, I16F16::ONE).to_num(); let expected_q1 = core::f32::consts::FRAC_PI_4; assert!( @@ -95,7 +79,6 @@ mod tests { "Q1: expected {expected_q1}, got {q1}" ); - // atan2(1, -1) should be 3π/4 ≈ 2.356 let q2: f32 = atan2(I16F16::ONE, -I16F16::ONE).to_num(); let expected_q2 = core::f32::consts::PI - core::f32::consts::FRAC_PI_4; assert!( @@ -103,7 +86,6 @@ mod tests { "Q2: expected {expected_q2}, got {q2}" ); - // atan2(-1, -1) should be -3π/4 ≈ -2.356 let q3: f32 = atan2(-I16F16::ONE, -I16F16::ONE).to_num(); let expected_q3 = -core::f32::consts::PI + core::f32::consts::FRAC_PI_4; assert!( @@ -111,7 +93,6 @@ mod tests { "Q3: expected {expected_q3}, got {q3}" ); - // atan2(-1, 1) should be -π/4 ≈ -0.785 let q4: f32 = atan2(-I16F16::ONE, I16F16::ONE).to_num(); let expected_q4 = -core::f32::consts::FRAC_PI_4; assert!( @@ -122,33 +103,27 @@ mod tests { #[test] fn atan2_axis_values() { - // Test values along the axes - // atan2(0, 1) = 0 assert!(approx_eq(atan2(I16F16::ZERO, I16F16::ONE), 0.0)); - // atan2(0, -1) = π assert!(approx_eq( atan2(I16F16::ZERO, -I16F16::ONE), core::f32::consts::PI )); - // atan2(1, 0) = π/2 assert!(approx_eq( atan2(I16F16::ONE, I16F16::ZERO), core::f32::consts::FRAC_PI_2 )); - // atan2(-1, 0) = -π/2 assert!(approx_eq( atan2(-I16F16::ONE, I16F16::ZERO), -core::f32::consts::FRAC_PI_2 )); - // atan2(0, 0) = 0 (undefined but we return 0) + // atan2(0, 0) is undefined; the crate returns 0. assert!(approx_eq(atan2(I16F16::ZERO, I16F16::ZERO), 0.0)); } - // Tests asin domain validation rejects |x| > 1 #[test] fn asin_domain_check() { assert!(asin(I16F16::from_num(1.5)).is_err()); @@ -157,17 +132,17 @@ mod tests { } #[test] - fn acos_special_values() { - let result: f32 = acos(I16F16::ONE).unwrap().to_num(); + fn acos_special_values() -> TestResult { + let result: f32 = acos(I16F16::ONE)?.to_num(); assert!(result.abs() < 0.01); + Ok(()) } #[test] - fn sin_asin_roundtrip() { - // sin(asin(x)) ≈ x for x in [-1, 1] + fn sin_asin_roundtrip() -> TestResult { for i in -9..=9 { let x = I16F16::from_num(i) * I16F16::from_num(0.1); - let result = sin(asin(x).unwrap()); + let result = sin(asin(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -175,14 +150,14 @@ mod tests { "sin(asin({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] - fn cos_acos_roundtrip() { - // cos(acos(x)) ≈ x for x in [-1, 1] + fn cos_acos_roundtrip() -> TestResult { for i in -9..=9 { let x = I16F16::from_num(i) * I16F16::from_num(0.1); - let result = cos(acos(x).unwrap()); + let result = cos(acos(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -190,11 +165,11 @@ mod tests { "cos(acos({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] fn tan_atan_roundtrip() { - // tan(atan(x)) ≈ x for various x for i in -10..=10 { let x = I16F16::from_num(i) * I16F16::from_num(0.5); let result = tan(atan(x)); @@ -208,63 +183,59 @@ mod tests { } #[test] - fn asin_boundary_values() { - // asin at exactly ±1 should return ±π/2 - let result_pos: f32 = asin(I16F16::ONE).unwrap().to_num(); + fn asin_boundary_values() -> TestResult { + let result_pos: f32 = asin(I16F16::ONE)?.to_num(); let expected_pos = core::f32::consts::FRAC_PI_2; assert!( (result_pos - expected_pos).abs() < 0.01, "asin(1) = {result_pos}, expected {expected_pos}" ); - let result_neg: f32 = asin(-I16F16::ONE).unwrap().to_num(); + let result_neg: f32 = asin(-I16F16::ONE)?.to_num(); let expected_neg = -core::f32::consts::FRAC_PI_2; assert!( (result_neg - expected_neg).abs() < 0.01, "asin(-1) = {result_neg}, expected {expected_neg}" ); - // Near boundary (0.999) let near_one = I16F16::from_num(0.999); assert!(asin(near_one).is_ok()); - let result: f32 = asin(near_one).unwrap().to_num(); + let result: f32 = asin(near_one)?.to_num(); assert!( result > 1.4 && result < 1.6, "asin(0.999) should be near π/2" ); + Ok(()) } #[test] - fn acos_boundary_values() { - // acos at exactly ±1 should return 0 or π - let result_pos: f32 = acos(I16F16::ONE).unwrap().to_num(); + fn acos_boundary_values() -> TestResult { + let result_pos: f32 = acos(I16F16::ONE)?.to_num(); assert!( result_pos.abs() < 0.01, "acos(1) = {result_pos}, expected 0" ); - let result_neg: f32 = acos(-I16F16::ONE).unwrap().to_num(); + let result_neg: f32 = acos(-I16F16::ONE)?.to_num(); let expected_neg = core::f32::consts::PI; assert!( (result_neg - expected_neg).abs() < 0.01, "acos(-1) = {result_neg}, expected {expected_neg}" ); + Ok(()) } #[test] fn sin_cos_large_angles() { - // Test angle reduction with large angles - let large_angle = I16F16::from_num(100.0); // Well beyond 2π + let large_angle = I16F16::from_num(100.0); let (s, c) = sin_cos(large_angle); - // Should still satisfy sin²+cos² = 1 let sum_sq: f32 = (s * s + c * c).to_num(); assert!( (sum_sq - 1.0).abs() < 0.05, "sin²(100) + cos²(100) = {sum_sq}, expected ~1.0" ); - // Test negative large angle let neg_large_angle = I16F16::from_num(-100.0); let (s2, c2) = sin_cos(neg_large_angle); let sum_sq2: f32 = (s2 * s2 + c2 * c2).to_num(); @@ -276,7 +247,6 @@ mod tests { #[test] fn atan_large_values() { - // atan of large values should approach ±π/2 let large = I16F16::from_num(1000.0); let result: f32 = atan(large).to_num(); let expected = core::f32::consts::FRAC_PI_2; @@ -295,83 +265,73 @@ mod tests { } #[test] - fn asin_near_negative_one() { - // Test asin for values very close to -1, exercising the boundary path - // where sqrt_term is very small and x is negative + fn asin_near_negative_one() -> TestResult { let near_neg_one = I16F16::from_num(-0.9999); let result = asin(near_neg_one); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); - // asin(-0.9999) should be very close to -π/2 + let val: f32 = result?.to_num(); let expected = -core::f32::consts::FRAC_PI_2; assert!( (val - expected).abs() < 0.1, "asin(-0.9999) = {val}, expected ~{expected}" ); - // Test exactly at -1 let neg_one = -I16F16::ONE; let result_exact = asin(neg_one); assert!(result_exact.is_ok()); - let val_exact: f32 = result_exact.unwrap().to_num(); + let val_exact: f32 = result_exact?.to_num(); assert!( (val_exact - expected).abs() < 0.01, "asin(-1) = {val_exact}, expected ~{expected}" ); + Ok(()) } #[test] - fn asin_boundary_very_close() { - // Test with higher precision type to hit the small sqrt_term boundary path - // The threshold is 2^(-15) ≈ 0.0000305, so we need sqrt(1-x²) < 0.0000305 - // This requires 1-x² < 2^(-30), meaning x must be within 2^(-31) of ±1 - // Using I32F32 for higher precision - - // Create a value extremely close to -1 but not exactly -1 - // -1 + 2^(-32) (smallest representable distance) + fn asin_boundary_very_close() -> TestResult { + // The boundary path needs √(1 - x²) < 2^-15, so x within 2^-31 of -1. let near_neg_one = -I32F32::ONE + I32F32::from_bits(1); let result = asin(near_neg_one); assert!(result.is_ok()); - let val: f64 = result.unwrap().to_num(); + let val: f64 = result?.to_num(); let expected = -core::f64::consts::FRAC_PI_2; assert!( (val - expected).abs() < 0.01, "asin(~-1) = {val}, expected ~{expected}" ); - // Also test positive boundary let near_pos_one = I32F32::ONE - I32F32::from_bits(1); let result_pos = asin(near_pos_one); assert!(result_pos.is_ok()); - let val_pos: f64 = result_pos.unwrap().to_num(); + let val_pos: f64 = result_pos?.to_num(); let expected_pos = core::f64::consts::FRAC_PI_2; assert!( (val_pos - expected_pos).abs() < 0.01, "asin(~1) = {val_pos}, expected ~{expected_pos}" ); + Ok(()) } #[test] fn asin_domain_error_message() { - // Test that domain errors are properly returned for out-of-range values let too_large = I16F16::from_num(1.5); - let err = asin(too_large).unwrap_err(); - assert!(matches!(err, fixed_analytics::Error::DomainError { .. })); + assert!(matches!( + asin(too_large), + Err(fixed_analytics::Error::DomainError { .. }) + )); let too_small = I16F16::from_num(-1.5); - let err2 = asin(too_small).unwrap_err(); - assert!(matches!(err2, fixed_analytics::Error::DomainError { .. })); + assert!(matches!( + asin(too_small), + Err(fixed_analytics::Error::DomainError { .. }) + )); } #[test] fn sin_cos_very_large_positive_angle() { - // Test with very large positive angle to exercise the reduction loop - // The loop iterates while reduced > pi, subtracting 2π each time - // 200.0 radians is about 32 times 2π, so the loop will iterate many times let very_large = I16F16::from_num(200.0); let (s, c) = sin_cos(very_large); - // Should still satisfy sin²+cos² = 1 let sum_sq: f32 = (s * s + c * c).to_num(); assert!( (sum_sq - 1.0).abs() < 0.1, @@ -381,12 +341,9 @@ mod tests { #[test] fn sin_cos_very_large_negative_angle() { - // Test with very large negative angle to exercise the second reduction loop - // The loop iterates while reduced < -pi, adding 2π each time let very_large_neg = I16F16::from_num(-200.0); let (s, c) = sin_cos(very_large_neg); - // Should still satisfy sin²+cos² = 1 let sum_sq: f32 = (s * s + c * c).to_num(); assert!( (sum_sq - 1.0).abs() < 0.1, @@ -396,11 +353,10 @@ mod tests { #[test] fn sin_cos_extreme_angles() { - // Test near I16F16::MAX - previously would produce garbage with iterative reduction + // Iterative reduction once returned garbage near MAX. let extreme = I16F16::MAX - I16F16::from_num(1.0); let (s, c) = sin_cos(extreme); - // Must satisfy Pythagorean identity regardless of input magnitude let sum_sq: f32 = (s * s + c * c).to_num(); assert!( (sum_sq - 1.0).abs() < 0.05, @@ -410,7 +366,6 @@ mod tests { sum_sq ); - // Test negative extreme let neg_extreme = I16F16::MIN + I16F16::from_num(1.0); let (s2, c2) = sin_cos(neg_extreme); let sum_sq2: f32 = (s2 * s2 + c2 * c2).to_num(); @@ -425,24 +380,22 @@ mod tests { #[test] fn sin_cos_known_large_values() { - // 1000 radians ≈ 159.15 * 2π, remainder ≈ 0.96 rad let large = I16F16::from_num(1000.0); let (s, c) = sin_cos(large); - // Compare against f64 reference - let expected_sin = 1000.0_f64.sin() as f32; - let expected_cos = 1000.0_f64.cos() as f32; + let expected_sin = 1000.0_f64.sin(); + let expected_cos = 1000.0_f64.cos(); assert!( - (s.to_num::() - expected_sin).abs() < 0.01, + (s.to_num::() - expected_sin).abs() < 0.01, "sin(1000) = {}, expected {}", - s.to_num::(), + s.to_num::(), expected_sin ); assert!( - (c.to_num::() - expected_cos).abs() < 0.01, + (c.to_num::() - expected_cos).abs() < 0.01, "cos(1000) = {}, expected {}", - c.to_num::(), + c.to_num::(), expected_cos ); } @@ -457,8 +410,7 @@ mod tests { #[test] fn tan_i16f16_near_positive_pole() { - // tan saturates to MAX when approaching π/2 from below - // Saturation occurs within ~3 ULPs of π/2 (≈ 0.00005) + // tan saturates within about 3 ulps of π/2. let far_from_pole = I16F16::from_num(FRAC_PI_2 - 0.0001); let near_pole = I16F16::from_num(FRAC_PI_2 - 0.00003); @@ -474,17 +426,13 @@ mod tests { #[test] fn tan_i16f16_near_negative_pole() { - // tan saturates to MIN when approaching -π/2 from below (more negative) - // The saturation is asymmetric due to fixed-point representation of -π/2 + // Saturation is asymmetric, since -π/2 rounds differently from π/2. let far_from_pole = I16F16::from_num(-FRAC_PI_2 - 0.0001); let near_pole = I16F16::from_num(-FRAC_PI_2 - 0.00005); - // Approaching from below (more negative), tan goes to +∞ then wraps to MIN - // Test the basic property: extreme values near the pole let at_neg_pole_offset = tan(far_from_pole); let closer_to_neg_pole = tan(near_pole); - // Values should have large absolute magnitude near the pole assert!( at_neg_pole_offset.to_num::().abs() > 1000.0, "tan near -π/2 should have large magnitude" @@ -497,14 +445,12 @@ mod tests { #[test] fn tan_i16f16_at_3pi_over_2() { - // tan also has a pole at 3π/2 use core::f64::consts::PI; let pole_3pi2 = 3.0 * PI / 2.0; let far_from_pole = I16F16::from_num(pole_3pi2 - 0.0001); let near_pole = I16F16::from_num(pole_3pi2 - 0.00005); - // Near 3π/2, tan approaches +∞ from below assert!( !is_max_16(tan(far_from_pole)), "tan(3π/2 - 0.0001) should not saturate" @@ -517,15 +463,13 @@ mod tests { } } -/// Angle reduction near the type bounds and inverse-function accuracy. #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod reduction_and_accuracy { use crate::unit::support::Lcg; + use crate::unit::support::TestResult; use fixed::types::{I24F8, I32F32, I64F64}; use fixed_analytics::{CordicNumber, acos, asin, atan, atan2, sin_cos}; - /// Reference reduction with the type's own 2π, in f64. fn reduced_f64(angle: T) -> f64 { let two_pi: f64 = (T::pi() + T::pi()).to_num(); let a: f64 = angle.to_num(); @@ -590,14 +534,14 @@ mod reduction_and_accuracy { } #[test] - fn inverse_functions_track_f64_at_i64f64() { + fn inverse_functions_track_f64_at_i64f64() -> TestResult { let mut rng = Lcg(0xA7A); for _ in 0..2000 { let v = rng.range(-1.0, 1.0); let x = I64F64::from_num(v); - let asin_v: f64 = asin(x).unwrap().to_num(); + let asin_v: f64 = asin(x)?.to_num(); assert!((asin_v - v.asin()).abs() < 2e-15, "asin({v}) = {asin_v}"); - let acos_v: f64 = acos(x).unwrap().to_num(); + let acos_v: f64 = acos(x)?.to_num(); assert!((acos_v - v.acos()).abs() < 2e-15, "acos({v}) = {acos_v}"); let big = (rng.range(-30.0, 30.0)).exp2() * v.signum(); let atan_v: f64 = atan(I64F64::from_num(big)).to_num(); @@ -613,5 +557,6 @@ mod reduction_and_accuracy { "atan2({y}, {x2}) = {atan2_v}" ); } + Ok(()) } } diff --git a/tests/unit/ops/exponential.rs b/tests/unit/ops/exponential.rs index 62f54c1..7946e5b 100644 --- a/tests/unit/ops/exponential.rs +++ b/tests/unit/ops/exponential.rs @@ -1,8 +1,6 @@ -//! Tests for exponential and logarithmic functions - #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod tests { + use crate::unit::support::TestResult; use fixed::types::I16F16; use fixed_analytics::{exp, ln, log2, log10, pow2}; @@ -26,10 +24,11 @@ mod tests { } #[test] - fn ln_special_values() { - assert!(approx_eq(ln(I16F16::ONE).unwrap(), 0.0, TOLERANCE)); - let ln_e: f32 = ln(I16F16::E).unwrap().to_num(); + fn ln_special_values() -> TestResult { + assert!(approx_eq(ln(I16F16::ONE)?, 0.0, TOLERANCE)); + let ln_e: f32 = ln(I16F16::E)?.to_num(); assert!((ln_e - 1.0).abs() < 0.25, "ln(e) = {ln_e}, expected ~1.0"); + Ok(()) } #[test] @@ -40,38 +39,27 @@ mod tests { } #[test] - fn log2_powers_of_two() { - assert!(approx_eq( - log2(I16F16::from_num(1.0)).unwrap(), - 0.0, - TOLERANCE - )); - assert!(approx_eq(log2(I16F16::from_num(2.0)).unwrap(), 1.0, 0.25)); - assert!(approx_eq(log2(I16F16::from_num(4.0)).unwrap(), 2.0, 0.3)); - assert!(approx_eq(log2(I16F16::from_num(8.0)).unwrap(), 3.0, 0.4)); + fn log2_powers_of_two() -> TestResult { + assert!(approx_eq(log2(I16F16::from_num(1.0))?, 0.0, TOLERANCE)); + assert!(approx_eq(log2(I16F16::from_num(2.0))?, 1.0, 0.25)); + assert!(approx_eq(log2(I16F16::from_num(4.0))?, 2.0, 0.3)); + assert!(approx_eq(log2(I16F16::from_num(8.0))?, 3.0, 0.4)); + Ok(()) } #[test] - fn log10_powers_of_ten() { - assert!(approx_eq( - log10(I16F16::from_num(1.0)).unwrap(), - 0.0, - TOLERANCE - )); - assert!(approx_eq(log10(I16F16::from_num(10.0)).unwrap(), 1.0, 0.25)); - assert!(approx_eq( - log10(I16F16::from_num(100.0)).unwrap(), - 2.0, - 0.35 - )); + fn log10_powers_of_ten() -> TestResult { + assert!(approx_eq(log10(I16F16::from_num(1.0))?, 0.0, TOLERANCE)); + assert!(approx_eq(log10(I16F16::from_num(10.0))?, 1.0, 0.25)); + assert!(approx_eq(log10(I16F16::from_num(100.0))?, 2.0, 0.35)); + Ok(()) } #[test] - fn exp_ln_inverse() { - // exp(ln(x)) ≈ x (with limited precision due to CORDIC) + fn exp_ln_inverse() -> TestResult { for i in 1..5 { let x = I16F16::from_num(i); - let result = exp(ln(x).unwrap()); + let result = exp(ln(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -79,14 +67,14 @@ mod tests { "exp(ln({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] - fn pow2_log2_inverse() { - // pow2(log2(x)) ≈ x for positive x + fn pow2_log2_inverse() -> TestResult { for i in 1..8 { let x = I16F16::from_num(i); - let result = pow2(log2(x).unwrap()); + let result = pow2(log2(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -94,59 +82,54 @@ mod tests { "pow2(log2({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] - fn ln_exp_inverse() { - // ln(exp(x)) ≈ x for small x (where exp doesn't overflow) + fn ln_exp_inverse() -> TestResult { for i in -3..=2 { let x = I16F16::from_num(i); let result = ln(exp(x)); assert!(result.is_ok()); let x_f32: f32 = x.to_num(); - let result_f32: f32 = result.unwrap().to_num(); + let result_f32: f32 = result?.to_num(); assert!( (result_f32 - x_f32).abs() < 0.3, "ln(exp({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] - fn ln_near_zero() { - // ln approaches -infinity as x approaches 0 - // Test small positive values + fn ln_near_zero() -> TestResult { let small = I16F16::from_num(0.01); let result = ln(small); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); - // ln(0.01) ≈ -4.605 + let val: f32 = result?.to_num(); assert!(val < -3.0, "ln(0.01) = {val}, expected < -3.0"); - // Very small value let very_small = I16F16::from_num(0.001); let result2 = ln(very_small); assert!(result2.is_ok()); - let val2: f32 = result2.unwrap().to_num(); - // ln(0.001) ≈ -6.908 + let val2: f32 = result2?.to_num(); assert!(val2 < -5.0, "ln(0.001) = {val2}, expected < -5.0"); + Ok(()) } #[test] - fn ln_at_one() { - // ln(1) should be exactly 0 + fn ln_at_one() -> TestResult { let result = ln(I16F16::ONE); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); + let val: f32 = result?.to_num(); assert!(val.abs() < 0.01, "ln(1) = {val}, expected 0"); + Ok(()) } #[test] fn exp_large_negative() { - // exp of large negative values should approach 0 let neg_large = I16F16::from_num(-10.0); let result: f32 = exp(neg_large).to_num(); - // exp(-10) ≈ 0.0000454 assert!( (0.0..0.01).contains(&result), "exp(-10) = {result}, expected ~0" @@ -155,7 +138,6 @@ mod tests { #[test] fn exp_zero() { - // exp(0) should be exactly 1 let result: f32 = exp(I16F16::ZERO).to_num(); assert!( (result - 1.0).abs() < 0.001, @@ -164,83 +146,70 @@ mod tests { } #[test] - fn log2_at_one() { - // log2(1) should be exactly 0 + fn log2_at_one() -> TestResult { let result = log2(I16F16::ONE); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); + let val: f32 = result?.to_num(); assert!(val.abs() < 0.01, "log2(1) = {val}, expected 0"); + Ok(()) } #[test] - fn log10_at_one() { - // log10(1) should be exactly 0 + fn log10_at_one() -> TestResult { let result = log10(I16F16::ONE); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); + let val: f32 = result?.to_num(); assert!(val.abs() < 0.01, "log10(1) = {val}, expected 0"); + Ok(()) } #[test] fn exp_large_positive() { - // exp of large positive values should exercise the argument reduction loop let large = I16F16::from_num(5.0); let result: f32 = exp(large).to_num(); - // exp(5) ≈ 148.41 assert!( result > 100.0 && result < 200.0, "exp(5) = {result}, expected ~148" ); - // Test even larger value to ensure multiple reduction iterations let larger = I16F16::from_num(8.0); let result2: f32 = exp(larger).to_num(); - // exp(8) ≈ 2981 assert!(result2 > 2000.0, "exp(8) = {result2}, expected > 2000"); } #[test] - fn ln_large_values() { - // ln of large values should exercise the argument reduction loop + fn ln_large_values() -> TestResult { let large = I16F16::from_num(1000.0); let result = ln(large); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); - // ln(1000) ≈ 6.9 + let val: f32 = result?.to_num(); assert!((val - 6.9).abs() < 0.5, "ln(1000) = {val}, expected ~6.9"); - // Test even larger value let larger = I16F16::from_num(10000.0); let result2 = ln(larger); assert!(result2.is_ok()); - let val2: f32 = result2.unwrap().to_num(); - // ln(10000) ≈ 9.2 + let val2: f32 = result2?.to_num(); assert!( (val2 - 9.2).abs() < 0.5, "ln(10000) = {val2}, expected ~9.2" ); + Ok(()) } #[test] - fn ln_very_small_values() { - // Test very small values to ensure the multiply-by-2 loop iterates + fn ln_very_small_values() -> TestResult { let tiny = I16F16::from_num(0.0001); let result = ln(tiny); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); - // ln(0.0001) ≈ -9.2 + let val: f32 = result?.to_num(); assert!(val < -8.0, "ln(0.0001) = {val}, expected < -8.0"); + Ok(()) } #[test] fn exp_overflow_to_max() { - // exp of very large positive values should return max when scale > max_shift - // For I16F16: total_bits = 32, max_shift = 31 - // Need scale > 31, i.e., x > 31 * ln(2) ≈ 21.5 - // exp(25) triggers scale=36 which exceeds max_shift=31 let very_large = I16F16::from_num(25.0); let result: f32 = exp(very_large).to_num(); - // Should return max value let max: f32 = I16F16::MAX.to_num(); assert!( (result - max).abs() < 1.0, @@ -250,13 +219,8 @@ mod tests { #[test] fn exp_underflow_to_zero() { - // exp of very large negative values should return zero when -scale > max_shift - // For I16F16: total_bits = 32, max_shift = 31 - // Need -scale > 31, i.e., x < -31 * ln(2) ≈ -21.5 - // exp(-25) triggers scale=-36 which exceeds -max_shift=-31 let very_negative = I16F16::from_num(-25.0); let result: f32 = exp(very_negative).to_num(); - // Should return zero assert!(result == 0.0, "exp(-25) = {result}, expected 0"); } @@ -264,39 +228,29 @@ mod tests { use super::*; use fixed::types::I32F32; - /// Check if I16F16 value is saturated to MAX (within 0.01%) fn is_max_16(val: I16F16) -> bool { val.to_num::() >= I16F16::MAX.to_num::() * 0.9999 } - /// Check if I16F16 value is exactly zero fn is_zero_16(val: I16F16) -> bool { val == I16F16::ZERO } - /// Check if I32F32 value is saturated to MAX (within 0.01%) fn is_max_32(val: I32F32) -> bool { val.to_num::() >= I32F32::MAX.to_num::() * 0.9999 } - /// Check if I32F32 value is exactly zero fn is_zero_32(val: I32F32) -> bool { val == I32F32::ZERO } - // ===== exp saturation thresholds ===== - // The result is zero once e^x rounds to nearest below half an ulp. - // I16F16: saturates to MAX at x >= 10.4, to zero at x <= -11.8 - // I32F32: saturates to MAX at x >= 21.49, to zero at x <= -22.9 - + // Zero once the result rounds below half an ulp. #[test] fn exp_i16f16_upper_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_16(exp(I16F16::from_num(10.39))), "exp(10.39) should not saturate" ); - // At threshold: should saturate assert!( is_max_16(exp(I16F16::from_num(10.4))), "exp(10.4) should saturate to MAX" @@ -305,12 +259,10 @@ mod tests { #[test] fn exp_i16f16_lower_threshold() { - // Above threshold: rounds to one ulp, should NOT be zero assert!( !is_zero_16(exp(I16F16::from_num(-11.7))), "exp(-11.7) should not be zero" ); - // At threshold: rounds below half an ulp, should be zero assert!( is_zero_16(exp(I16F16::from_num(-11.8))), "exp(-11.8) should be zero" @@ -319,15 +271,11 @@ mod tests { #[test] fn exp_deep_negative_band_is_zero() { - // Exercises the floor-adjustment path where the truncated - // quotient sits exactly at -max_shift and the adjustment - // pushes it past: trunc(-21.6/ln2) = -31 = -max_shift for - // I16F16, then r < 0 decrements scale to -32. + // trunc(-21.6/ln2) = -31 for I16F16, and r < 0 then steps scale to -32. assert!( is_zero_16(exp(I16F16::from_num(-21.6))), "exp(-21.6) should be zero" ); - // Same band for I32F32: trunc(-43.7/ln2) = -63 = -max_shift. assert!( is_zero_32(exp(I32F32::from_num(-43.7))), "exp(-43.7) should be zero" @@ -336,12 +284,10 @@ mod tests { #[test] fn exp_i32f32_upper_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_32(exp(I32F32::from_num(21.48))), "exp(21.48) should not saturate" ); - // At threshold: should saturate assert!( is_max_32(exp(I32F32::from_num(21.49))), "exp(21.49) should saturate to MAX" @@ -350,31 +296,22 @@ mod tests { #[test] fn exp_i32f32_lower_threshold() { - // Above threshold: rounds to one ulp, should NOT be zero assert!( !is_zero_32(exp(I32F32::from_num(-22.8))), "exp(-22.8) should not be zero" ); - // At threshold: rounds below half an ulp, should be zero assert!( is_zero_32(exp(I32F32::from_num(-22.9))), "exp(-22.9) should be zero" ); } - // ===== pow2 saturation thresholds ===== - // The result is zero once 2^x rounds to nearest below half an ulp. - // I16F16: saturates to MAX at x >= 15.0, to zero at x <= -17.1 - // I32F32: saturates to MAX at x >= 31.0, to zero at x <= -33.1 - #[test] fn pow2_i16f16_upper_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_16(pow2(I16F16::from_num(14.9))), "pow2(14.9) should not saturate" ); - // At threshold: should saturate assert!( is_max_16(pow2(I16F16::from_num(15.0))), "pow2(15.0) should saturate to MAX" @@ -383,12 +320,10 @@ mod tests { #[test] fn pow2_i16f16_lower_threshold() { - // Above threshold: rounds to one ulp, should NOT be zero assert!( !is_zero_16(pow2(I16F16::from_num(-16.9))), "pow2(-16.9) should not be zero" ); - // At threshold: rounds below half an ulp, should be zero assert!( is_zero_16(pow2(I16F16::from_num(-17.1))), "pow2(-17.1) should be zero" @@ -397,12 +332,10 @@ mod tests { #[test] fn pow2_i32f32_upper_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_32(pow2(I32F32::from_num(30.9))), "pow2(30.9) should not saturate" ); - // At threshold: should saturate assert!( is_max_32(pow2(I32F32::from_num(31.0))), "pow2(31.0) should saturate to MAX" @@ -411,12 +344,10 @@ mod tests { #[test] fn pow2_i32f32_lower_threshold() { - // Above threshold: rounds to one ulp, should NOT be zero assert!( !is_zero_32(pow2(I32F32::from_num(-32.9))), "pow2(-32.9) should not be zero" ); - // At threshold: rounds below half an ulp, should be zero assert!( is_zero_32(pow2(I32F32::from_num(-33.1))), "pow2(-33.1) should be zero" @@ -425,11 +356,10 @@ mod tests { } } -/// Argument reduction, saturation, and the bit-scan normalisation in `ln`. #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod reduction { use crate::unit::support::Lcg; + use crate::unit::support::TestResult; use fixed::types::{I4F60, I16F16, I24F8, I32F32, I64F64}; use fixed_analytics::{CordicNumber, exp, ln, log2, pow, pow2, sqrt}; @@ -460,8 +390,7 @@ mod reduction { #[test] fn exp_estimate_corrected_downward() { - // Just below a negative multiple of ln2 the reciprocal estimate - // rounds up to the multiple; r < 0 then steps k down. + // Just below a negative multiple of ln2 the estimate rounds up; r < 0 steps k down. let x16 = I16F16::LN_2.mul_int(-3) - I16F16::from_bits(1); let got16: f64 = exp(x16).to_num(); assert!((got16 - 0.125).abs() < 3e-5, "exp({x16}) = {got16}"); @@ -522,7 +451,7 @@ mod reduction { } #[test] - fn ln_normalisation_branches() { + fn ln_normalisation_branches() -> TestResult { let ulp = I64F64::from_bits(1); let cases: [(I64F64, &str); 8] = [ (I64F64::from_num(2), "exactly two"), @@ -538,7 +467,7 @@ mod reduction { (I64F64::from_bits(1), "one ulp"), ]; for (x, label) in cases { - let got: f64 = ln(x).unwrap().to_num(); + let got: f64 = ln(x)?.to_num(); let want: f64 = x.to_num::().ln(); assert!( (got - want).abs() < 1e-15 * want.abs().max(1.0), @@ -547,68 +476,53 @@ mod reduction { } for k in 1..=60 { let want = f64::from(k) * core::f64::consts::LN_2; - let got_log2: f64 = log2(I64F64::from_num(1u64 << k)).unwrap().to_num(); + let got_log2: f64 = log2(I64F64::from_num(1u64 << k))?.to_num(); assert!( (got_log2 - f64::from(k)).abs() < 1e-14, "log2(2^{k}) = {got_log2}" ); - let got_ln: f64 = ln(I64F64::from_num(1u64 << k)).unwrap().to_num(); + let got_ln: f64 = ln(I64F64::from_num(1u64 << k))?.to_num(); assert!( (got_ln - want).abs() < 1e-14, "ln(2^{k}) = {got_ln}, want {want}" ); } assert!(ln(I16F16::from_num(4) + I16F16::from_bits(1)).is_ok()); + Ok(()) } #[test] - fn pow_special_cases_and_domain() { - assert_eq!(pow(I16F16::ZERO, I16F16::ZERO).unwrap(), I16F16::ONE); - assert_eq!(pow(I16F16::from_num(7), I16F16::ZERO).unwrap(), I16F16::ONE); - assert_eq!( - pow(I16F16::ZERO, I16F16::from_num(2)).unwrap(), - I16F16::ZERO - ); - assert_eq!(pow(I16F16::ONE, I16F16::from_num(-9)).unwrap(), I16F16::ONE); + fn pow_special_cases_and_domain() -> TestResult { + assert_eq!(pow(I16F16::ZERO, I16F16::ZERO)?, I16F16::ONE); + assert_eq!(pow(I16F16::from_num(7), I16F16::ZERO)?, I16F16::ONE); + assert_eq!(pow(I16F16::ZERO, I16F16::from_num(2))?, I16F16::ZERO); + assert_eq!(pow(I16F16::ONE, I16F16::from_num(-9))?, I16F16::ONE); assert!(pow(I16F16::ZERO, -I16F16::ONE).is_err()); assert!(pow(-I16F16::ONE, I16F16::from_num(0.5)).is_err()); assert!(pow(-I16F16::from_num(2), I16F16::from_num(-1.5)).is_err()); assert_eq!( - pow(I16F16::from_num(200), I16F16::from_num(3)).unwrap(), + pow(I16F16::from_num(200), I16F16::from_num(3))?, I16F16::MAX ); assert_eq!( - pow(-I16F16::from_num(200), I16F16::from_num(3)).unwrap(), + pow(-I16F16::from_num(200), I16F16::from_num(3))?, -I16F16::MAX ); - assert_eq!( - pow(-I16F16::from_num(3), I16F16::ZERO).unwrap(), - I16F16::ONE - ); - assert_eq!( - pow(-I16F16::ONE, I16F16::from_num(7)).unwrap(), - -I16F16::ONE - ); - assert_eq!( - pow(-I16F16::ONE, I16F16::from_num(-8)).unwrap(), - I16F16::ONE - ); + assert_eq!(pow(-I16F16::from_num(3), I16F16::ZERO)?, I16F16::ONE); + assert_eq!(pow(-I16F16::ONE, I16F16::from_num(7))?, -I16F16::ONE); + assert_eq!(pow(-I16F16::ONE, I16F16::from_num(-8))?, I16F16::ONE); for (b, e, want) in [ (-2.0, 3.0, -8.0), (-2.0, 2.0, 4.0), (-2.0, -1.0, -0.5), (-0.5, 3.0, -0.125), ] { - let got16: f64 = pow(I16F16::from_num(b), I16F16::from_num(e)) - .unwrap() - .to_num(); + let got16: f64 = pow(I16F16::from_num(b), I16F16::from_num(e))?.to_num(); assert!( (got16 - want).abs() < 0.05, "I16F16 {b}^{e} = {got16}, want {want}" ); - let got64: f64 = pow(I64F64::from_num(b), I64F64::from_num(e)) - .unwrap() - .to_num(); + let got64: f64 = pow(I64F64::from_num(b), I64F64::from_num(e))?.to_num(); assert!( (got64 - want).abs() < 1e-11, "I64F64 {b}^{e} = {got64}, want {want}" @@ -616,22 +530,21 @@ mod reduction { } // Parity is read from the raw bits, so it survives exponents beyond i32. let huge = I64F64::from_num(1u64 << 40); - assert_eq!(pow(-I64F64::ONE, huge).unwrap(), I64F64::ONE); - assert_eq!(pow(-I64F64::ONE, huge + I64F64::ONE).unwrap(), -I64F64::ONE); - let got: f64 = pow(I16F16::from_num(2), I16F16::from_num(10)) - .unwrap() - .to_num(); + assert_eq!(pow(-I64F64::ONE, huge)?, I64F64::ONE); + assert_eq!(pow(-I64F64::ONE, huge + I64F64::ONE)?, -I64F64::ONE); + let got: f64 = pow(I16F16::from_num(2), I16F16::from_num(10))?.to_num(); assert!((got - 1024.0).abs() < 0.5, "I16F16 2^10 = {got}"); + Ok(()) } #[test] - fn pow_tracks_f64_at_i64f64() { + fn pow_tracks_f64_at_i64f64() -> TestResult { let mut rng = Lcg(0x90); for _ in 0..2000 { let base = I64F64::from_num((rng.range(-20.0, 20.0)).exp2()); let exponent = I64F64::from_num(rng.range(-3.0, 3.0)); let (b, e): (f64, f64) = (base.to_num(), exponent.to_num()); - let got: f64 = pow(base, exponent).unwrap().to_num(); + let got: f64 = pow(base, exponent)?.to_num(); let want = b.powf(e); // Small results are quantised to a few raw ulps. let tol = 8.0f64.mul_add(2f64.powi(-64), 1e-14 * want); @@ -641,24 +554,26 @@ mod reduction { ); } let x = I64F64::from_num(2.5); - let root: f64 = pow(x, I64F64::from_num(0.5)).unwrap().to_num(); - let want: f64 = sqrt(x).unwrap().to_num(); + let root: f64 = pow(x, I64F64::from_num(0.5))?.to_num(); + let want: f64 = sqrt(x)?.to_num(); assert!((root - want).abs() < 1e-15, "2.5^0.5 = {root}, want {want}"); - let inv: f64 = pow(x, -I64F64::ONE).unwrap().to_num(); + let inv: f64 = pow(x, -I64F64::ONE)?.to_num(); assert!((inv - 0.4).abs() < 1e-15, "2.5^-1 = {inv}"); + Ok(()) } #[test] - fn ln_tracks_f64_at_i64f64() { + fn ln_tracks_f64_at_i64f64() -> TestResult { let mut rng = Lcg(0x1A); for _ in 0..3000 { let x = I64F64::from_num((rng.range(-60.0, 62.0)).exp2()); - let got: f64 = ln(x).unwrap().to_num(); + let got: f64 = ln(x)?.to_num(); let want = x.to_num::().ln(); assert!( (got - want).abs() < 1e-15 * want.abs().max(1.0), "ln({x}) = {got}, want {want}" ); } + Ok(()) } } diff --git a/tests/unit/ops/hyperbolic.rs b/tests/unit/ops/hyperbolic.rs index 0dc0f65..28e034c 100644 --- a/tests/unit/ops/hyperbolic.rs +++ b/tests/unit/ops/hyperbolic.rs @@ -1,8 +1,6 @@ -//! Tests for hyperbolic functions - #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod tests { + use crate::unit::support::TestResult; use fixed::types::{I16F16, I32F32}; use fixed_analytics::{acosh, acoth, asinh, atanh, cosh, coth, sinh, sinh_cosh, tanh}; @@ -29,7 +27,6 @@ mod tests { #[test] fn hyperbolic_identity() { - // cosh²(x) - sinh²(x) = 1 for i in -5..=5 { let x = I16F16::from_num(i) * I16F16::from_num(0.2); let (sh, ch) = sinh_cosh(x); @@ -61,7 +58,6 @@ mod tests { #[test] fn acoth_domain_check() { - // acoth requires |x| > 1 assert!(acoth(I16F16::from_num(0.5)).is_err()); assert!(acoth(I16F16::ONE).is_err()); assert!(acoth(-I16F16::ONE).is_err()); @@ -70,30 +66,27 @@ mod tests { } #[test] - fn acoth_values() { - // acoth(x) = atanh(1/x) - // acoth(2) = atanh(0.5) ≈ 0.5493 + fn acoth_values() -> TestResult { let result = acoth(I16F16::from_num(2.0)); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); + let val: f32 = result?.to_num(); assert!( (val - 0.5493).abs() < TOLERANCE, "acoth(2) expected ~0.5493, got {val}" ); - // acoth(-2) = atanh(-0.5) ≈ -0.5493 let result_neg = acoth(I16F16::from_num(-2.0)); assert!(result_neg.is_ok()); - let val_neg: f32 = result_neg.unwrap().to_num(); + let val_neg: f32 = result_neg?.to_num(); assert!( (val_neg + 0.5493).abs() < TOLERANCE, "acoth(-2) expected ~-0.5493, got {val_neg}" ); + Ok(()) } #[test] fn sinh_asinh_roundtrip() { - // sinh(asinh(x)) ≈ x for various x for i in -10..=10 { let x = I16F16::from_num(i) * I16F16::from_num(0.3); let result = sinh(asinh(x)); @@ -107,12 +100,11 @@ mod tests { } #[test] - fn cosh_acosh_roundtrip() { - // cosh(acosh(x)) ≈ x for x >= 1 + fn cosh_acosh_roundtrip() -> TestResult { for i in 1..=10 { let x = I16F16::from_num(i) * I16F16::from_num(0.5); if x >= I16F16::ONE { - let result = cosh(acosh(x).unwrap()); + let result = cosh(acosh(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -121,14 +113,14 @@ mod tests { ); } } + Ok(()) } #[test] - fn tanh_atanh_roundtrip() { - // tanh(atanh(x)) ≈ x for x in (-1, 1) + fn tanh_atanh_roundtrip() -> TestResult { for i in -9..=9 { let x = I16F16::from_num(i) * I16F16::from_num(0.1); - let result = tanh(atanh(x).unwrap()); + let result = tanh(atanh(x)?); let x_f32: f32 = x.to_num(); let result_f32: f32 = result.to_num(); assert!( @@ -136,45 +128,41 @@ mod tests { "tanh(atanh({x_f32})) = {result_f32}, expected {x_f32}" ); } + Ok(()) } #[test] - fn atanh_near_boundary() { - // atanh approaches infinity as |x| approaches 1 - // Test values close to but not at the boundary + fn atanh_near_boundary() -> TestResult { let near_one = I16F16::from_num(0.99); let result = atanh(near_one); assert!(result.is_ok()); - let val: f32 = result.unwrap().to_num(); - // atanh(0.99) ≈ 2.647 + let val: f32 = result?.to_num(); assert!(val > 2.0, "atanh(0.99) = {val}, expected > 2.0"); let neg_near_one = I16F16::from_num(-0.99); let result_neg = atanh(neg_near_one); assert!(result_neg.is_ok()); - let val_neg: f32 = result_neg.unwrap().to_num(); + let val_neg: f32 = result_neg?.to_num(); assert!(val_neg < -2.0, "atanh(-0.99) = {val_neg}, expected < -2.0"); + Ok(()) } #[test] - fn acosh_at_boundary() { - // acosh(1) should be exactly 0 - let result: f32 = acosh(I16F16::ONE).unwrap().to_num(); + fn acosh_at_boundary() -> TestResult { + let result: f32 = acosh(I16F16::ONE)?.to_num(); assert!(result.abs() < 0.01, "acosh(1) = {result}, expected 0"); - // acosh near boundary (1.01) let near_one = I16F16::from_num(1.01); - let result_near: f32 = acosh(near_one).unwrap().to_num(); - // acosh(1.01) ≈ 0.141 + let result_near: f32 = acosh(near_one)?.to_num(); assert!( result_near > 0.0 && result_near < 0.3, "acosh(1.01) = {result_near}, expected ~0.14" ); + Ok(()) } #[test] fn tanh_large_values() { - // tanh should approach ±1 for large values let large = I16F16::from_num(10.0); let result: f32 = tanh(large).to_num(); assert!( @@ -192,11 +180,9 @@ mod tests { #[test] fn sinh_cosh_large_values() { - // Test argument reduction for large values let large = I16F16::from_num(5.0); let (s, c) = sinh_cosh(large); - // Verify cosh²-sinh² = 1 identity let diff: f32 = (c * c - s * s).to_num(); assert!( (diff - 1.0).abs() < 0.2, @@ -206,51 +192,43 @@ mod tests { #[test] fn coth_at_zero() { - // coth(0) is undefined (pole), should return DomainError let result = coth(I16F16::ZERO); assert!(result.is_err(), "coth(0) should return Err"); } #[test] - fn coth_nonzero_values() { - // coth(x) = cosh(x)/sinh(x) - // coth(1) ≈ 1.3130 - let result: f32 = coth(I16F16::ONE).unwrap().to_num(); + fn coth_nonzero_values() -> TestResult { + let result: f32 = coth(I16F16::ONE)?.to_num(); assert!( (result - 1.3130).abs() < TOLERANCE, "coth(1) = {result}, expected ~1.3130" ); - // coth(-1) ≈ -1.3130 - let result_neg: f32 = coth(-I16F16::ONE).unwrap().to_num(); + let result_neg: f32 = coth(-I16F16::ONE)?.to_num(); assert!( (result_neg + 1.3130).abs() < TOLERANCE, "coth(-1) = {result_neg}, expected ~-1.3130" ); + Ok(()) } #[test] fn sinh_cosh_small_values_high_precision() { - // Test Taylor series approximation for high-precision types (≥24 frac bits) - // Uses fifth/sixth-order Taylor series - let small = I32F32::from_num(0.03); // Below 0.05 threshold + let small = I32F32::from_num(0.03); let (s, c) = sinh_cosh(small); - // sinh(0.03) ≈ 0.03 (Taylor: sinh(x) ≈ x) let s_f32: f32 = s.to_num(); assert!( (s_f32 - 0.03).abs() < 0.01, "sinh(0.03) = {s_f32}, expected ~0.03" ); - // cosh(0.03) ≈ 1.00045 (Taylor: 1 + x²/2) let c_f32: f32 = c.to_num(); assert!( (c_f32 - 1.00045).abs() < 0.01, "cosh(0.03) = {c_f32}, expected ~1.00045" ); - // Test negative small value for high precision let small_neg = I32F32::from_num(-0.03); let (s_neg, c_neg) = sinh_cosh(small_neg); let s_neg_f32: f32 = s_neg.to_num(); @@ -264,10 +242,8 @@ mod tests { "cosh(-0.03) = {c_neg_f32}, expected ~1.00045" ); - // Additional test with even smaller value to ensure full Taylor path coverage let tiny = core::hint::black_box(I32F32::from_num(0.01)); let (s_tiny, c_tiny) = sinh_cosh(tiny); - // Use black_box to prevent optimization let s_tiny_f32: f32 = core::hint::black_box(s_tiny).to_num(); let c_tiny_f32: f32 = core::hint::black_box(c_tiny).to_num(); assert!( @@ -282,27 +258,21 @@ mod tests { #[test] fn sinh_cosh_small_values() { - // Test Taylor series approximation for very small values (|x| < 0.1) - // sinh(x) ≈ x for small x - // cosh(x) ≈ 1 + x²/2 for small x let small = I16F16::from_num(0.05); let (s, c) = sinh_cosh(small); - // sinh(0.05) ≈ 0.05 (Taylor: sinh(x) ≈ x) let s_f32: f32 = s.to_num(); assert!( (s_f32 - 0.05).abs() < 0.01, "sinh(0.05) = {s_f32}, expected ~0.05" ); - // cosh(0.05) ≈ 1.00125 (Taylor: 1 + x²/2 = 1 + 0.00125) let c_f32: f32 = c.to_num(); assert!( (c_f32 - 1.00125).abs() < 0.01, "cosh(0.05) = {c_f32}, expected ~1.00125" ); - // Also test negative small value let small_neg = I16F16::from_num(-0.05); let (s_neg, c_neg) = sinh_cosh(small_neg); let s_neg_f32: f32 = s_neg.to_num(); @@ -321,38 +291,28 @@ mod tests { use super::*; use fixed::types::I32F32; - /// Check if I16F16 value is saturated to MAX (within 0.01%) fn is_max_16(val: I16F16) -> bool { val.to_num::() >= I16F16::MAX.to_num::() * 0.9999 } - /// Check if I16F16 value is saturated to MIN (within 0.01%) fn is_min_16(val: I16F16) -> bool { val.to_num::() <= I16F16::MIN.to_num::() * 0.9999 } - /// Check if I32F32 value is saturated to MAX (within 0.01%) fn is_max_32(val: I32F32) -> bool { val.to_num::() >= I32F32::MAX.to_num::() * 0.9999 } - /// Check if I32F32 value is saturated to MIN (within 0.01%) fn is_min_32(val: I32F32) -> bool { val.to_num::() <= I32F32::MIN.to_num::() * 0.9999 } - // ===== sinh saturation thresholds ===== - // I16F16: saturates to MAX at x >= 11.1, to MIN at x <= -11.1 - // I32F32: saturates to MAX at x >= 22.2, to MIN at x <= -22.2 - #[test] fn sinh_i16f16_positive_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_16(sinh(I16F16::from_num(11.0))), "sinh(11.0) should not saturate" ); - // At threshold: should saturate assert!( is_max_16(sinh(I16F16::from_num(11.1))), "sinh(11.1) should saturate to MAX" @@ -361,12 +321,10 @@ mod tests { #[test] fn sinh_i16f16_negative_threshold() { - // Above threshold: should NOT saturate assert!( !is_min_16(sinh(I16F16::from_num(-11.0))), "sinh(-11.0) should not saturate to MIN" ); - // At threshold: should saturate assert!( is_min_16(sinh(I16F16::from_num(-11.1))), "sinh(-11.1) should saturate to MIN" @@ -375,12 +333,10 @@ mod tests { #[test] fn sinh_i32f32_positive_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_32(sinh(I32F32::from_num(22.1))), "sinh(22.1) should not saturate" ); - // At threshold: should saturate assert!( is_max_32(sinh(I32F32::from_num(22.2))), "sinh(22.2) should saturate to MAX" @@ -389,31 +345,22 @@ mod tests { #[test] fn sinh_i32f32_negative_threshold() { - // Above threshold: should NOT saturate assert!( !is_min_32(sinh(I32F32::from_num(-22.1))), "sinh(-22.1) should not saturate to MIN" ); - // At threshold: should saturate assert!( is_min_32(sinh(I32F32::from_num(-22.2))), "sinh(-22.2) should saturate to MIN" ); } - // ===== cosh saturation thresholds ===== - // I16F16: saturates to MAX at |x| >= 11.1 - // I32F32: saturates to MAX at |x| >= 22.2 - // (cosh is even, so positive and negative thresholds are symmetric) - #[test] fn cosh_i16f16_positive_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_16(cosh(I16F16::from_num(11.0))), "cosh(11.0) should not saturate" ); - // At threshold: should saturate assert!( is_max_16(cosh(I16F16::from_num(11.1))), "cosh(11.1) should saturate to MAX" @@ -422,12 +369,10 @@ mod tests { #[test] fn cosh_i16f16_negative_threshold() { - // Above threshold (less negative): should NOT saturate assert!( !is_max_16(cosh(I16F16::from_num(-11.0))), "cosh(-11.0) should not saturate" ); - // At threshold: should saturate assert!( is_max_16(cosh(I16F16::from_num(-11.1))), "cosh(-11.1) should saturate to MAX" @@ -436,12 +381,10 @@ mod tests { #[test] fn cosh_i32f32_positive_threshold() { - // Below threshold: should NOT saturate assert!( !is_max_32(cosh(I32F32::from_num(22.1))), "cosh(22.1) should not saturate" ); - // At threshold: should saturate assert!( is_max_32(cosh(I32F32::from_num(22.2))), "cosh(22.2) should saturate to MAX" @@ -450,24 +393,19 @@ mod tests { #[test] fn cosh_i32f32_negative_threshold() { - // Above threshold (less negative): should NOT saturate assert!( !is_max_32(cosh(I32F32::from_num(-22.1))), "cosh(-22.1) should not saturate" ); - // At threshold: should saturate assert!( is_max_32(cosh(I32F32::from_num(-22.2))), "cosh(-22.2) should saturate to MAX" ); } - // ===== sinh/cosh threshold consistency ===== - #[test] fn sinh_cosh_thresholds_match() { - // sinh and cosh should have the same saturation threshold - // (both grow as e^x/2 for large |x|) + // Both grow as e^x/2, so they saturate together. assert!( is_max_16(sinh(I16F16::from_num(11.1))), "sinh(11.1) should saturate" @@ -488,11 +426,10 @@ mod tests { } } -/// Integer-bit-poor types, and the logarithmic inverse forms. #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod wide_and_narrow { use crate::unit::support::Lcg; + use crate::unit::support::TestResult; use fixed::types::{I4F12, I4F60, I16F16, I32F32, I64F64}; use fixed_analytics::{acosh, asinh, sinh_cosh}; @@ -521,7 +458,7 @@ mod wide_and_narrow { } #[test] - fn asinh_acosh_track_f64_at_i64f64() { + fn asinh_acosh_track_f64_at_i64f64() -> TestResult { let mut rng = Lcg(0xA5); for _ in 0..2000 { let sign = if rng.unit() < 0.5 { -1.0 } else { 1.0 }; @@ -534,7 +471,7 @@ mod wide_and_narrow { ); // acosh is ill-conditioned near 1; start where f64 is trustworthy. let w = I64F64::from_num(1.0 + (rng.range(-10.0, 61.0)).exp2()); - let got_w: f64 = acosh(w).unwrap().to_num(); + let got_w: f64 = acosh(w)?.to_num(); let want_w = w.to_num::().acosh(); assert!( (got_w - want_w).abs() < 1e-13 * want_w.abs().max(1.0), @@ -543,13 +480,14 @@ mod wide_and_narrow { } // Both branches around the 1.5 threshold of acosh. for w in [1.01f64, 1.25, 1.49, 1.5, 1.51, 2.0] { - let got: f64 = acosh(I64F64::from_num(w)).unwrap().to_num(); + let got: f64 = acosh(I64F64::from_num(w))?.to_num(); assert!((got - w.acosh()).abs() < 1e-9, "acosh({w}) = {got}"); } + Ok(()) } #[test] - fn asinh_acosh_near_the_type_bounds() { + fn asinh_acosh_near_the_type_bounds() -> TestResult { // The product would overflow, so the logarithm is split. for (x, want) in [ (I16F16::MAX, (2.0 * I16F16::MAX.to_num::()).ln()), @@ -561,17 +499,17 @@ mod wide_and_narrow { } for x in [I16F16::MAX, I16F16::from_num(20_000)] { let want = (2.0f64 * x.to_num::()).ln(); - let got: f64 = acosh(x).unwrap().to_num(); + let got: f64 = acosh(x)?.to_num(); assert!((got - want).abs() < 2e-4, "acosh({x}) = {got}, want {want}"); } let want_32 = (2.0f64 * I32F32::MAX.to_num::()).ln(); - let got_32: f64 = acosh(I32F32::MAX).unwrap().to_num(); + let got_32: f64 = acosh(I32F32::MAX)?.to_num(); assert!( (got_32 - want_32).abs() < 1e-8, "acosh(I32F32::MAX) = {got_32}" ); let want_64 = (2.0f64 * I64F64::MAX.to_num::()).ln(); - let got_64: f64 = acosh(I64F64::MAX).unwrap().to_num(); + let got_64: f64 = acosh(I64F64::MAX)?.to_num(); assert!( (got_64 - want_64).abs() < 1e-12, "acosh(I64F64::MAX) = {got_64}" @@ -582,6 +520,7 @@ mod wide_and_narrow { "asinh(MAX) = {got_asinh}, want {want_64}" ); assert_eq!(asinh(-I64F64::from_num(3)), -asinh(I64F64::from_num(3))); + Ok(()) } #[test] diff --git a/tests/unit/ops/mod.rs b/tests/unit/ops/mod.rs index 8a461a1..f8f5bc8 100644 --- a/tests/unit/ops/mod.rs +++ b/tests/unit/ops/mod.rs @@ -1,5 +1,3 @@ -//! Tests for mathematical operations - mod algebraic; mod circular; mod exponential; diff --git a/tests/unit/smoke.rs b/tests/unit/smoke.rs index 8b13594..2a28094 100644 --- a/tests/unit/smoke.rs +++ b/tests/unit/smoke.rs @@ -1,8 +1,7 @@ -//! Smoke tests and multi-type tests for the library API - #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod tests { + use crate::unit::support::TestResult; + use core::hint::black_box; use fixed::types::I16F16; use fixed_analytics::{ acos, acosh, acoth, asin, asinh, atan, atan2, atanh, cos, cosh, coth, exp, ln, log2, log10, @@ -12,69 +11,69 @@ mod tests { #[test] fn smoke_test_trig() { let angle = I16F16::from_num(0.5); - let _ = sin(angle); - let _ = cos(angle); - let _ = tan(angle); - let _ = sin_cos(angle); + black_box(sin(angle)); + black_box(cos(angle)); + black_box(tan(angle)); + black_box(sin_cos(angle)); } #[test] - fn smoke_test_inverse_trig() { + fn smoke_test_inverse_trig() -> TestResult { let x = I16F16::from_num(0.5); - let _ = asin(x); - let _ = acos(x); - let _ = atan(x); - let _ = atan2(x, I16F16::ONE); + black_box(asin(x)?); + black_box(acos(x)?); + black_box(atan(x)); + black_box(atan2(x, I16F16::ONE)); + Ok(()) } #[test] - fn smoke_test_hyperbolic() { + fn smoke_test_hyperbolic() -> TestResult { let x = I16F16::from_num(0.5); - let _ = sinh(x); - let _ = cosh(x); - let _ = tanh(x); - let _ = coth(x); - let _ = sinh_cosh(x); + black_box(sinh(x)); + black_box(cosh(x)); + black_box(tanh(x)); + black_box(coth(x)?); + black_box(sinh_cosh(x)); + Ok(()) } #[test] - fn smoke_test_inverse_hyperbolic() { + fn smoke_test_inverse_hyperbolic() -> TestResult { let x = I16F16::from_num(0.5); - let _ = asinh(x); - let _ = atanh(x); + black_box(asinh(x)); + black_box(atanh(x)?); let x_large = I16F16::from_num(1.5); - let _ = acosh(x_large); - let _ = acoth(x_large); + black_box(acosh(x_large)?); + black_box(acoth(x_large)?); + Ok(()) } #[test] - fn smoke_test_exponential() { + fn smoke_test_exponential() -> TestResult { let x = I16F16::from_num(0.5); - let _ = exp(x); - let _ = ln(x); - let _ = log2(x); - let _ = log10(x); + black_box(exp(x)); + black_box(ln(x)?); + black_box(log2(x)?); + black_box(log10(x)?); + Ok(()) } #[test] - fn smoke_test_algebraic() { + fn smoke_test_algebraic() -> TestResult { let x = I16F16::from_num(2.0); - let _ = sqrt(x).unwrap(); + black_box(sqrt(x)?); + Ok(()) } } -// ========================================================================== -// Multi-type tests: verify functions work with different fixed-point types -// ========================================================================== - #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod multi_type { + use crate::unit::support::TestResult; use fixed::types::{I8F24, I32F32}; use fixed_analytics::{acos, asin, atan, exp, ln, sin_cos, sinh_cosh, sqrt}; - // I32F32 tests - higher precision (32 fractional bits) #[test] fn trig_i32f32() { let angle = I32F32::from_num(0.5); @@ -88,7 +87,6 @@ mod multi_type { #[test] fn trig_identity_i32f32() { - // sin²(x) + cos²(x) = 1 with higher precision for i in -10..=10 { let angle = I32F32::from_num(i) * I32F32::from_num(0.3); let (s, c) = sin_cos(angle); @@ -104,13 +102,12 @@ mod multi_type { } #[test] - fn inverse_trig_i32f32() { + fn inverse_trig_i32f32() -> TestResult { let x = I32F32::from_num(0.5); - let asin_val = asin(x).unwrap(); - let acos_val = acos(x).unwrap(); + let asin_val = asin(x)?; + let acos_val = acos(x)?; let atan_val = atan(x); - // Verify results are in expected range let asin_f64: f64 = asin_val.to_num(); assert!(asin_f64 > 0.5 && asin_f64 < 0.6, "asin(0.5) ≈ 0.524"); @@ -119,13 +116,13 @@ mod multi_type { let atan_f64: f64 = atan_val.to_num(); assert!(atan_f64 > 0.45 && atan_f64 < 0.5, "atan(0.5) ≈ 0.464"); + Ok(()) } #[test] fn hyperbolic_i32f32() { let x = I32F32::from_num(0.5); let (s, c) = sinh_cosh(x); - // cosh²(x) - sinh²(x) = 1 let diff: f64 = (c * c - s * s).to_num(); assert!( (diff - 1.0).abs() < 0.01, @@ -134,30 +131,30 @@ mod multi_type { } #[test] - fn exp_ln_i32f32() { + fn exp_ln_i32f32() -> TestResult { let x = I32F32::from_num(2.0); - let result = exp(ln(x).unwrap()); + let result = exp(ln(x)?); let result_f64: f64 = result.to_num(); assert!( (result_f64 - 2.0).abs() < 0.3, "exp(ln(2)) = {result_f64} (I32F32), expected ~2.0" ); + Ok(()) } #[test] - fn sqrt_i32f32() { + fn sqrt_i32f32() -> TestResult { let x = I32F32::from_num(4.0); - let result: f64 = sqrt(x).unwrap().to_num(); + let result: f64 = sqrt(x)?.to_num(); assert!( (result - 2.0).abs() < 0.001, "sqrt(4) = {result} (I32F32), expected 2.0" ); + Ok(()) } - // I8F24 tests - small integer part, high fractional precision #[test] fn trig_i8f24() { - // I8F24 can represent values from about -128 to 127 with 24 fractional bits let angle = I8F24::from_num(0.5); let (s, c) = sin_cos(angle); let sum_sq: f32 = (s * s + c * c).to_num(); @@ -179,40 +176,39 @@ mod multi_type { } #[test] - fn sqrt_i8f24() { + fn sqrt_i8f24() -> TestResult { let x = I8F24::from_num(2.0); - let result: f32 = sqrt(x).unwrap().to_num(); + let result: f32 = sqrt(x)?.to_num(); assert!( (result - 1.414).abs() < 0.01, "sqrt(2) = {result} (I8F24), expected ~1.414" ); + Ok(()) } #[test] - fn ln_i8f24() { + fn ln_i8f24() -> TestResult { let x = I8F24::from_num(2.0); - let result: f32 = ln(x).unwrap().to_num(); + let result: f32 = ln(x)?.to_num(); assert!( (result - 0.693).abs() < 0.2, "ln(2) = {result} (I8F24), expected ~0.693" ); + Ok(()) } } -// I8F8 tests - lower precision, 16-bit total #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod i8f8 { + use crate::unit::support::TestResult; use fixed::types::I8F8; use fixed_analytics::{sin_cos, sqrt}; #[test] fn basic_trig() { - // I8F8 has limited precision but should still work let angle = I8F8::from_num(0.5); let (s, c) = sin_cos(angle); let sum_sq: f32 = (s * s + c * c).to_num(); - // Lower precision type, allow more tolerance assert!( (sum_sq - 1.0).abs() < 0.2, "sin²+cos² = {sum_sq} (I8F8), expected ~1.0" @@ -220,12 +216,13 @@ mod i8f8 { } #[test] - fn basic_sqrt() { + fn basic_sqrt() -> TestResult { let x = I8F8::from_num(4.0); - let result: f32 = sqrt(x).unwrap().to_num(); + let result: f32 = sqrt(x)?.to_num(); assert!( (result - 2.0).abs() < 0.1, "sqrt(4) = {result} (I8F8), expected ~2.0" ); + Ok(()) } } diff --git a/tests/unit/support.rs b/tests/unit/support.rs index 9ad7257..3202f9a 100644 --- a/tests/unit/support.rs +++ b/tests/unit/support.rs @@ -1,6 +1,5 @@ -//! Shared helpers for the unit tests. +pub type TestResult = Result<(), Box>; -/// Deterministic generator for reproducible sweeps. pub struct Lcg(pub u64); impl Lcg { @@ -13,16 +12,10 @@ impl Lcg { (x ^ (x >> 29)).wrapping_mul(0x9E37_79B9_7F4A_7C15) ^ (x >> 32) } - /// Uniform in [0, 1). - #[allow( - clippy::cast_precision_loss, - reason = "53 bits of randomness is plenty" - )] - pub const fn unit(&mut self) -> f64 { - (self.next_u64() >> 11) as f64 / (1u64 << 53) as f64 + pub fn unit(&mut self) -> f64 { + fixed::types::U0F64::from_bits(self.next_u64() >> 11 << 11).to_num() } - /// Uniform in [lo, hi). pub fn range(&mut self, lo: f64, hi: f64) -> f64 { (hi - lo).mul_add(self.unit(), lo) } diff --git a/tests/unit/tables/chebyshev.rs b/tests/unit/tables/chebyshev.rs index 566826d..fc5c51b 100644 --- a/tests/unit/tables/chebyshev.rs +++ b/tests/unit/tables/chebyshev.rs @@ -1,26 +1,14 @@ -//! Tests for Chebyshev polynomial coefficient tables - #[cfg(test)] -#[allow( - clippy::cast_precision_loss, - reason = "test code uses f64 casts for verification" -)] mod tests { use fixed_analytics::tables::chebyshev::{COS_Q_HI, COS_Q_LO, SIN_P_HI, SIN_P_LO}; - const SCALE: f64 = (1_u64 << 63) as f64; - fn i1f63_to_f64(bits: i64) -> f64 { - (bits as f64) / SCALE + fixed::types::I1F63::from_bits(bits).to_num() } - // The constant term (last element, c₀) of each polynomial should match - // the Taylor series leading coefficient to high precision. This guards - // against coefficient ordering mistakes and regeneration drift. - + // c₀ must match the Taylor coefficient, which catches ordering mistakes and drift. #[test] fn sin_constant_term_is_neg_one_sixth() { - // (sin(x)-x)/x³ → -1/6 at x=0 let expected = -1.0 / 6.0; let lo = i1f63_to_f64(*SIN_P_LO.last().unwrap_or(&0)); let hi = i1f63_to_f64(*SIN_P_HI.last().unwrap_or(&0)); @@ -36,7 +24,6 @@ mod tests { #[test] fn cos_constant_term_is_neg_one_half() { - // (cos(x)-1)/x² → -1/2 at x=0 let expected = -0.5; let lo = i1f63_to_f64(*COS_Q_LO.last().unwrap_or(&0)); let hi = i1f63_to_f64(*COS_Q_HI.last().unwrap_or(&0)); @@ -74,20 +61,13 @@ mod tests { } } -/// `horner` wraps on the strength of every intermediate staying below 1 in -/// magnitude: check the bound and the bit-identity with saturating Horner. +/// `horner` wraps safely only while every intermediate stays below 1. #[cfg(test)] -#[allow( - clippy::cast_precision_loss, - clippy::unwrap_used, - reason = "test code uses f64 for verification" -)] mod wrapping_invariant { - use fixed::types::{I16F16, I32F32, I64F64}; + use fixed::types::{I1F63, I16F16, I32F32, I64F64}; use fixed_analytics::CordicNumber; use fixed_analytics::tables::chebyshev::{COS_Q_HI, COS_Q_LO, SIN_P_HI, SIN_P_LO, horner}; - const SCALE: f64 = (1_u64 << 63) as f64; const U_MAX: f64 = core::f64::consts::FRAC_PI_4 * core::f64::consts::FRAC_PI_4; fn max_intermediate(coeffs: &[i64]) -> f64 { @@ -95,10 +75,10 @@ mod wrapping_invariant { for i in 0..=1000 { let u = U_MAX * f64::from(i) / 1000.0; let mut iter = coeffs.iter(); - let mut acc = *iter.next().unwrap_or(&0) as f64 / SCALE; + let mut acc: f64 = iter.next().map_or(0.0, |&c| I1F63::from_bits(c).to_num()); worst = worst.max(acc.abs()); for &c in iter { - acc = u.mul_add(acc, c as f64 / SCALE); + acc = u.mul_add(acc, I1F63::from_bits(c).to_num()); worst = worst.max(acc.abs()); worst = worst.max((u * acc).abs()); } @@ -116,7 +96,7 @@ mod wrapping_invariant { fn saturating_horner(coeffs: &[i64; N], x: T) -> T { let mut iter = coeffs.iter(); - let mut result = T::from_i1f63(*iter.next().unwrap()); + let mut result = iter.next().map_or_else(T::zero, |&c| T::from_i1f63(c)); for &c in iter { result = T::from_i1f63(c).saturating_add(x.saturating_mul(result)); } diff --git a/tests/unit/tables/circular.rs b/tests/unit/tables/circular.rs index a5b1aff..f506709 100644 --- a/tests/unit/tables/circular.rs +++ b/tests/unit/tables/circular.rs @@ -1,11 +1,4 @@ -//! Tests for circular CORDIC lookup tables - #[cfg(test)] -#[allow( - clippy::indexing_slicing, - clippy::cast_precision_loss, - reason = "test code uses direct indexing and f64 casts" -)] mod tests { use fixed_analytics::tables::circular::ATAN_TABLE; @@ -16,19 +9,15 @@ mod tests { #[test] fn atan_table_spot_check() { - // ATAN_TABLE[0] = atan(2^0) = atan(1) = π/4 ≈ 0.7853981633974483 - // In I1F63 format: π/4 * 2^63 - // Use u64 for 2^63 to avoid i64 overflow - let scale = (1_u64 << 63) as f64; - let atan_1: f64 = (ATAN_TABLE[0] as f64) / scale; + let [t0, t1, ..] = ATAN_TABLE; + let atan_1: f64 = fixed::types::I1F63::from_bits(t0).to_num(); let expected = core::f64::consts::FRAC_PI_4; assert!( (atan_1 - expected).abs() < 1e-15, "atan(1) = {atan_1}, expected {expected}" ); - // ATAN_TABLE[1] = atan(2^(-1)) = atan(0.5) ≈ 0.4636476090008061 - let atan_half: f64 = (ATAN_TABLE[1] as f64) / scale; + let atan_half: f64 = fixed::types::I1F63::from_bits(t1).to_num(); let expected_half = 0.5_f64.atan(); assert!( (atan_half - expected_half).abs() < 1e-15, @@ -38,31 +27,27 @@ mod tests { #[test] fn atan_table_decreasing_until_convergence() { - // Each entry should be smaller than or equal to the previous - // (last few entries may be equal due to precision limits) - for i in 1..ATAN_TABLE.len() { + for (i, (prev, cur)) in ATAN_TABLE.iter().zip(ATAN_TABLE.iter().skip(1)).enumerate() { assert!( - ATAN_TABLE[i] <= ATAN_TABLE[i - 1], - "ATAN_TABLE[{}] = {} should be <= ATAN_TABLE[{}] = {}", - i, - ATAN_TABLE[i], - i - 1, - ATAN_TABLE[i - 1] + cur <= prev, + "ATAN_TABLE[{}] = {cur} should be <= ATAN_TABLE[{i}] = {prev}", + i + 1 ); } } #[test] fn atan_table_strictly_decreasing_early() { - // First 60 entries should be strictly decreasing - for i in 1..60 { + for (i, (prev, cur)) in ATAN_TABLE + .iter() + .zip(ATAN_TABLE.iter().skip(1)) + .enumerate() + .take(59) + { assert!( - ATAN_TABLE[i] < ATAN_TABLE[i - 1], - "ATAN_TABLE[{}] = {} should be < ATAN_TABLE[{}] = {}", - i, - ATAN_TABLE[i], - i - 1, - ATAN_TABLE[i - 1] + cur < prev, + "ATAN_TABLE[{}] = {cur} should be < ATAN_TABLE[{i}] = {prev}", + i + 1 ); } } diff --git a/tests/unit/tables/hyperbolic.rs b/tests/unit/tables/hyperbolic.rs index 6d6cf53..4ad5006 100644 --- a/tests/unit/tables/hyperbolic.rs +++ b/tests/unit/tables/hyperbolic.rs @@ -1,15 +1,7 @@ -//! Tests for hyperbolic CORDIC lookup tables - #[cfg(test)] -#[allow( - clippy::indexing_slicing, - clippy::cast_precision_loss, - reason = "test code uses direct indexing and f64 casts" -)] mod tests { use fixed_analytics::tables::hyperbolic::{ATANH_HALF, ATANH_TABLE, needs_repeat}; - /// Repeat indices for hyperbolic CORDIC convergence (used only in tests). const REPEAT_INDICES: [u32; 5] = [4, 13, 40, 121, 364]; #[test] @@ -19,42 +11,40 @@ mod tests { #[test] fn atanh_table_decreasing_until_convergence() { - // Each entry should be smaller than or equal to the previous - // (last few entries may be equal due to precision limits) - for i in 1..ATANH_TABLE.len() { + for (i, (prev, cur)) in ATANH_TABLE + .iter() + .zip(ATANH_TABLE.iter().skip(1)) + .enumerate() + { assert!( - ATANH_TABLE[i] <= ATANH_TABLE[i - 1], - "ATANH_TABLE[{}] = {} should be <= ATANH_TABLE[{}] = {}", - i, - ATANH_TABLE[i], - i - 1, - ATANH_TABLE[i - 1] + cur <= prev, + "ATANH_TABLE[{}] = {cur} should be <= ATANH_TABLE[{i}] = {prev}", + i + 1 ); } } #[test] fn atanh_table_strictly_decreasing_early() { - // First 60 entries should be strictly decreasing - for i in 1..60 { + for (i, (prev, cur)) in ATANH_TABLE + .iter() + .zip(ATANH_TABLE.iter().skip(1)) + .enumerate() + .take(59) + { assert!( - ATANH_TABLE[i] < ATANH_TABLE[i - 1], - "ATANH_TABLE[{}] = {} should be < ATANH_TABLE[{}] = {}", - i, - ATANH_TABLE[i], - i - 1, - ATANH_TABLE[i - 1] + cur < prev, + "ATANH_TABLE[{}] = {cur} should be < ATANH_TABLE[{i}] = {prev}", + i + 1 ); } } #[test] fn needs_repeat_correct_indices() { - // Verify needs_repeat returns true for expected indices for &idx in &REPEAT_INDICES { assert!(needs_repeat(idx), "needs_repeat({idx}) should be true"); } - // Verify needs_repeat returns false for non-repeat indices assert!(!needs_repeat(5)); assert!(!needs_repeat(100)); assert!(!needs_repeat(0)); @@ -63,25 +53,20 @@ mod tests { #[test] fn atanh_half_matches_table() { - // ATANH_HALF should equal ATANH_TABLE[0] - assert_eq!(ATANH_HALF, ATANH_TABLE[0]); + assert_eq!(Some(&ATANH_HALF), ATANH_TABLE.first()); } #[test] fn atanh_table_spot_check() { - // ATANH_TABLE[0] = atanh(2^(-1)) = atanh(0.5) ≈ 0.5493061443340548 - // In I1F63 format: atanh(0.5) * 2^63 - // Use u64 for 2^63 to avoid i64 overflow - let scale = (1_u64 << 63) as f64; - let atanh_half: f64 = (ATANH_TABLE[0] as f64) / scale; + let [t0, t1, ..] = ATANH_TABLE; + let atanh_half: f64 = fixed::types::I1F63::from_bits(t0).to_num(); let expected = 0.5_f64.atanh(); assert!( (atanh_half - expected).abs() < 1e-15, "atanh(0.5) = {atanh_half}, expected {expected}" ); - // ATANH_TABLE[1] = atanh(2^(-2)) = atanh(0.25) ≈ 0.2554128118829953 - let atanh_quarter: f64 = (ATANH_TABLE[1] as f64) / scale; + let atanh_quarter: f64 = fixed::types::I1F63::from_bits(t1).to_num(); let expected_quarter = 0.25_f64.atanh(); assert!( (atanh_quarter - expected_quarter).abs() < 1e-15, diff --git a/tests/unit/tables/mod.rs b/tests/unit/tables/mod.rs index a0b07c8..63b3277 100644 --- a/tests/unit/tables/mod.rs +++ b/tests/unit/tables/mod.rs @@ -1,5 +1,3 @@ -//! Tests for lookup tables - mod chebyshev; mod circular; mod hyperbolic; diff --git a/tests/unit/traits.rs b/tests/unit/traits.rs index a0c556e..d121d06 100644 --- a/tests/unit/traits.rs +++ b/tests/unit/traits.rs @@ -1,17 +1,14 @@ -//! Tests for `CordicNumber` trait implementations - #[cfg(test)] mod tests { use fixed::types::{I4F12, I4F60, I8F8, I8F24, I16F16, I20F12, I24F8, I32F32, I48F16, I64F64}; use fixed_analytics::CordicNumber; #[test] - #[allow(clippy::approx_constant, reason = "testing pi approximation")] fn basic_operations_i16f16() { let x = I16F16::from_num(2.5); assert_eq!(I16F16::zero(), I16F16::ZERO); assert_eq!(I16F16::one(), I16F16::ONE); - assert!(I16F16::pi() > I16F16::from_num(3.14)); + assert!(I16F16::pi() > I16F16::from_num(314) / 100); assert!(I16F16::pi() < I16F16::from_num(3.15)); assert!(!x.is_negative()); assert!((-x).is_negative()); @@ -19,33 +16,23 @@ mod tests { #[test] fn generic_impl_works_for_various_types() { - // Test that the generic impl works for various type configurations - - // I8F24: 8 integer bits, 24 fractional bits (32-bit total) let _: I8F24 = I8F24::pi(); let _: I8F24 = I8F24::frac_pi_2(); - // I24F8: 24 integer bits, 8 fractional bits let _: I24F8 = I24F8::pi(); let _: I24F8 = I24F8::frac_pi_2(); - // I4F12: 4 integer bits, 12 fractional bits (16-bit total) let _: I4F12 = I4F12::pi(); - // I20F12: 20 integer bits, 12 fractional bits (32-bit total) let _: I20F12 = I20F12::pi(); - // I48F16: 48 integer bits, 16 fractional bits (64-bit total) let _: I48F16 = I48F16::pi(); - // I4F60: 4 integer bits, 60 fractional bits (64-bit total) - high precision let _: I4F60 = I4F60::pi(); } #[test] fn from_i1f63_across_types() { - // Test that from_i1f63 works correctly across different types - // 0.5 in I1F63 format let half_bits: i64 = 0x4000_0000_0000_0000; let i8f8_half: f32 = I8F8::from_i1f63(half_bits).to_num(); @@ -60,7 +47,6 @@ mod tests { let i64f64_half: f64 = I64F64::from_i1f63(half_bits).to_num(); assert!((i64f64_half - 0.5).abs() < 1e-15); - // Also test a non-standard type let i24f8_half: f32 = I24F8::from_i1f63(half_bits).to_num(); assert!((i24f8_half - 0.5).abs() < 0.01); @@ -75,7 +61,6 @@ mod tests { assert_eq!(I32F32::frac_bits(), 32); assert_eq!(I64F64::frac_bits(), 64); - // Non-standard types assert_eq!(I8F24::frac_bits(), 24); assert_eq!(I24F8::frac_bits(), 8); assert_eq!(I4F12::frac_bits(), 12); @@ -84,8 +69,7 @@ mod tests { #[test] fn div_overflow_saturates_to_min() { - // Dividing a negative number by zero should saturate to MIN - // (signs disagree: negative / "positive zero"). + // Signs disagree, since zero counts as positive. let neg = I16F16::from_num(-1); let zero = I16F16::ZERO; assert_eq!(neg.div(zero), I16F16::MIN); @@ -100,8 +84,6 @@ mod tests { #[test] fn frac_pi_4_values() { - // Test the frac_pi_4() default implementation - // π/4 ≈ 0.7854 let pi_4_16: f32 = I16F16::frac_pi_4().to_num(); assert!( (pi_4_16 - core::f32::consts::FRAC_PI_4).abs() < 0.001, @@ -122,9 +104,7 @@ mod tests { } } -/// The integer and wrapping helpers added for the fast paths. #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod fast_path_helpers { use fixed::types::{I4F4, I8F8, I16F16, I32F32, I64F64}; use fixed_analytics::CordicNumber; diff --git a/tests/unit/verification.rs b/tests/unit/verification.rs index f41ae9a..1a32d72 100644 --- a/tests/unit/verification.rs +++ b/tests/unit/verification.rs @@ -1,22 +1,6 @@ -//! Mathematical correctness verification tests. -//! -//! These tests verify the CORDIC implementations against known mathematical -//! properties: reference values, identities, inverse roundtrips, monotonicity, -//! and output bounds. - -// Test-specific lints - these are acceptable in test code -#![allow( - clippy::unwrap_used, - clippy::cast_possible_truncation, - clippy::suboptimal_flops, - clippy::cast_lossless, - clippy::manual_range_contains, - reason = "test code uses these patterns for conciseness and clarity" -)] - #[cfg(test)] mod reference_comparison { - //! Compare against f64 reference implementations across sampled inputs. + use crate::unit::support::TestResult; use fixed::types::I16F16; use fixed_analytics::{ @@ -24,25 +8,21 @@ mod reference_comparison { sqrt, tan, tanh, }; - /// Deterministic pseudo-random bit generator for reproducible sampling. fn sample_bits(seed: u64, index: u64) -> i32 { - // Simple LCG-style mixing let mut x = seed.wrapping_add(index.wrapping_mul(0x9E37_79B9_7F4A_7C15)); x = x.wrapping_mul(0x5851_F42D_4C95_7F2D); x = x ^ (x >> 32); - x as i32 + u32::try_from(x & 0xFFFF_FFFF).map_or(0, u32::cast_signed) } const SAMPLES: u64 = 2000; const SEED: u64 = 0xDEAD_BEEF_CAFE_BABE; - // Error tolerance: I16F16 has 16 fractional bits, so ~4-5 decimal digits - // Allow 2^-12 ≈ 0.000244 absolute error for well-conditioned functions const TRIG_TOL: f64 = 0.001; const SQRT_TOL: f64 = 0.0005; - const EXP_TOL: f64 = 0.01; // exp amplifies errors + const EXP_TOL: f64 = 0.01; const LOG_TOL: f64 = 0.01; - const HYPER_TOL: f64 = 0.02; // hyperbolic has lower precision due to iteration repeats + const HYPER_TOL: f64 = 0.02; #[test] fn sin_cos_vs_f64() { @@ -54,7 +34,6 @@ mod reference_comparison { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip extreme values where f64 comparison becomes meaningless if x_f64.abs() > 100.0 { continue; } @@ -80,9 +59,6 @@ mod reference_comparison { c.to_num::() ); } - - // Uncomment to see actual max errors during development: - // println!("Max sin error: {max_sin_err}, max cos error: {max_cos_err}"); } #[test] @@ -92,12 +68,10 @@ mod reference_comparison { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip near poles and extreme values let frac_pi_2 = core::f64::consts::FRAC_PI_2; if x_f64.abs() > 50.0 { continue; } - // Skip near ±π/2 where tan explodes let reduced = x_f64 % core::f64::consts::PI; if (reduced.abs() - frac_pi_2).abs() < 0.1 { continue; @@ -106,7 +80,6 @@ mod reference_comparison { let result = tan(x); let expected = x_f64.tan(); - // Use relative error for large values let err = if expected.abs() > 1.0 { (result.to_num::() - expected).abs() / expected.abs() } else { @@ -114,7 +87,7 @@ mod reference_comparison { }; assert!( - err < TRIG_TOL * 2.0, // tan accumulates more error + err < TRIG_TOL * 2.0, "tan({x_f64}): got {}, expected {expected}", result.to_num::() ); @@ -154,7 +127,6 @@ mod reference_comparison { let y_f64: f64 = y.to_num(); let x_f64: f64 = x.to_num(); - // Skip extreme values if y_f64.abs() > 1000.0 || x_f64.abs() > 1000.0 { continue; } @@ -173,11 +145,8 @@ mod reference_comparison { #[test] fn asin_acos_vs_f64() { - // Test over domain [-1, 1] for i in 0..500 { - // Map to [-1, 1] - let t = (i as f64) / 499.0 * 2.0 - 1.0; - // Avoid exact ±1 where precision is worst + let t = f64::mul_add(f64::from(i) / 499.0, 2.0, -1.0); let t = t * 0.999; let x = I16F16::from_num(t); let x_f64: f64 = x.to_num(); @@ -205,19 +174,17 @@ mod reference_comparison { } #[test] - fn sqrt_vs_f64() { + fn sqrt_vs_f64() -> TestResult { for i in 0..SAMPLES { let bits = sample_bits(SEED, i); - // Only test non-negative values let bits = bits & 0x7FFF_FFFF; let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - let result = sqrt(x).unwrap(); + let result = sqrt(x)?; let expected = x_f64.sqrt(); let err = (result.to_num::() - expected).abs(); - // Use relative error for larger values let tol = if expected > 1.0 { SQRT_TOL * expected } else { @@ -230,6 +197,7 @@ mod reference_comparison { result.to_num::() ); } + Ok(()) } #[test] @@ -239,8 +207,7 @@ mod reference_comparison { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip values where exp would overflow I16F16 (max ~32768) - // ln(32768) ≈ 10.4 + // exp overflows I16F16 above ln(32768) ≈ 10.4. if x_f64 > 10.0 || x_f64 < -20.0 { continue; } @@ -248,7 +215,6 @@ mod reference_comparison { let result = exp(x); let expected = x_f64.exp(); - // Use relative error since exp spans many orders of magnitude let err = if expected > 0.01 { (result.to_num::() - expected).abs() / expected } else { @@ -267,13 +233,12 @@ mod reference_comparison { fn ln_vs_f64() { for i in 0..SAMPLES { let bits = sample_bits(SEED, i); - // Only positive values let bits = (bits & 0x7FFF_FFFF).max(1); let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); if x_f64 < 0.001 { - continue; // ln near 0 is extreme + continue; } if let Ok(result) = ln(x) { @@ -330,7 +295,6 @@ mod reference_comparison { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip large values where overflow occurs if x_f64.abs() > 8.0 { continue; } @@ -339,7 +303,6 @@ mod reference_comparison { let expected_sinh = x_f64.sinh(); let expected_cosh = x_f64.cosh(); - // Use relative error for large results let sinh_err = if expected_sinh.abs() > 1.0 { (sh.to_num::() - expected_sinh).abs() / expected_sinh.abs() } else { @@ -390,7 +353,6 @@ mod reference_comparison { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // asinh involves sqrt which loses precision for large values if x_f64.abs() > 20.0 { continue; } @@ -409,10 +371,8 @@ mod reference_comparison { #[test] fn acosh_vs_f64() { - // Test over domain [1, moderate] - precision degrades for large values for i in 0..500 { - // Map to [1, 20] - let t = 1.0 + (i as f64) / 499.0 * 19.0; + let t = f64::mul_add(f64::from(i) / 499.0, 19.0, 1.0); let x = I16F16::from_num(t); let x_f64: f64 = x.to_num(); @@ -430,10 +390,8 @@ mod reference_comparison { #[test] fn atanh_vs_f64() { - // Test over domain (-1, 1) for i in 0..500 { - // Map to (-0.99, 0.99) - let t = (i as f64) / 499.0 * 1.98 - 0.99; + let t = f64::mul_add(f64::from(i) / 499.0, 1.98, -0.99); let x = I16F16::from_num(t); let x_f64: f64 = x.to_num(); @@ -452,8 +410,6 @@ mod reference_comparison { #[cfg(test)] mod identities { - //! Verify fundamental mathematical identities. - use fixed::types::I16F16; use fixed_analytics::{exp, ln, sin_cos, sinh_cosh}; @@ -461,7 +417,7 @@ mod identities { let mut x = seed.wrapping_add(index.wrapping_mul(0x9E37_79B9_7F4A_7C15)); x = x.wrapping_mul(0x5851_F42D_4C95_7F2D); x = x ^ (x >> 32); - x as i32 + u32::try_from(x & 0xFFFF_FFFF).map_or(0, u32::cast_signed) } const SAMPLES: u64 = 1000; @@ -469,7 +425,6 @@ mod identities { #[test] fn pythagorean_identity() { - // sin²(x) + cos²(x) = 1 for i in 0..SAMPLES { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); @@ -482,7 +437,7 @@ mod identities { let (s, c) = sin_cos(x); let s_f64: f64 = s.to_num(); let c_f64: f64 = c.to_num(); - let sum_sq = s_f64 * s_f64 + c_f64 * c_f64; + let sum_sq = c_f64.mul_add(c_f64, s_f64 * s_f64); assert!( (sum_sq - 1.0).abs() < 0.01, @@ -493,13 +448,11 @@ mod identities { #[test] fn hyperbolic_identity() { - // cosh²(x) - sinh²(x) = 1 for i in 0..SAMPLES { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip large values where overflow occurs if x_f64.abs() > 8.0 { continue; } @@ -507,7 +460,7 @@ mod identities { let (sh, ch) = sinh_cosh(x); let sh_f64: f64 = sh.to_num(); let ch_f64: f64 = ch.to_num(); - let diff = ch_f64 * ch_f64 - sh_f64 * sh_f64; + let diff = sh_f64.mul_add(-sh_f64, ch_f64 * ch_f64); assert!( (diff - 1.0).abs() < 0.02, @@ -518,13 +471,11 @@ mod identities { #[test] fn angle_addition() { - // sin(a + b) = sin(a)cos(b) + cos(a)sin(b) for i in 0..200 { for j in 0..200 { let a_bits = sample_bits(SEED, i); let b_bits = sample_bits(SEED ^ 0xFFFF, j); - // Use smaller values to stay in good precision range let a = I16F16::from_bits(a_bits >> 8); let b = I16F16::from_bits(b_bits >> 8); @@ -532,9 +483,11 @@ mod identities { let (sin_b, cos_b) = sin_cos(b); let (sin_ab, _) = sin_cos(a + b); - // sin(a+b) should equal sin(a)cos(b) + cos(a)sin(b) - let expected = sin_a.to_num::() * cos_b.to_num::() - + cos_a.to_num::() * sin_b.to_num::(); + let expected = f64::mul_add( + cos_a.to_num::(), + sin_b.to_num::(), + sin_a.to_num::() * cos_b.to_num::(), + ); let actual: f64 = sin_ab.to_num(); assert!( @@ -549,15 +502,11 @@ mod identities { #[test] fn double_angle() { - // sin(2x) = 2*sin(x)*cos(x) - // cos(2x) = cos²(x) - sin²(x) for i in 0..SAMPLES { let bits = sample_bits(SEED, i); - // Keep x small to avoid angle reduction errors compounding - let x = I16F16::from_bits(bits >> 8); // Much smaller range + let x = I16F16::from_bits(bits >> 8); let x_f64: f64 = x.to_num(); - // Skip large angles where errors compound if x_f64.abs() > 50.0 { continue; } @@ -568,7 +517,9 @@ mod identities { let (sin_2x, cos_2x) = sin_cos(two_x); let expected_sin_2x = 2.0 * sin_x.to_num::() * cos_x.to_num::(); - let expected_cos_2x = cos_x.to_num::().powi(2) - sin_x.to_num::().powi(2); + let expected_cos_2x = sin_x + .to_num::() + .mul_add(-sin_x.to_num::(), cos_x.to_num::().powi(2)); let sin_err = (sin_2x.to_num::() - expected_sin_2x).abs(); let cos_err = (cos_2x.to_num::() - expected_cos_2x).abs(); @@ -590,21 +541,17 @@ mod identities { #[test] fn exp_addition() { - // exp(a + b) = exp(a) * exp(b) for i in 0..200 { for j in 0..200 { let a_bits = sample_bits(SEED, i); let b_bits = sample_bits(SEED ^ 0xFFFF, j); - // Keep values small to avoid overflow - // Map to roughly [-2, 2] let a = I16F16::from_bits(a_bits >> 14); let b = I16F16::from_bits(b_bits >> 14); let a_f64: f64 = a.to_num(); let b_f64: f64 = b.to_num(); - // Skip if sum would overflow if (a_f64 + b_f64).abs() > 8.0 { continue; } @@ -628,13 +575,11 @@ mod identities { #[test] fn ln_multiplication() { - // ln(a * b) = ln(a) + ln(b) for i in 0..200 { for j in 0..200 { let a_bits = sample_bits(SEED, i); let b_bits = sample_bits(SEED ^ 0xFFFF, j); - // Keep values in range [0.5, 128] for reasonable precision let a_bits = (a_bits & 0x007F_FFFF).max(0x0000_8000); let b_bits = (b_bits & 0x007F_FFFF).max(0x0000_8000); let a = I16F16::from_bits(a_bits); @@ -643,7 +588,6 @@ mod identities { let a_f64: f64 = a.to_num(); let b_f64: f64 = b.to_num(); - // Skip if product would overflow or underflow significantly let product = a_f64 * b_f64; if product > 10000.0 || product < 0.1 { continue; @@ -667,8 +611,6 @@ mod identities { #[cfg(test)] mod roundtrips { - //! Verify that inverse functions are actually inverses. - use fixed::types::I16F16; use fixed_analytics::{ acos, acosh, asin, asinh, atan, atanh, cos, cosh, exp, ln, sin, sinh, tan, tanh, @@ -678,9 +620,8 @@ mod roundtrips { #[test] fn sin_asin_roundtrip() { - // For x in [-1, 1]: sin(asin(x)) ≈ x for i in 0..200 { - let t = (i as f64) / 199.0 * 1.98 - 0.99; // (-0.99, 0.99) + let t = f64::mul_add(f64::from(i) / 199.0, 1.98, -0.99); let x = I16F16::from_num(t); if let Ok(asin_x) = asin(x) { @@ -693,10 +634,9 @@ mod roundtrips { #[test] fn asin_sin_roundtrip() { - // For x in [-π/2, π/2]: asin(sin(x)) ≈ x let half_pi = core::f64::consts::FRAC_PI_2 * 0.95; for i in 0..200 { - let t = (i as f64) / 199.0 * 2.0 * half_pi - half_pi; + let t = f64::mul_add(f64::from(i) / 199.0 * 2.0, half_pi, -half_pi); let x = I16F16::from_num(t); let sin_x = sin(x); @@ -709,9 +649,8 @@ mod roundtrips { #[test] fn cos_acos_roundtrip() { - // For x in [-1, 1]: cos(acos(x)) ≈ x for i in 0..200 { - let t = (i as f64) / 199.0 * 1.98 - 0.99; + let t = f64::mul_add(f64::from(i) / 199.0, 1.98, -0.99); let x = I16F16::from_num(t); if let Ok(acos_x) = acos(x) { @@ -724,10 +663,9 @@ mod roundtrips { #[test] fn tan_atan_roundtrip() { - // atan(tan(x)) ≈ x for x in (-π/2, π/2) let half_pi = core::f64::consts::FRAC_PI_2 * 0.9; for i in 0..200 { - let t = (i as f64) / 199.0 * 2.0 * half_pi - half_pi; + let t = f64::mul_add(f64::from(i) / 199.0 * 2.0, half_pi, -half_pi); let x = I16F16::from_num(t); let tan_x = tan(x); let roundtrip = atan(tan_x); @@ -739,9 +677,8 @@ mod roundtrips { #[test] fn exp_ln_roundtrip() { - // For x > 0: exp(ln(x)) ≈ x for i in 1..200 { - let t = (i as f64) / 20.0; // (0.05, 10) + let t = f64::from(i) / 20.0; let x = I16F16::from_num(t); if let Ok(ln_x) = ln(x) { @@ -758,9 +695,8 @@ mod roundtrips { #[test] fn ln_exp_roundtrip() { - // ln(exp(x)) ≈ x for small x for i in 0..200 { - let t = (i as f64) / 20.0 - 5.0; // (-5, 5) + let t = f64::from(i) / 20.0 - 5.0; let x = I16F16::from_num(t); let exp_x = exp(x); @@ -773,9 +709,8 @@ mod roundtrips { #[test] fn sinh_asinh_roundtrip() { - // Hyperbolic functions have lower precision, use smaller range for i in 0..200 { - let t = (i as f64) / 50.0 - 2.0; // (-2, 2) + let t = f64::from(i) / 50.0 - 2.0; let x = I16F16::from_num(t); let sinh_x = sinh(x); let roundtrip = asinh(sinh_x); @@ -791,10 +726,9 @@ mod roundtrips { #[test] fn cosh_acosh_roundtrip() { - // For x >= 0: acosh(cosh(x)) ≈ x (cosh is not injective, so test non-negative) - // Use smaller range for better precision + // cosh is not injective, so only x ≥ 0. for i in 0..200 { - let t = (i as f64) / 50.0; // (0, 4) + let t = f64::from(i) / 50.0; let x = I16F16::from_num(t); let cosh_x = cosh(x); @@ -811,9 +745,8 @@ mod roundtrips { #[test] fn tanh_atanh_roundtrip() { - // For x in (-1, 1): tanh(atanh(x)) ≈ x for i in 0..200 { - let t = (i as f64) / 199.0 * 1.9 - 0.95; // (-0.95, 0.95) + let t = f64::mul_add(f64::from(i) / 199.0, 1.9, -0.95); let x = I16F16::from_num(t); if let Ok(atanh_x) = atanh(x) { @@ -830,19 +763,18 @@ mod roundtrips { } #[cfg(test)] -#[allow(clippy::unwrap_used, reason = "test code uses unwrap for conciseness")] mod monotonicity { - //! Verify that monotonic functions are actually monotonic. + use crate::unit::support::TestResult; use fixed::types::I16F16; use fixed_analytics::{asin, atan, exp, ln, sin, sqrt, tanh}; #[test] - fn sqrt_is_increasing() { + fn sqrt_is_increasing() -> TestResult { let mut prev = I16F16::ZERO; for i in 0..1000 { let x = I16F16::from_bits(i * 1000); - let y = sqrt(x).unwrap(); + let y = sqrt(x)?; assert!( y >= prev, "sqrt({}) = {} < sqrt({}) = {}", @@ -853,14 +785,14 @@ mod monotonicity { ); prev = y; } + Ok(()) } #[test] fn exp_is_increasing() { let mut prev = I16F16::ZERO; for i in 0..500 { - // Range from -5 to 5 - let t = (i as f64) / 50.0 - 5.0; + let t = f64::from(i) / 50.0 - 5.0; let x = I16F16::from_num(t); let y = exp(x); assert!( @@ -894,11 +826,10 @@ mod monotonicity { #[test] fn sin_increasing_on_neg_half_pi_to_half_pi() { - // sin is increasing on [-π/2, π/2] let half_pi = core::f64::consts::FRAC_PI_2; let mut prev = I16F16::MIN; for i in 0..200 { - let t = (i as f64) / 199.0 * 2.0 * half_pi - half_pi; + let t = f64::mul_add(f64::from(i) / 199.0 * 2.0, half_pi, -half_pi); let x = I16F16::from_num(t); let y = sin(x); assert!( @@ -915,8 +846,7 @@ mod monotonicity { fn atan_is_increasing() { let mut prev = I16F16::MIN; for i in 0..500 { - // Range from -100 to 100 - let t = (i as f64) / 2.5 - 100.0; + let t = f64::from(i) / 2.5 - 100.0; let x = I16F16::from_num(t); let y = atan(x); assert!( @@ -933,7 +863,7 @@ mod monotonicity { fn asin_is_increasing() { let mut prev = I16F16::MIN; for i in 0..200 { - let t = (i as f64) / 199.0 * 1.98 - 0.99; // (-0.99, 0.99) + let t = f64::mul_add(f64::from(i) / 199.0, 1.98, -0.99); let x = I16F16::from_num(t); if let Ok(y) = asin(x) { assert!( @@ -951,8 +881,7 @@ mod monotonicity { fn tanh_is_increasing() { let mut prev = I16F16::MIN; for i in 0..500 { - // Use smaller range to avoid saturation at ±1 - let t = (i as f64) / 50.0 - 5.0; // (-5, 5) + let t = f64::from(i) / 50.0 - 5.0; let x = I16F16::from_num(t); let y = tanh(x); assert!( @@ -968,8 +897,6 @@ mod monotonicity { #[cfg(test)] mod bounds { - //! Verify output bounds for functions with known ranges. - use fixed::types::I16F16; use fixed_analytics::{atan, atan2, cos, exp, sin, sqrt, tanh}; @@ -977,7 +904,7 @@ mod bounds { let mut x = seed.wrapping_add(index.wrapping_mul(0x9E37_79B9_7F4A_7C15)); x = x.wrapping_mul(0x5851_F42D_4C95_7F2D); x = x ^ (x >> 32); - x as i32 + u32::try_from(x & 0xFFFF_FFFF).map_or(0, u32::cast_signed) } const SAMPLES: u64 = 2000; @@ -985,7 +912,6 @@ mod bounds { #[test] fn sin_in_bounds() { - // sin(x) ∈ [-1, 1] for i in 0..SAMPLES { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); @@ -993,7 +919,7 @@ mod bounds { let y_f64: f64 = y.to_num(); assert!( - y_f64 >= -1.01 && y_f64 <= 1.01, + (-1.01..=1.01).contains(&y_f64), "sin({}) = {} out of bounds", x.to_num::(), y_f64 @@ -1003,7 +929,6 @@ mod bounds { #[test] fn cos_in_bounds() { - // cos(x) ∈ [-1, 1] for i in 0..SAMPLES { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); @@ -1011,7 +936,7 @@ mod bounds { let y_f64: f64 = y.to_num(); assert!( - y_f64 >= -1.01 && y_f64 <= 1.01, + (-1.01..=1.01).contains(&y_f64), "cos({}) = {} out of bounds", x.to_num::(), y_f64 @@ -1021,7 +946,6 @@ mod bounds { #[test] fn tanh_in_bounds() { - // tanh(x) ∈ (-1, 1) for i in 0..SAMPLES { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); @@ -1043,7 +967,6 @@ mod bounds { let bits = sample_bits(SEED, i); let x = I16F16::from_bits(bits); - // sqrt returns Err for negative inputs, Ok for non-negative match sqrt(x) { Ok(y) => { assert!( @@ -1071,7 +994,6 @@ mod bounds { let x = I16F16::from_bits(bits); let x_f64: f64 = x.to_num(); - // Skip extreme values if x_f64 > 10.0 || x_f64 < -20.0 { continue; } @@ -1088,7 +1010,6 @@ mod bounds { #[test] fn atan_in_bounds() { - // atan(x) ∈ (-π/2, π/2) let half_pi = core::f64::consts::FRAC_PI_2; for i in 0..SAMPLES { let bits = sample_bits(SEED, i); @@ -1107,8 +1028,6 @@ mod bounds { #[test] fn atan2_in_bounds() { - // atan2(y, x) ∈ [-π, π] - // Allow slightly wider tolerance for extreme input values let pi = core::f64::consts::PI; for i in 0..SAMPLES { let y_bits = sample_bits(SEED, i); @@ -1116,7 +1035,6 @@ mod bounds { let y = I16F16::from_bits(y_bits); let x = I16F16::from_bits(x_bits); - // Skip extreme values near I16F16::MAX that cause precision issues if y.to_num::().abs() > 30000.0 || x.to_num::().abs() > 30000.0 { continue; } From 11f0b50fa96d46f66fb3613054f595f6c78da26b Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 17:40:56 +0000 Subject: [PATCH 4/9] CI jobs set RUSTUP_TOOLCHAIN so the committed toolchain pin does not override the toolchain each job tests --- .github/workflows/ci.yml | 20 ++++++++++++++++++++ .github/workflows/publish.yml | 2 ++ 2 files changed, 22 insertions(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 3eee593..d583133 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -11,11 +11,15 @@ on: env: CARGO_TERM_COLOR: always RUSTFLAGS: -Dwarnings +# Each job sets RUSTUP_TOOLCHAIN to the toolchain it tests, since it overrides +# rust-toolchain.toml, which pins the toolchain for check.sh. jobs: test: name: Test runs-on: ${{ matrix.os }} + env: + RUSTUP_TOOLCHAIN: ${{ matrix.rust }} strategy: fail-fast: false matrix: @@ -38,6 +42,8 @@ jobs: no-std: name: no_std runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable @@ -46,6 +52,8 @@ jobs: lint: name: Lint runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable @@ -61,6 +69,8 @@ jobs: msrv: name: MSRV runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: '1.95' steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@1.95 @@ -70,6 +80,8 @@ jobs: coverage: name: Coverage runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable @@ -85,6 +97,8 @@ jobs: security: name: Security runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: rustsec/audit-check@v1.4.1 @@ -94,6 +108,8 @@ jobs: no-panic: name: Verify No Panics runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable @@ -117,6 +133,8 @@ jobs: bench-check: name: Benchmarks Compile runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable @@ -126,6 +144,8 @@ jobs: accuracy: name: Accuracy Gate runs-on: ubuntu-latest + env: + RUSTUP_TOOLCHAIN: stable steps: - uses: actions/checkout@v4 diff --git a/.github/workflows/publish.yml b/.github/workflows/publish.yml index bcab3e2..bc901ce 100644 --- a/.github/workflows/publish.yml +++ b/.github/workflows/publish.yml @@ -20,6 +20,8 @@ on: env: CARGO_TERM_COLOR: always + # Overrides rust-toolchain.toml, which pins the toolchain for check.sh. + RUSTUP_TOOLCHAIN: stable concurrency: group: publish From ea0732dceb28e26d11f03785d900afa52dda677f Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 20:56:36 +0000 Subject: [PATCH 5/9] Package: leave out the four derivation pages, which docs.rs doesn't show --- Cargo.toml | 3 +++ 1 file changed, 3 insertions(+) diff --git a/Cargo.toml b/Cargo.toml index 793984f..b0470db 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -11,6 +11,9 @@ readme = "README.md" license = "MIT" keywords = ["fixed-point", "trigonometry", "math", "no_std", "deterministic"] categories = ["mathematics", "no-std", "algorithms", "embedded"] +# Derivation notes for maintainers, which the code comments cite; docs.rs +# shows only the pages lib.rs and the modules include. +exclude = ["docs/exp.md", "docs/sin_cos.md", "docs/sinh_cosh.md", "docs/sqrt.md"] [package.metadata.docs.rs] features = ["std"] From d37874863adfb6c7e55d1ddfd6062132cbf4e343 Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 21:05:01 +0000 Subject: [PATCH 6/9] Tests use the crate's Error as their error type, so they build without the std feature, as the publish workflow's cargo test --no-default-features requires. --- tests/unit/support.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/unit/support.rs b/tests/unit/support.rs index 3202f9a..14534bb 100644 --- a/tests/unit/support.rs +++ b/tests/unit/support.rs @@ -1,4 +1,4 @@ -pub type TestResult = Result<(), Box>; +pub type TestResult = Result<(), fixed_analytics::Error>; pub struct Lcg(pub u64); From 0ed564af7e0bf5bf780286a0a2dbdbde38231cb2 Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 21:05:01 +0000 Subject: [PATCH 7/9] CI runs the tests without default features too, which the publish workflow requires before a release. --- .github/workflows/ci.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index d583133..76aa271 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -48,6 +48,7 @@ jobs: - uses: actions/checkout@v4 - uses: dtolnay/rust-toolchain@stable - run: cargo check --no-default-features + - run: cargo test --no-default-features lint: name: Lint From 969d73fe3bd788ceeb6b10d5313ce7a4794bcf6b Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 21:19:01 +0000 Subject: [PATCH 8/9] Cargo.lock moves crossbeam-epoch, a dev-dependency through criterion, off the yanked 0.9.18. --- Cargo.lock | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Cargo.lock b/Cargo.lock index 77f975d..c1d8b83 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -177,9 +177,9 @@ dependencies = [ [[package]] name = "crossbeam-epoch" -version = "0.9.18" +version = "0.9.21" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "5b82ac4a3c2ca9c3460964f020e1402edd5753411d7737aa39c3714ad1b5420e" +checksum = "dc74980687109a3b14c72fd458107bf0baa1da1a1a805e178d15501ba9b86d9d" dependencies = [ "crossbeam-utils", ] From 77a0ce7e92d4692ff44b47423cbcfd6ae3a541be Mon Sep 17 00:00:00 2001 From: GeEom Date: Sun, 25 Oct 2026 22:18:25 +0000 Subject: [PATCH 9/9] v3.1.1: the package leaves out the toolchain and CI files; the release step moves to action-gh-release v2 --- .github/workflows/publish.yml | 2 +- Cargo.lock | 2 +- Cargo.toml | 13 +++++++++---- README.md | 4 ++-- 4 files changed, 13 insertions(+), 8 deletions(-) diff --git a/.github/workflows/publish.yml b/.github/workflows/publish.yml index bc901ce..702b5a3 100644 --- a/.github/workflows/publish.yml +++ b/.github/workflows/publish.yml @@ -101,7 +101,7 @@ jobs: CARGO_REGISTRY_TOKEN: ${{ secrets.CARGO_REGISTRY_TOKEN }} - name: Create GitHub release - uses: softprops/action-gh-release@v1 + uses: softprops/action-gh-release@v2 with: tag_name: v${{ needs.check.outputs.version }} generate_release_notes: true diff --git a/Cargo.lock b/Cargo.lock index c1d8b83..d09fe43 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -222,7 +222,7 @@ dependencies = [ [[package]] name = "fixed_analytics" -version = "3.1.0" +version = "3.1.1" dependencies = [ "criterion", "fixed", diff --git a/Cargo.toml b/Cargo.toml index b0470db..8be539c 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "fixed_analytics" -version = "3.1.0" +version = "3.1.1" edition = "2024" rust-version = "1.95" authors = ["David Gathercole"] @@ -11,9 +11,14 @@ readme = "README.md" license = "MIT" keywords = ["fixed-point", "trigonometry", "math", "no_std", "deterministic"] categories = ["mathematics", "no-std", "algorithms", "embedded"] -# Derivation notes for maintainers, which the code comments cite; docs.rs -# shows only the pages lib.rs and the modules include. -exclude = ["docs/exp.md", "docs/sin_cos.md", "docs/sinh_cosh.md", "docs/sqrt.md"] +# Derivation notes for maintainers, which the code comments cite (docs.rs +# shows only the pages lib.rs and the modules include), and the local gate +# and CI configuration, which would pin anyone building the packaged crate to +# this repo's toolchain. +exclude = [ + "docs/exp.md", "docs/sin_cos.md", "docs/sinh_cosh.md", "docs/sqrt.md", + "rust-toolchain.toml", "clippy.toml", "codecov.yml", ".github/", +] [package.metadata.docs.rs] features = ["std"] diff --git a/README.md b/README.md index d71f1f5..8a6b5e1 100644 --- a/README.md +++ b/README.md @@ -30,14 +30,14 @@ Requires Rust 1.95 or later. ```toml [dependencies] -fixed_analytics = "3.1.0" +fixed_analytics = "3.1.1" ``` For `no_std` environments: ```toml [dependencies] -fixed_analytics = { version = "3.1.0", default-features = false } +fixed_analytics = { version = "3.1.1", default-features = false } ``` ## Available Functions