diff --git a/.github/workflows/unittests-rust.yml b/.github/workflows/unittests-rust.yml index 253d59957..65bc707c9 100644 --- a/.github/workflows/unittests-rust.yml +++ b/.github/workflows/unittests-rust.yml @@ -26,6 +26,8 @@ jobs: - name: Run clippy run: | poe clippy + - name: Install gsl + run: sudo apt-get install libgsl0-dev - name: Run Rust unit tests run: | poe rtest diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 47f1d172d..ed70a0acc 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -47,7 +47,7 @@ repos: - id: fmt name: cargo fmt description: Format Rust files with cargo fmt. - entry: cargo fmt -- + entry: cargo fmt --all -- language: system files: ^crates/.*\.rs$ args: [] diff --git a/Cargo.lock b/Cargo.lock index b29daaaa5..d43fd3b0b 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -3,20 +3,41 @@ version = 3 [[package]] -name = "adler" -version = "1.0.2" +name = "GSL" +version = "7.0.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "f26201604c87b1e01bd3d98f8d5d9a8fcbb815e8cedb41ffccbeb4bf593a35fe" +checksum = "db3943d5a15b5c46e991124abee6a1bc89c7c9ffb25dbb8aeb4eab926fd9b307" +dependencies = [ + "GSL-sys", + "paste", +] + +[[package]] +name = "GSL-sys" +version = "3.0.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "4577670dcc0720995dc39f04c438595eaae8ccc27f4aafd3e572dd408d01bd9d" +dependencies = [ + "libc", + "pkg-config", +] + +[[package]] +name = "adler2" +version = "2.0.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "512761e0bb2578dd7380c6baaa0f4ce03e84f95e960231d1dec8bf4d7d6e2627" [[package]] name = "ahash" -version = "0.8.3" +version = "0.8.11" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2c99f64d1e06488f620f932677e24bc6e2897582980441ae90a671415bd7ec2f" +checksum = "e89da841a80418a9b391ebaea17f5c112ffaaa96f621d2c285b5174da76b9011" dependencies = [ "cfg-if", "once_cell", "version_check", + "zerocopy", ] [[package]] @@ -30,15 +51,15 @@ dependencies = [ [[package]] name = "allocator-api2" -version = "0.2.16" +version = "0.2.21" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "0942ffc6dcaadf03badf6e6a2d0228460359d5e34b57ccdc720b7382dfbd5ec5" +checksum = "683d7910e743518b0e34f1186f92494becacb047c7b6bf616c96772180fef923" [[package]] name = "anstyle" -version = "1.0.8" +version = "1.0.10" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "1bec1de6f59aedf83baf9ff929c98f2ad654b97c9510f4e70cf6f661d49fd5b1" +checksum = "55cc3b69f167a1ef2e161439aa98aed94e6028e5f9a59be9a6ffb47aef1651f9" [[package]] name = "arraydeque" @@ -63,15 +84,15 @@ dependencies = [ [[package]] name = "autocfg" -version = "1.1.0" +version = "1.4.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "d468802bab17cbc0cc575e9b053f41e72aa36bfa6b7f55e3529ffa43161b97fa" +checksum = "ace50bade8e6234aa140d9a2f552bbee1db4d353f69b8217bc503490fc1a9f26" [[package]] name = "bitflags" -version = "2.6.0" +version = "2.8.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "b048fb63fd8b5923fc5aa7b340d8e156aec7ec02f0c78fa8a6ddc2613f6f71de" +checksum = "8f68f53c83ab957f72c32642f3868eec03eb974d1fb82e453128456482613d36" [[package]] name = "block-buffer" @@ -84,9 +105,9 @@ dependencies = [ [[package]] name = "bstr" -version = "1.10.0" +version = "1.11.3" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "40723b8fb387abc38f4f4a37c09073622e41dd12327033091ef8950659e6dc0c" +checksum = "531a9155a481e2ee699d4f98f43c0ca4ff8ee1bfd55c31e9e98fb29d2b176fe0" dependencies = [ "memchr", "serde", @@ -106,9 +127,9 @@ checksum = "baf1de4339761588bc0619e3cbc0120ee582ebb74b53b4efbf79117bd2da40fd" [[package]] name = "cpufeatures" -version = "0.2.13" +version = "0.2.17" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "51e852e6dc9a5bed1fae92dd2375037bf2b768725bf3be87811edee3249d09ad" +checksum = "59ed5838eebb26a2bb2e58f6d5b5316989ae9d08bab10e0e6d103e656d1b0280" dependencies = [ "libc", ] @@ -124,9 +145,9 @@ dependencies = [ [[package]] name = "crossbeam-deque" -version = "0.8.5" +version = "0.8.6" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "613f8cc01fe9cf1a3eb3d7f488fd2fa8388403e97039e2f73692932e291a770d" +checksum = "9dd111b7b7f7d55b72c0a6ae361660ee5853c9af73f70c3c2ef6858b950e2e51" dependencies = [ "crossbeam-epoch", "crossbeam-utils", @@ -143,9 +164,9 @@ dependencies = [ [[package]] name = "crossbeam-utils" -version = "0.8.20" +version = "0.8.21" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "22ec99545bb0ed0ea7bb9b8e1e9122ea386ff8a48c0922e43f36d45ab09e0e80" +checksum = "d0a5c400df2834b80a4c3327b3aad3a4c4cd4de0629063962b03235697506a28" [[package]] name = "crypto-common" @@ -167,7 +188,7 @@ dependencies = [ "ndarray-npy", "predicates", "tar", - "thiserror", + "thiserror 1.0.69", "yaml-rust2", ] @@ -205,52 +226,53 @@ dependencies = [ name = "ekore" version = "0.0.1" dependencies = [ - "float-cmp", + "GSL", + "float-cmp 0.9.0", "num", ] [[package]] name = "encoding_rs" -version = "0.8.34" +version = "0.8.35" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "b45de904aa0b010bce2ab45264d0631681847fa7b6f2eaa7dab7619943bc4f59" +checksum = "75030f3c4f45dafd7586dd6780965a8c7e8e285a5ecb86713e63a79c5b2766f3" dependencies = [ "cfg-if", ] [[package]] name = "errno" -version = "0.3.9" +version = "0.3.10" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "534c5cf6194dfab3db3242765c03bbe257cf92f22b38f6bc0c58d59108a820ba" +checksum = "33d852cb9b869c2a9b3df2f71a3074817f01e1844f839a144f5fcef059a4eb5d" dependencies = [ "libc", - "windows-sys 0.52.0", + "windows-sys", ] [[package]] name = "fastrand" -version = "2.1.0" +version = "2.3.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "9fc0510504f03c51ada170672ac806f1f105a88aa97a5281117e1ddc3368e51a" +checksum = "37909eebbb50d72f9059c3b6d82c0463f2ff062c9e95845c43a6c9c0355411be" [[package]] name = "filetime" -version = "0.2.24" +version = "0.2.25" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "bf401df4a4e3872c4fe8151134cf483738e74b67fc934d6532c882b3d24a4550" +checksum = "35c0522e981e68cbfa8c3f978441a5f34b30b96e146b33cd3359176b50fe8586" dependencies = [ "cfg-if", "libc", "libredox", - "windows-sys 0.59.0", + "windows-sys", ] [[package]] name = "flate2" -version = "1.0.31" +version = "1.1.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "7f211bbe8e69bbd0cfdea405084f128ae8b4aaa6b0b522fc8f2b009084797920" +checksum = "11faaf5a5236997af9848be0bef4db95824b1d534ebc64d0f0c6cf3e67bd38dc" dependencies = [ "crc32fast", "miniz_oxide", @@ -265,6 +287,15 @@ dependencies = [ "num-traits", ] +[[package]] +name = "float-cmp" +version = "0.10.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b09cf3155332e944990140d967ff5eceb70df778b34f77d8075db46e4704e6d8" +dependencies = [ + "num-traits", +] + [[package]] name = "generic-array" version = "0.14.7" @@ -275,11 +306,23 @@ dependencies = [ "version_check", ] +[[package]] +name = "getrandom" +version = "0.3.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "43a49c392881ce6d5c3b8cb70f98717b7c07aabbdff06687b9030dbfbe2725f8" +dependencies = [ + "cfg-if", + "libc", + "wasi", + "windows-targets", +] + [[package]] name = "globset" -version = "0.4.14" +version = "0.4.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "57da3b9b5b85bd66f31093f8c408b90a74431672542466497dcbdfdc02034be1" +checksum = "15f1ce686646e7f1e19bf7d5533fe443a45dbfb990e00629110797578b42fb19" dependencies = [ "aho-corasick", "bstr", @@ -301,9 +344,9 @@ dependencies = [ [[package]] name = "hashbrown" -version = "0.14.0" +version = "0.14.5" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2c6201b9ff9fd90a5a3bac2e56a830d0caa509576f0e503818ee82c181b3437a" +checksum = "e5274423e17b7c9fc20b6e7e208532f9b19825d82dfd615708b70edd83df41f1" dependencies = [ "ahash", "allocator-api2", @@ -320,9 +363,9 @@ dependencies = [ [[package]] name = "ignore" -version = "0.4.22" +version = "0.4.23" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "b46810df39e66e925525d6e38ce1e7f6e1d208f72dc39757880fcb66e2c58af1" +checksum = "6d89fd380afde86567dfba715db065673989d6253f42b88179abd3eae47bda4b" dependencies = [ "crossbeam-deque", "globset", @@ -336,9 +379,9 @@ dependencies = [ [[package]] name = "libc" -version = "0.2.155" +version = "0.2.170" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "97b3888a4aecf77e811145cadf6eef5901f4782c53886191b2f693f24761847c" +checksum = "875b3680cb2f8f71bdcf9a30f38d48282f5d3c95cbf9b3fa57269bb5d5c06828" [[package]] name = "libredox" @@ -353,15 +396,15 @@ dependencies = [ [[package]] name = "linux-raw-sys" -version = "0.4.14" +version = "0.4.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "78b3ae25bc7c8c38cec158d1f2757ee79e9b3740fbc7ccf0e59e4b08d793fa89" +checksum = "d26c52dbd32dccf2d10cac7725f8eae5296885fb5703b261f7d0a0739ec807ab" [[package]] name = "log" -version = "0.4.22" +version = "0.4.26" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a7a70ba024b9dc04c27ea2f0c0548feb474ec5c54bba33a7f72f873a39d07b24" +checksum = "30bde2b3dc3671ae49d8e2e9f044c7c005836e7a023ee57cffa25ab82764bb9e" [[package]] name = "lz4_flex" @@ -390,11 +433,11 @@ checksum = "78ca9ab1a0babb1e7d5695e3530886289c18cf2f87ec19a575a0abdce112e3a3" [[package]] name = "miniz_oxide" -version = "0.7.4" +version = "0.8.5" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "b8a240ddb74feaf34a79a7add65a741f3167852fba007066dcac1ca548d89c08" +checksum = "8e3e04debbb59698c15bacbb6d93584a8c0ca9cc3213cb423d31f760d8843ce5" dependencies = [ - "adler", + "adler2", ] [[package]] @@ -432,9 +475,9 @@ checksum = "61807f77802ff30975e01f4f071c8ba10c022052f98b3294119f3e615d13e5be" [[package]] name = "num" -version = "0.4.1" +version = "0.4.3" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "b05180d69e3da0e530ba2a1dae5110317e49e3b7f3d41be227dc5f92e49ee7af" +checksum = "35bd024e8b2ff75562e5f34e7f4905839deb4b22955ef5e73d2fea1b9813cb23" dependencies = [ "num-bigint", "num-complex", @@ -446,39 +489,37 @@ dependencies = [ [[package]] name = "num-bigint" -version = "0.4.3" +version = "0.4.6" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "f93ab6289c7b344a8a9f60f88d80aa20032336fe78da341afc91c8a2341fc75f" +checksum = "a5e44f723f1133c9deac646763579fdb3ac745e418f2a7af9cd0c431da1f20b9" dependencies = [ - "autocfg", "num-integer", "num-traits", ] [[package]] name = "num-complex" -version = "0.4.3" +version = "0.4.6" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "02e0d21255c828d6f128a1e41534206671e8c3ea0c62f32291e808dc82cff17d" +checksum = "73f88a1307638156682bada9d7604135552957b7818057dcef22705b4d509495" dependencies = [ "num-traits", ] [[package]] name = "num-integer" -version = "0.1.45" +version = "0.1.46" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "225d3389fb3509a24c93f5c29eb6bde2586b98d9f016636dff58d7c6f7569cd9" +checksum = "7969661fd2958a5cb096e56c8e1ad0444ac2bbcd0061bd28660485a44879858f" dependencies = [ - "autocfg", "num-traits", ] [[package]] name = "num-iter" -version = "0.1.43" +version = "0.1.45" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "7d03e6c028c5dc5cac6e2dec0efda81fc887605bb3d884578bb6d6bf7514e252" +checksum = "1429034a0490724d0075ebb2bc9e875d6503c3cf69e235a8941aa757d83ef5bf" dependencies = [ "autocfg", "num-integer", @@ -487,11 +528,10 @@ dependencies = [ [[package]] name = "num-rational" -version = "0.4.1" +version = "0.4.2" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "0638a1c9d0a3c0914158145bc76cff373a75a627e6ecbfb71cbe6f453a5a19b0" +checksum = "f83d14da390562dca69fc84082e73e548e1ad308d24accdedd2720017cb37824" dependencies = [ - "autocfg", "num-bigint", "num-integer", "num-traits", @@ -499,35 +539,41 @@ dependencies = [ [[package]] name = "num-traits" -version = "0.2.16" +version = "0.2.19" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "f30b0abd723be7e2ffca1272140fac1a2f084c77ec3e123c192b66af1ee9e6c2" +checksum = "071dfc062690e90b734c0b2273ce72ad0ffa95f0c74596bc250dcfd960262841" dependencies = [ "autocfg", ] [[package]] name = "once_cell" -version = "1.18.0" +version = "1.20.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "945462a4b81e43c4e3ba96bd7b49d834c6f61198356aa858733bc4acf3cbe62e" + +[[package]] +name = "paste" +version = "1.0.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "dd8b5dd2ae5ed71462c540258bedcb51965123ad7e7ccf4b9a8cafaa4a63576d" +checksum = "57c0d7b74b563b49d38dae00a0c37d4d6de9b432382b2892f0574ddcae73fd0a" [[package]] name = "pest" -version = "2.7.11" +version = "2.7.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "cd53dff83f26735fdc1ca837098ccf133605d794cdae66acfc2bfac3ec809d95" +checksum = "8b7cafe60d6cf8e62e1b9b2ea516a089c008945bb5a275416789e7db0bc199dc" dependencies = [ "memchr", - "thiserror", + "thiserror 2.0.11", "ucd-trie", ] [[package]] name = "pest_derive" -version = "2.7.11" +version = "2.7.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2a548d2beca6773b1c244554d36fcf8548a8a58e74156968211567250e48e49a" +checksum = "816518421cfc6887a0d62bf441b6ffb4536fcc926395a69e1a85852d4363f57e" dependencies = [ "pest", "pest_generator", @@ -535,9 +581,9 @@ dependencies = [ [[package]] name = "pest_generator" -version = "2.7.11" +version = "2.7.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "3c93a82e8d145725dcbaf44e5ea887c8a869efdcc28706df2d08c69e17077183" +checksum = "7d1396fd3a870fc7838768d171b4616d5c91f6cc25e377b673d714567d99377b" dependencies = [ "pest", "pest_meta", @@ -548,24 +594,30 @@ dependencies = [ [[package]] name = "pest_meta" -version = "2.7.11" +version = "2.7.15" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a941429fea7e08bedec25e4f6785b6ffaacc6b755da98df5ef3e7dcf4a124c4f" +checksum = "e1e58089ea25d717bfd31fb534e4f3afcc2cc569c70de3e239778991ea3b7dea" dependencies = [ "once_cell", "pest", "sha2", ] +[[package]] +name = "pkg-config" +version = "0.3.31" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "953ec861398dccce10c670dfeaf3ec4911ca479e9c02154b3a215178c5f566f2" + [[package]] name = "predicates" -version = "3.1.2" +version = "3.1.3" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "7e9086cc7640c29a356d1a29fd134380bee9d8f79a17410aa76e7ad295f42c97" +checksum = "a5d19ee57562043d37e82899fade9a22ebab7be9cef5026b07fda9cdd4293573" dependencies = [ "anstyle", "difflib", - "float-cmp", + "float-cmp 0.10.0", "normalize-line-endings", "predicates-core", "regex", @@ -573,15 +625,15 @@ dependencies = [ [[package]] name = "predicates-core" -version = "1.0.8" +version = "1.0.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "ae8177bee8e75d6846599c6b9ff679ed51e882816914eec639944d7c9aa11931" +checksum = "727e462b119fe9c93fd0eb1429a5f7647394014cf3c04ab2c0350eeb09095ffa" [[package]] name = "predicates-tree" -version = "1.0.11" +version = "1.0.12" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "41b740d195ed3166cd147c8047ec98db0e22ec019eb8eeb76d343b795304fb13" +checksum = "72dd2d6d381dfb73a193c7fca536518d7caee39fc8503f74e7dc0be0531b425c" dependencies = [ "predicates-core", "termtree", @@ -589,9 +641,9 @@ dependencies = [ [[package]] name = "proc-macro2" -version = "1.0.86" +version = "1.0.93" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "5e719e8df665df0d1c8fbfd238015744736151d4445ec0836b8e628aae103b77" +checksum = "60946a68e5f9d28b0dc1c21bb8a97ee7d018a8b322fa57838ba31cc878e22d99" dependencies = [ "unicode-ident", ] @@ -611,9 +663,9 @@ dependencies = [ [[package]] name = "quote" -version = "1.0.36" +version = "1.0.38" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "0fa76aaf39101c457836aec0ce2316dbdc3ab723cdda1c6bd4e6ad4208acaca7" +checksum = "0e4dccaaaf89514f546c693ddc140f729f958c247918a13380cccc6078391acc" dependencies = [ "proc-macro2", ] @@ -626,18 +678,18 @@ checksum = "60a357793950651c4ed0f3f52338f53b2f809f32d83a07f72909fa13e4c6c1e3" [[package]] name = "redox_syscall" -version = "0.5.3" +version = "0.5.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2a908a6e00f1fdd0dfd9c0eb08ce85126f6d8bbda50017e74bc4a4b7d4a926a4" +checksum = "82b568323e98e49e2a0899dcee453dd679fae22d69adf9b11dd508d1549b7e2f" dependencies = [ "bitflags", ] [[package]] name = "regex" -version = "1.10.6" +version = "1.11.1" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "4219d74c6b67a3654a9fbebc4b419e22126d13d2f3c4a07ee0cb61ff79a79619" +checksum = "b544ef1b4eac5dc2db33ea63606ae9ffcfac26c1416a2806ae0bf5f56b201191" dependencies = [ "aho-corasick", "memchr", @@ -647,9 +699,9 @@ dependencies = [ [[package]] name = "regex-automata" -version = "0.4.7" +version = "0.4.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "38caf58cc5ef2fed281f89292ef23f6365465ed9a41b7a7754eb4e26496c92df" +checksum = "809e8dc61f6de73b46c85f4c96486310fe304c434cfa43669d7b40f711150908" dependencies = [ "aho-corasick", "memchr", @@ -658,21 +710,21 @@ dependencies = [ [[package]] name = "regex-syntax" -version = "0.8.4" +version = "0.8.5" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "7a66a03ae7c801facd77a29370b4faec201768915ac14a721ba36f20bc9c209b" +checksum = "2b15c43186be67a4fd63bee50d0303afffcef381492ebe2c5d87f324e1b8815c" [[package]] name = "rustix" -version = "0.38.34" +version = "0.38.44" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "70dc5ec042f7a43c4a73241207cecc9873a06d45debb38b329f8541d85c2730f" +checksum = "fdb5bc1ae2baa591800df16c9ca78619bf65c0488b41b96ccec5d11220d8c154" dependencies = [ "bitflags", "errno", "libc", "linux-raw-sys", - "windows-sys 0.52.0", + "windows-sys", ] [[package]] @@ -686,18 +738,18 @@ dependencies = [ [[package]] name = "serde" -version = "1.0.208" +version = "1.0.218" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "cff085d2cb684faa248efb494c39b68e522822ac0de72ccf08109abde717cfb2" +checksum = "e8dfc9d19bdbf6d17e22319da49161d5d0108e4188e8b680aef6299eed22df60" dependencies = [ "serde_derive", ] [[package]] name = "serde_derive" -version = "1.0.208" +version = "1.0.218" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "24008e81ff7613ed8e5ba0cfaf24e2c2f1e5b8a0495711e44fcd4882fca62bcf" +checksum = "f09503e191f4e797cb8aac08e9a4a4695c5edf6a2e70e376d961ddd5c969f82b" dependencies = [ "proc-macro2", "quote", @@ -723,9 +775,9 @@ checksum = "a2eb9349b6444b326872e140eb1cf5e7c522154d69e7a0ffb0fb81c06b37543f" [[package]] name = "syn" -version = "2.0.74" +version = "2.0.98" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "1fceb41e3d546d0bd83421d3409b1460cc7444cd389341a4c880fe7a042cb3d7" +checksum = "36147f1a48ae0ec2b5b3bc5b537d267457555a10dc06f3dbc8cb11ba3006d3b1" dependencies = [ "proc-macro2", "quote", @@ -734,9 +786,9 @@ dependencies = [ [[package]] name = "tar" -version = "0.4.41" +version = "0.4.44" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "cb797dad5fb5b76fcf519e702f4a589483b5ef06567f160c392832c1f5e44909" +checksum = "1d863878d212c87a19c1a610eb53bb01fe12951c0501cf5a0d65f724914a667a" dependencies = [ "filetime", "libc", @@ -745,36 +797,58 @@ dependencies = [ [[package]] name = "tempfile" -version = "3.10.1" +version = "3.17.1" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "85b77fafb263dd9d05cbeac119526425676db3784113aa9295c88498cbf8bff1" +checksum = "22e5a0acb1f3f55f65cc4a866c361b2fb2a0ff6366785ae6fbb5f85df07ba230" dependencies = [ "cfg-if", "fastrand", + "getrandom", + "once_cell", "rustix", - "windows-sys 0.52.0", + "windows-sys", ] [[package]] name = "termtree" -version = "0.4.1" +version = "0.5.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "8f50febec83f5ee1df3015341d8bd429f2d1cc62bcba7ea2076759d315084683" + +[[package]] +name = "thiserror" +version = "1.0.69" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "3369f5ac52d5eb6ab48c6b4ffdc8efbcad6b89c765749064ba298f2c68a16a76" +checksum = "b6aaf5339b578ea85b50e080feb250a3e8ae8cfcdff9a461c9ec2904bc923f52" +dependencies = [ + "thiserror-impl 1.0.69", +] [[package]] name = "thiserror" -version = "1.0.63" +version = "2.0.11" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d452f284b73e6d76dd36758a0c8684b1d5be31f92b89d07fd5822175732206fc" +dependencies = [ + "thiserror-impl 2.0.11", +] + +[[package]] +name = "thiserror-impl" +version = "1.0.69" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "c0342370b38b6a11b6cc11d6a805569958d54cfa061a29969c3b5ce2ea405724" +checksum = "4fee6c4efc90059e10f81e6d42c60a18f76588c3d74cb83a0b242a2b6c7504c1" dependencies = [ - "thiserror-impl", + "proc-macro2", + "quote", + "syn", ] [[package]] name = "thiserror-impl" -version = "1.0.63" +version = "2.0.11" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a4558b58466b9ad7ca0f102865eccc95938dca1a74a856f2b57b6629050da261" +checksum = "26afc1baea8a989337eeb52b6e72a039780ce45c3edfcc9c5b9d112feeb173c2" dependencies = [ "proc-macro2", "quote", @@ -793,27 +867,27 @@ dependencies = [ [[package]] name = "typenum" -version = "1.17.0" +version = "1.18.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "42ff0bf0c66b8238c6f3b578df37d0b7848e55df8577b3f74f92a69acceeb825" +checksum = "1dccffe3ce07af9386bfd29e80c0ab1a8205a2fc34e4bcd40364df902cfa8f3f" [[package]] name = "ucd-trie" -version = "0.1.6" +version = "0.1.7" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "ed646292ffc8188ef8ea4d1e0e0150fb15a5c2e12ad9b8fc191ae7a8a7f3c4b9" +checksum = "2896d95c02a80c6d6a5d6e953d479f5ddf2dfdb6a244441010e373ac0fb88971" [[package]] name = "unicode-ident" -version = "1.0.12" +version = "1.0.17" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "3354b9ac3fae1ff6755cb6db53683adb661634f67557942dea4facebec0fee4b" +checksum = "00e2473a93778eb0bad35909dff6a10d28e63f792f16ed15e404fca9d5eeedbe" [[package]] name = "version_check" -version = "0.9.4" +version = "0.9.5" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "49874b5167b65d7193b8aba1567f5c7d93d001cafc34600cee003eda787e483f" +checksum = "0b928f33d975fc6ad9f86c8f283853ad26bdd5b10b7f1542aa2fa15e2289105a" [[package]] name = "walkdir" @@ -826,21 +900,21 @@ dependencies = [ ] [[package]] -name = "winapi-util" -version = "0.1.9" +name = "wasi" +version = "0.13.3+wasi-0.2.2" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "cf221c93e13a30d793f7645a0e7762c55d169dbb0a49671918a2319d289b10bb" +checksum = "26816d2e1a4a36a2940b96c5296ce403917633dff8f3440e9b236ed6f6bacad2" dependencies = [ - "windows-sys 0.59.0", + "wit-bindgen-rt", ] [[package]] -name = "windows-sys" -version = "0.52.0" +name = "winapi-util" +version = "0.1.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "282be5f36a8ce781fad8c8ae18fa3f9beff57ec1b52cb3de0789201425d9a33d" +checksum = "cf221c93e13a30d793f7645a0e7762c55d169dbb0a49671918a2319d289b10bb" dependencies = [ - "windows-targets", + "windows-sys", ] [[package]] @@ -916,11 +990,20 @@ version = "0.52.6" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "589f6da84c646204747d1270a2a5661ea66ed1cced2631d546fdfb155959f9ec" +[[package]] +name = "wit-bindgen-rt" +version = "0.33.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "3268f3d866458b787f390cf61f4bbb563b922d091359f9608842999eaee3943c" +dependencies = [ + "bitflags", +] + [[package]] name = "xattr" -version = "1.3.1" +version = "1.4.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "8da84f1a25939b27f6820d92aed108f83ff920fdf11a7b19366c27c4cda81d4f" +checksum = "e105d177a3871454f754b33bb0ee637ecaaac997446375fd3e5d43a2ed00c909" dependencies = [ "libc", "linux-raw-sys", @@ -938,6 +1021,26 @@ dependencies = [ "hashlink", ] +[[package]] +name = "zerocopy" +version = "0.7.35" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1b9b4fd18abc82b8136838da5d50bae7bdea537c574d8dc1a34ed098d6c166f0" +dependencies = [ + "zerocopy-derive", +] + +[[package]] +name = "zerocopy-derive" +version = "0.7.35" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "fa4f8080344d4671fb4e831a13ad1e68092748387dfc4f55e356242fae12ce3e" +dependencies = [ + "proc-macro2", + "quote", + "syn", +] + [[package]] name = "zip" version = "0.5.13" @@ -947,5 +1050,5 @@ dependencies = [ "byteorder", "crc32fast", "flate2", - "thiserror", + "thiserror 1.0.69", ] diff --git a/crates/eko/src/lib.rs b/crates/eko/src/lib.rs index b79cc48b2..242f404b5 100644 --- a/crates/eko/src/lib.rs +++ b/crates/eko/src/lib.rs @@ -74,7 +74,7 @@ fn unravel_qed_ns(res: Vec>>, order_qcd: usize, order_qed: usiz target } -/// Intergration kernel inside quad. +/// Integration kernel inside quad. /// /// # Safety /// This is the connection from Python, so we don't know what is on the other side. @@ -97,6 +97,9 @@ pub unsafe extern "C" fn rust_quad_ker(u: f64, rargs: *mut c_void) -> f64 { re: Vec::::new(), im: Vec::::new(), }; + let n3lo_ad_variation = std::slice::from_raw_parts(args.n3lo_ad_variation, 7) + .try_into() + .unwrap(); if args.is_ome { if is_singlet { @@ -125,7 +128,13 @@ pub unsafe extern "C" fn rust_quad_ker(u: f64, rargs: *mut c_void) -> f64 { let gamma_singlet_qed = ekore::anomalous_dimensions::unpolarized::spacelike::gamma_singlet_qed; raw = unravel_qed( - gamma_singlet_qed(args.order_qcd, args.order_qed, &mut c, args.nf), + gamma_singlet_qed( + args.order_qcd, + args.order_qed, + &mut c, + args.nf, + n3lo_ad_variation, + ), args.order_qcd, args.order_qed, ); @@ -135,7 +144,12 @@ pub unsafe extern "C" fn rust_quad_ker(u: f64, rargs: *mut c_void) -> f64 { false => ekore::anomalous_dimensions::unpolarized::spacelike::gamma_singlet_qcd, }; raw = unravel( - gamma_singlet_qcd(args.order_qcd, &mut c, args.nf), + gamma_singlet_qcd( + args.order_qcd, + &mut c, + args.nf, + n3lo_ad_variation[0..4].try_into().unwrap(), + ), args.order_qcd, ); } @@ -144,14 +158,27 @@ pub unsafe extern "C" fn rust_quad_ker(u: f64, rargs: *mut c_void) -> f64 { let gamma_valence_qed = ekore::anomalous_dimensions::unpolarized::spacelike::gamma_valence_qed; raw = unravel_qed( - gamma_valence_qed(args.order_qcd, args.order_qed, &mut c, args.nf), + gamma_valence_qed( + args.order_qcd, + args.order_qed, + &mut c, + args.nf, + n3lo_ad_variation[4..7].try_into().unwrap(), + ), args.order_qcd, args.order_qed, ); } else { let gamma_ns_qed = ekore::anomalous_dimensions::unpolarized::spacelike::gamma_ns_qed; raw = unravel_qed_ns( - gamma_ns_qed(args.order_qcd, args.order_qed, args.mode0, &mut c, args.nf), + gamma_ns_qed( + args.order_qcd, + args.order_qed, + args.mode0, + &mut c, + args.nf, + n3lo_ad_variation[4..7].try_into().unwrap(), + ), args.order_qcd, args.order_qed, ); @@ -162,7 +189,13 @@ pub unsafe extern "C" fn rust_quad_ker(u: f64, rargs: *mut c_void) -> f64 { true => ekore::anomalous_dimensions::polarized::spacelike::gamma_ns_qcd, false => ekore::anomalous_dimensions::unpolarized::spacelike::gamma_ns_qcd, }; - let res = gamma_ns_qcd(args.order_qcd, args.mode0, &mut c, args.nf); + let res = gamma_ns_qcd( + args.order_qcd, + args.mode0, + &mut c, + args.nf, + n3lo_ad_variation[4..7].try_into().unwrap(), + ); for el in res.iter().take(args.order_qcd) { raw.re.push(el.re); raw.im.push(el.im); @@ -282,6 +315,8 @@ pub struct QuadArgs { pub a_half_x: u8, pub a_half_y: u8, pub alphaem_running: bool, + // additional param required for N3LO + pub n3lo_ad_variation: *const u8, } /// Empty placeholder function for python callback. @@ -367,5 +402,6 @@ pub unsafe extern "C" fn empty_args() -> QuadArgs { a_half_x: 0, a_half_y: 0, alphaem_running: false, + n3lo_ad_variation: [0, 0, 0, 0, 0, 0, 0].as_ptr(), } } diff --git a/crates/ekore/Cargo.toml b/crates/ekore/Cargo.toml index 9d5acd2a9..dd7a9717d 100644 --- a/crates/ekore/Cargo.toml +++ b/crates/ekore/Cargo.toml @@ -22,3 +22,4 @@ num = "0.4.1" [dev-dependencies] float-cmp = "0.9.0" +GSL = "7.0" diff --git a/crates/ekore/refs.bib b/crates/ekore/refs.bib index 13d85e8b5..667f22713 100644 --- a/crates/ekore/refs.bib +++ b/crates/ekore/refs.bib @@ -1,3 +1,5 @@ +# Item keys MUST be valid rust identifiers, i.e. they MUST NOT contain e.g. `:` (which is usually the only relevant character). +# By convention we just drop the `:` before the year. @article{Moch2004pa, author = "Moch, S. and Vermaseren, J. A. M. and Vogt, A.", title = "{The Three loop splitting functions in QCD: The @@ -156,3 +158,68 @@ @article{deFlorian2016gvk pages = "056", year = "2016" } +@article{Moch2017uml, + author = "Moch, S. and Ruijl, B. and Ueda, T. and Vermaseren, J. A. M. and Vogt, A.", + title = "{Four-Loop Non-Singlet Splitting Functions in the Planar Limit and Beyond}", + eprint = "1707.08315", + archivePrefix = "arXiv", + primaryClass = "hep-ph", + reportNumber = "DESY-17-106, NIKHEF-2017-034, LTH-1139", + doi = "10.1007/JHEP10(2017)041", + journal = "JHEP", + volume = "10", + pages = "041", + year = "2017" +} +@article{Falcioni2024qpd, + author = "Falcioni, G. and Herzog, F. and Moch, S. and Pelloni, A. and Vogt, A.", + title = "{Four-loop splitting functions in QCD \textendash{} the gluon-gluon case \textendash{}}", + eprint = "2410.08089", + archivePrefix = "arXiv", + primaryClass = "hep-ph", + reportNumber = "ZU-TH 47/24, DESY-24-144, LTH 1384", + doi = "10.1016/j.physletb.2024.139194", + journal = "Phys. Lett. B", + volume = "860", + pages = "139194", + year = "2025" +} +@article{Falcioni2024xyt, + author = "Falcioni, G. and Herzog, F. and Moch, S. and Pelloni, A. and Vogt, A.", + title = "{Four-loop splitting functions in QCD \textendash{} The quark-to-gluon case}", + eprint = "2404.09701", + archivePrefix = "arXiv", + primaryClass = "hep-ph", + reportNumber = "ZU-TH 20/24, DESY-24-053, LTH 1367", + doi = "10.1016/j.physletb.2024.138906", + journal = "Phys. Lett. B", + volume = "856", + pages = "138906", + year = "2024" +} +@article{Falcioni2023luc, + author = "Falcioni, G. and Herzog, F. and Moch, S. and Vogt, A.", + title = "{Four-loop splitting functions in QCD \textendash{} The quark-quark case}", + eprint = "2302.07593", + archivePrefix = "arXiv", + primaryClass = "hep-ph", + reportNumber = "DESY 23--022, LTH 1333", + doi = "10.1016/j.physletb.2023.137944", + journal = "Phys. Lett. B", + volume = "842", + pages = "137944", + year = "2023" +} +@article{Falcioni2023vqq, + author = "Falcioni, G. and Herzog, F. and Moch, S. and Vogt, A.", + title = "{Four-loop splitting functions in QCD \textendash{} The gluon-to-quark case}", + eprint = "2307.04158", + archivePrefix = "arXiv", + primaryClass = "hep-ph", + reportNumber = "DESY 23-096, LTH 1345", + doi = "10.1016/j.physletb.2023.138215", + journal = "Phys. Lett. B", + volume = "846", + pages = "138215", + year = "2023" +} diff --git a/crates/ekore/src/anomalous_dimensions/polarized/spacelike.rs b/crates/ekore/src/anomalous_dimensions/polarized/spacelike.rs index bb394deda..dd636b0ad 100644 --- a/crates/ekore/src/anomalous_dimensions/polarized/spacelike.rs +++ b/crates/ekore/src/anomalous_dimensions/polarized/spacelike.rs @@ -9,7 +9,13 @@ pub mod as2; // pub mod as3; /// Compute the tower of the non-singlet anomalous dimensions. -pub fn gamma_ns_qcd(order_qcd: usize, mode: u16, c: &mut Cache, nf: u8) -> Vec> { +pub fn gamma_ns_qcd( + order_qcd: usize, + mode: u16, + c: &mut Cache, + nf: u8, + _n3lo_variation: [u8; 3], +) -> Vec> { let mut gamma_ns = vec![Complex::::zero(); order_qcd]; gamma_ns[0] = as1::gamma_ns(c, nf); // NLO and beyond @@ -36,7 +42,12 @@ pub fn gamma_ns_qcd(order_qcd: usize, mode: u16, c: &mut Cache, nf: u8) -> Vec Vec<[[Complex; 2]; 2]> { +pub fn gamma_singlet_qcd( + order_qcd: usize, + c: &mut Cache, + nf: u8, + _n3lo_variation: [u8; 4], +) -> Vec<[[Complex; 2]; 2]> { let mut gamma_S = vec![ [ [Complex::::zero(), Complex::::zero()], diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike.rs index 08de553fe..ff4da18cc 100644 --- a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike.rs +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike.rs @@ -12,9 +12,20 @@ pub mod as1; pub mod as1aem1; pub mod as2; pub mod as3; +pub mod as4; /// Compute the tower of the non-singlet anomalous dimensions. -pub fn gamma_ns_qcd(order_qcd: usize, mode: u16, c: &mut Cache, nf: u8) -> Vec> { +/// +/// `n3lo_variation = (ns_p, ns_m, ns_v)` is a list indicating which variation should +/// be used. `variation = 1,2` is the upper/lower bound, while any other value +/// returns the central (averaged) value. +pub fn gamma_ns_qcd( + order_qcd: usize, + mode: u16, + c: &mut Cache, + nf: u8, + n3lo_variation: [u8; 3], +) -> Vec> { let mut gamma_ns = vec![Complex::::zero(); order_qcd]; gamma_ns[0] = as1::gamma_ns(c, nf); // NLO and beyond @@ -23,7 +34,7 @@ pub fn gamma_ns_qcd(order_qcd: usize, mode: u16, c: &mut Cache, nf: u8) -> Vec as2::gamma_nsp(c, nf), // To fill the full valence vector in NNLO we need to add gamma_ns^1 explicitly here PID_NSM | PID_NSV => as2::gamma_nsm(c, nf), - _ => panic!("Unkown non-singlet sector element"), + _ => panic!("Unknown non-singlet sector element"), }; gamma_ns[1] = gamma_ns_1 } @@ -33,15 +44,34 @@ pub fn gamma_ns_qcd(order_qcd: usize, mode: u16, c: &mut Cache, nf: u8) -> Vec as3::gamma_nsp(c, nf), PID_NSM => as3::gamma_nsm(c, nf), PID_NSV => as3::gamma_nsv(c, nf), - _ => panic!("Unkown non-singlet sector element"), + _ => panic!("Unknown non-singlet sector element"), }; gamma_ns[2] = gamma_ns_2 } + // N3LO and beyond + if order_qcd >= 4 { + let gamma_ns_3 = match mode { + PID_NSP => as4::gamma_nsp(c, nf, n3lo_variation[0]), + PID_NSM => as4::gamma_nsm(c, nf, n3lo_variation[1]), + PID_NSV => as4::gamma_nsv(c, nf, n3lo_variation[2]), + _ => panic!("Unknown non-singlet sector element"), + }; + gamma_ns[3] = gamma_ns_3 + } gamma_ns } /// Compute the tower of the singlet anomalous dimension matrices. -pub fn gamma_singlet_qcd(order_qcd: usize, c: &mut Cache, nf: u8) -> Vec<[[Complex; 2]; 2]> { +/// +/// `n3lo_variation = (gg, gq, qg, qq)` is a list indicating which variation should +/// be used. `variation = 1,2` is the upper/lower bound, while any other value +/// returns the central (averaged) value. +pub fn gamma_singlet_qcd( + order_qcd: usize, + c: &mut Cache, + nf: u8, + n3lo_variation: [u8; 4], +) -> Vec<[[Complex; 2]; 2]> { let mut gamma_S = vec![ [ [Complex::::zero(), Complex::::zero()], @@ -58,6 +88,10 @@ pub fn gamma_singlet_qcd(order_qcd: usize, c: &mut Cache, nf: u8) -> Vec<[[Compl if order_qcd >= 3 { gamma_S[2] = as3::gamma_singlet(c, nf); } + // N3LO and beyond + if order_qcd >= 4 { + gamma_S[3] = as4::gamma_singlet(c, nf, n3lo_variation); + } gamma_S } @@ -68,6 +102,7 @@ pub fn gamma_ns_qed( mode: u16, c: &mut Cache, nf: u8, + n3lo_variation: [u8; 3], ) -> Vec>> { let col = vec![Complex::::zero(); order_qed + 1]; let mut gamma_ns = vec![col; order_qcd + 1]; @@ -78,7 +113,7 @@ pub fn gamma_ns_qed( PID_NSM_U | PID_NSM_D => PID_NSM, _ => mode, }; - let gamma_qcd = gamma_ns_qcd(order_qcd, qcd_mode, c, nf); + let gamma_qcd = gamma_ns_qcd(order_qcd, qcd_mode, c, nf, n3lo_variation); for j in 0..order_qcd { gamma_ns[1 + j][0] = gamma_qcd[j]; } @@ -86,7 +121,7 @@ pub fn gamma_ns_qed( gamma_ns[0][1] = match mode { PID_NSP_U | PID_NSM_U => EU2 * aem1::gamma_ns(c, nf), PID_NSP_D | PID_NSM_D => ED2 * aem1::gamma_ns(c, nf), - _ => panic!("Unkown non-singlet sector element"), + _ => panic!("Unknown non-singlet sector element"), }; if order_qed >= 2 { gamma_ns[0][2] = match mode { @@ -94,7 +129,7 @@ pub fn gamma_ns_qed( PID_NSP_D => ED2 * aem2::gamma_nspd(c, nf), PID_NSM_U => EU2 * aem2::gamma_nsmu(c, nf), PID_NSM_D => ED2 * aem2::gamma_nsmd(c, nf), - _ => panic!("Unkown non-singlet sector element"), + _ => panic!("Unknown non-singlet sector element"), }; } // QCDxQED corrections @@ -103,7 +138,7 @@ pub fn gamma_ns_qed( PID_NSP_D => ED2 * as1aem1::gamma_nsp(c, nf), PID_NSM_U => EU2 * as1aem1::gamma_nsm(c, nf), PID_NSM_D => ED2 * as1aem1::gamma_nsm(c, nf), - _ => panic!("Unkown non-singlet sector element"), + _ => panic!("Unknown non-singlet sector element"), }; gamma_ns } @@ -114,6 +149,7 @@ pub fn gamma_singlet_qed( order_qed: usize, c: &mut Cache, nf: u8, + n3lo_variation: [u8; 7], ) -> Vec; 4]; 4]>> { let col = vec![ [[ @@ -127,8 +163,14 @@ pub fn gamma_singlet_qed( let mut gamma_s = vec![col; order_qcd + 1]; // QCD corrections - let gamma_qcd_s = gamma_singlet_qcd(order_qcd, c, nf); - let gamma_qcd_nsp = gamma_ns_qcd(order_qcd, PID_NSP, c, nf); + let gamma_qcd_s = gamma_singlet_qcd(order_qcd, c, nf, n3lo_variation[0..4].try_into().unwrap()); + let gamma_qcd_nsp = gamma_ns_qcd( + order_qcd, + PID_NSP, + c, + nf, + n3lo_variation[4..7].try_into().unwrap(), + ); for j in 0..order_qcd { let s = gamma_qcd_s[j]; gamma_s[1 + j][0] = [ @@ -174,13 +216,14 @@ pub fn gamma_valence_qed( order_qed: usize, c: &mut Cache, nf: u8, + n3lo_variation: [u8; 3], ) -> Vec; 2]; 2]>> { let col = vec![[[Complex::::zero(), Complex::::zero(),]; 2]; order_qed + 1]; let mut gamma_v = vec![col; order_qcd + 1]; // QCD corrections - let gamma_qcd_nsv = gamma_ns_qcd(order_qcd, PID_NSV, c, nf); - let gamma_qcd_nsm = gamma_ns_qcd(order_qcd, PID_NSM, c, nf); + let gamma_qcd_nsv = gamma_ns_qcd(order_qcd, PID_NSV, c, nf, n3lo_variation); + let gamma_qcd_nsm = gamma_ns_qcd(order_qcd, PID_NSM, c, nf, n3lo_variation); for j in 0..order_qcd { gamma_v[1 + j][0] = [ [gamma_qcd_nsv[j], Complex::::zero()], @@ -210,16 +253,17 @@ mod tests { const NF: u8 = 3; const N: Complex = cmplx!(1., 0.); let mut c = Cache::new(N); + let n3lo_variation: [u8; 3] = [0, 0, 0]; assert_approx_eq_cmplx!( f64, - gamma_ns_qcd(3, PID_NSP, &mut c, NF)[0], + gamma_ns_qcd(3, PID_NSP, &mut c, NF, n3lo_variation)[0], cmplx!(0., 0.), epsilon = 1e-14 ); assert_approx_eq_cmplx_1d!( f64, - gamma_ns_qcd(2, PID_NSM, &mut c, NF), + gamma_ns_qcd(2, PID_NSM, &mut c, NF, n3lo_variation), [cmplx!(0., 0.); 2], 2, epsilon = 2e-6 @@ -227,7 +271,7 @@ mod tests { assert_approx_eq_cmplx_1d!( f64, - gamma_ns_qcd(3, PID_NSM, &mut c, NF), + gamma_ns_qcd(3, PID_NSM, &mut c, NF, n3lo_variation), [cmplx!(0., 0.); 3], 3, epsilon = 2e-4 @@ -235,7 +279,7 @@ mod tests { assert_approx_eq_cmplx_1d!( f64, - gamma_ns_qcd(3, PID_NSV, &mut c, NF), + gamma_ns_qcd(3, PID_NSV, &mut c, NF, n3lo_variation), [cmplx!(0., 0.); 3], 3, epsilon = 8e-4 @@ -247,12 +291,13 @@ mod tests { const NF: u8 = 3; const N: Complex = cmplx!(1., 0.); let mut c = Cache::new(N); + let n3lo_variation: [u8; 3] = [0, 0, 0]; // ns- for pid in [PID_NSM_U, PID_NSM_D] { assert_approx_eq_cmplx_2d!( f64, - gamma_ns_qed(1, 1, pid, &mut c, NF), + gamma_ns_qed(1, 1, pid, &mut c, NF, n3lo_variation), [[cmplx!(0., 0.); 2]; 2], 2, epsilon = 1e-5 @@ -264,12 +309,12 @@ mod tests { // as^0 a^0 must be trivial assert_approx_eq_cmplx!( f64, - gamma_ns_qed(1, 1, pid, &mut c, NF)[0][0], + gamma_ns_qed(1, 1, pid, &mut c, NF, n3lo_variation)[0][0], cmplx!(0., 0.) ); assert_approx_eq_cmplx!( f64, - gamma_ns_qed(1, 1, pid, &mut c, NF)[0][1], + gamma_ns_qed(1, 1, pid, &mut c, NF, n3lo_variation)[0][1], cmplx!(0., 0.), epsilon = 1e-5 ); @@ -281,8 +326,9 @@ mod tests { const NF: u8 = 3; const N: Complex = cmplx!(2., 0.); let mut c = Cache::new(N); + let n3lo_variation: [u8; 3] = [0, 0, 0]; - let g = gamma_valence_qed(3, 2, &mut c, NF); + let g = gamma_valence_qed(3, 2, &mut c, NF, n3lo_variation); // as^0 a^0 must be trivial assert_approx_eq_cmplx_2d!(f64, g[0][0], [[cmplx!(0., 0.); 2]; 2], 2); // reference from Python side @@ -303,7 +349,8 @@ mod tests { const NF: u8 = 3; const N: Complex = cmplx!(2., 0.); let mut c = Cache::new(N); - let g = gamma_singlet_qed(3, 2, &mut c, NF); + let n3lo_variation: [u8; 7] = [0, 0, 0, 0, 0, 0, 0]; + let g = gamma_singlet_qed(3, 2, &mut c, NF, n3lo_variation); // as^0 a^0 must be trivial assert_approx_eq_cmplx_2d!(f64, g[0][0], [[cmplx!(0., 0.); 4]; 4], 4); diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as3.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as3.rs index 66d61fc76..de4eac4fc 100644 --- a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as3.rs +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as3.rs @@ -1,6 +1,6 @@ //! |NNLO| |QCD|. -use ::num::complex::Complex; +use num::complex::Complex; use num::traits::Pow; use crate::cmplx; diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4.rs new file mode 100644 index 000000000..7f21e14e2 --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4.rs @@ -0,0 +1,67 @@ +//! |N3LO| |QCD|. + +mod ggg; +mod ggq; +mod gnsm; +mod gnsp; +mod gnsv; +mod gps; +mod gqg; + +use crate::harmonics::cache::Cache; +use num::complex::Complex; + +pub use ggg::gamma_gg; +pub use ggq::gamma_gq; +pub use gnsm::gamma_nsm; +pub use gnsp::gamma_nsp; +pub use gnsv::gamma_nsv; +pub use gps::gamma_ps; +pub use gqg::gamma_qg; + +/// Compute the singlet anomalous dimension matrix. +pub fn gamma_singlet(c: &mut Cache, nf: u8, variation: [u8; 4]) -> [[Complex; 2]; 2] { + let gamma_qq = gnsp::gamma_nsp(c, nf, variation[3]) + gps::gamma_ps(c, nf, variation[3]); + [ + [gamma_qq, gqg::gamma_qg(c, nf, variation[2])], + [ + ggq::gamma_gq(c, nf, variation[1]), + ggg::gamma_gg(c, nf, variation[0]), + ], + ] +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_momentum_conservation() { + let NF = 5; + let mut c = Cache::new(cmplx!(2., 0.)); + // Numbers are coming from the python implementation. + let quark_refs: [f64; 3] = [0.053441, 0.225674, -0.118792]; + let gluon_refs: [f64; 3] = [-0.0300842, 0.283004, -0.343172]; + for imod in [[0, 0, 0, 0], [1, 1, 1, 1], [2, 2, 2, 2]] { + let g_singlet = gamma_singlet(&mut c, NF, imod); + // quark conservation + assert_approx_eq_cmplx!( + f64, + g_singlet[0][0] + g_singlet[1][0], + cmplx!(quark_refs[imod[0] as usize], 0.), + rel = 2e-5 + ); + // gluon conservation + assert_approx_eq_cmplx!( + f64, + g_singlet[0][1] + g_singlet[1][1], + cmplx!(gluon_refs[imod[0] as usize], 0.), + rel = 6e-5 + ); + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggg.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggg.rs new file mode 100644 index 000000000..8f02ec240 --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggg.rs @@ -0,0 +1,213 @@ +/// Compute the singlet gluon-to-gluon anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::constants::{ZETA2, ZETA3}; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{ + lm11, lm11m1, lm11m2, lm12m1, lm12m2, lm13m1, lm13m2, lm14m1, +}; + +/// Compute the singlet gluon-to-gluon anomalous dimension. +/// +/// The routine is taken from [\[Falcioni:2024qpd\]][crate::bib::Falcioni2024qpd]. +/// +/// These are approximations for fixed `nf` = 3, 4 and 5 based on the +/// first 10 even moments together with small-x/large-x constraints. +/// The two sets providing the error estimate are called via `variation = 1` +/// and `variation = 2`. Any other value of `variation` invokes their average. +pub fn gamma_gg(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let nf_ = nf as f64; + let nf2 = nf_.pow(2); + let nf3 = nf_.pow(3); + + // The known large-x coefficients [except delta(1-x)] + let A4gluon = 40880.330 - 11714.246 * nf_ + 440.04876 * nf2 + 7.3627750 * nf3; + let mut B4gluon = 68587.64 - 18143.983 * nf_ + 423.81135 * nf2 + 9.0672154 * 0.1 * nf3; + + // The coefficient of delta(1-x), also called the virtual anomalous + // dimension. nf^0 and nf^1 are still approximate, but the error at + // nf^1 is far too small to be relevant in this context. + if variation == 1 { + B4gluon -= 0.2; + } else if variation == 2 { + B4gluon += 0.2; + } + + let Ccoeff = 8.5814120 * f64::pow(10., 4) - 1.3880515 * f64::pow(10., 4) * nf_ + + 1.3511111 * f64::pow(10., 2) * nf2; + let Dcoeff = + 5.4482808 * f64::pow(10., 4) - 4.3411337 * f64::pow(10., 3) * nf_ - 2.1333333 * 10. * nf2; + + let x1L4cff = 5.6460905 * 10. * nf_ - 3.6213992 * nf2; + let x1L3cff = 2.4755054 * f64::pow(10., 2) * nf_ - 4.0559671 * 10. * nf2 + 1.5802469 * nf3; + + // The known coefficients of 1/x*ln^a x terms, a = 3,2 + let bfkl0 = -8.3086173 * f64::pow(10., 3); + let bfkl1 = -1.0691199 * f64::pow(10., 5) - 9.9638304 * f64::pow(10., 2) * nf_; + + let x0L6cff = 1.44 * f64::pow(10., 2) - 2.7786008 * 10. * nf_ + 7.9012346 * 0.1 * nf2; + let x0L5cff = + -1.44 * f64::pow(10., 2) - 1.6208066 * f64::pow(10., 2) * nf_ + 1.4380247 * 10. * nf2; + let x0L4cff = 2.6165784 * f64::pow(10., 4) - 3.3447551 * f64::pow(10., 3) * nf_ + + 9.1522635 * 10. * nf2 + - 1.9753086 * 0.1 * nf3; + + // The resulting part of the function + let P3gg01 = bfkl0 * (-(6. / (-1. + n).powu(4))) + + bfkl1 * 2. / (-1. + n).powu(3) + + x0L6cff * 720. / n.powu(7) + + x0L5cff * -120. / n.powu(6) + + x0L4cff * 24. / n.powu(5) + + A4gluon * (-S1) + + B4gluon + + Ccoeff * lm11(c) + + Dcoeff * 1. / n + + x1L4cff * lm14m1(c) + + x1L3cff * lm13m1(c); + + // The selected approximations for nf = 3, 4, 5 + let P3ggApp1: Complex; + let P3ggApp2: Complex; + if nf == 3 { + P3ggApp1 = P3gg01 + - 421311.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 325557.0 * 1. / ((-1. + n) * n) + + 1679790.0 * (1. / (n + n.powu(2))) + - 1456863.0 * (1. / (2. + 3. * n + n.powu(2))) + + 3246307.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 2026324.0 * 2. / n.powu(3) + + 549188.0 * (-(6. / n.powu(4))) + + 8337.0 * lm11m1(c) + + 26718.0 * lm12m1(c) + - 27049.0 * lm13m2(c); + P3ggApp2 = P3gg01 + - 700113.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 2300581.0 * 1. / ((-1. + n) * n) + + 896407.0 * (1. / n - n / (2. + 3. * n + n.powu(2))) + - 162733.0 * (1. / (6. + 5. * n + n.powu(2))) + - 2661862.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 196759.0 * 2. / n.powu(3) + - 260607.0 * (-(6. / n.powu(4))) + + 84068.0 * lm11m1(c) + + 346318.0 * lm12m1(c) + + 315725.0 + * (-3. * S1.powu(2) + 6. * n * S1 * (ZETA2 - S2) + - 3. * (S2 + 2. * n * (S3 - ZETA3))) + / (3. * n.powu(2)); + } else if nf == 4 { + P3ggApp1 = P3gg01 + - 437084.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 361570.0 * 1. / ((-1. + n) * n) + + 1696070.0 * (1. / (n + n.powu(2))) + - 1457385.0 * (1. / (2. + 3. * n + n.powu(2))) + + 3195104.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 2009021.0 * 2. / n.powu(3) + + 544380.0 * (-(6. / n.powu(4))) + + 9938.0 * lm11m1(c) + + 24376.0 * lm12m1(c) + - 22143.0 * lm13m2(c); + P3ggApp2 = P3gg01 + - 706649.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 2274637.0 * 1. / ((-1. + n) * n) + + 836544.0 * (1. / n - n / (2. + 3. * n + n.powu(2))) + - 199929.0 * (1. / (6. + 5. * n + n.powu(2))) + - 2683760.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 168802.0 * 2. / n.powu(3) + - 250799.0 * (-(6. / n.powu(4))) + + 36967.0 * lm11m1(c) + + 24530.0 * lm12m1(c) + - 71470.0 * lm12m2(c); + } else if nf == 5 { + P3ggApp1 = P3gg01 + - 439426.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 293679.0 * 1. / ((-1. + n) * n) + + 1916281.0 * (1. / (n + n.powu(2))) + - 1615883.0 * (1. / (2. + 3. * n + n.powu(2))) + + 3648786.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 2166231.0 * 2. / n.powu(3) + + 594588.0 * (-(6. / n.powu(4))) + + 50406.0 * lm11m1(c) + + 24692.0 * lm12m1(c) + + 174067.0 * lm11m2(c); + P3ggApp2 = P3gg01 + - 705978.0 * (-(1. / (-1. + n).powu(2)) + 1. / n.powu(2)) + - 2192234.0 * 1. / ((-1. + n) * n) + + 1730508.0 * (1. / (2. + 3. * n + n.powu(2))) + + 353143.0 + * ((12. + 9. * n + n.powu(2)) + / (6. * n + 11. * n.powu(2) + 6. * n.powu(3) + n.powu(4))) + - 2602682.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 178960.0 * 2. / n.powu(3) + - 218133.0 * (-(6. / n.powu(4))) + + 2285.0 * lm11m1(c) + + 19295.0 * lm12m1(c) + - 13719.0 * lm12m2(c); + } else { + panic!("nf=6 is not available at N3LO"); + } + + // We return (for now) one of the two error-band representatives + // or the present best estimate, their average + let P3GGA = match variation { + 1 => P3ggApp1, + 2 => P3ggApp2, + _ => 0.5 * (P3ggApp1 + P3ggApp2), + }; + -P3GGA +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_reference_moments() { + fn gg3_moment(N: usize, nf: f64) -> f64 { + let nf2 = nf * nf; + let nf3 = nf2 * nf; + // From Eq. 5 of [Falcioni:2024qpd] + let mom_list = [ + 654.4627782205557 * nf - 245.6106197887179 * nf2 + 0.9249909688301847 * nf3, + 39876.123276008046 - 10103.4511350227 * nf + + 437.0988475397789 * nf2 + + 12.955565459350593 * nf3, + 53563.84353419538 - 14339.131035160317 * nf + + 652.7773306808972 * nf2 + + 16.654103652963503 * nf3, + 62279.7437813437 - 17150.696783851945 * nf + + 785.8806126875509 * nf2 + + 18.933103109772713 * nf3, + 68958.7532 - 19307.3854 * nf + 883.929802 * nf2 + 20.6112832 * nf3, + 74473.0024 - 21076.0320 * nf + 962.264417 * nf2 + 21.9511603 * nf3, + 79209.0111 - 22583.5268 * nf + 1027.80706 * nf2 + 23.0713754 * nf3, + 83378.4014 - 23901.3437 * nf + 1084.30677 * nf2 + 24.0362925 * nf3, + 87112.4096 - 25074.2309 * nf + 1134.04028 * nf2 + 24.8850403 * nf3, + 90499.2530 - 26132.2983 * nf + 1178.50283 * nf2 + 25.643327 * nf3, + ]; + mom_list[(N - 2) / 2] + } + for variation in [0, 1, 2] { + for NF in [3, 4, 5] { + for N in [2.0, 4.0, 6.0, 8.0, 10., 12., 14., 16., 18., 20.] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_gg(&mut c, NF, variation); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(gg3_moment(N as usize, NF as f64), 0.), + rel = 4e-4 + ); + } + } + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggq.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggq.rs new file mode 100644 index 000000000..7eedf7a8f --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/ggq.rs @@ -0,0 +1,189 @@ +/// Compute the singlet quark-to-gluon anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::constants::ZETA2; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{lm11, lm12, lm12m1, lm13, lm14, lm14m1, lm15, lm15m1}; + +/// Compute the singlet quark-to-gluon anomalous dimension. +/// +/// The routine is taken from [\[Falcioni:2024xyt\]][crate::bib::Falcioni2024xyt]. +/// +/// These are approximations for fixed `nf` = 3, 4 and 5 based on the +/// first 10 even moments together with small-x/large-x constraints. +/// The two sets providing the error estimate are called via `variation = 1` +/// and `variation = 2`. Any other value of `variation` invokes their average. +pub fn gamma_gq(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let nf_ = nf as f64; + let nf2 = nf_.pow(2); + + // Known large-x coefficients + let x1L5cff = 1.3443073 * 10. - 5.4869684 * 0.1 * nf_; + let x1L4cff = 3.7539831 * f64::pow(10., 2) - 3.4494742 * 10. * nf_ + 8.7791495 * 0.1 * nf2; + let y1L5cff = 2.2222222 * 10. - 5.4869684 * 0.1 * nf_; + let y1L4cff = 6.6242163 * f64::pow(10., 2) - 4.7992684 * 10. * nf_ + 8.7791495 * 0.1 * nf2; + + // Small-x, Casimir scaled from P_gg (approx. for bfkl1) + let bfkl0 = -8.3086173 * f64::pow(10., 3) / 2.25; + let bfkl1 = (-1.0691199 * f64::pow(10., 5) - nf_ * 9.9638304 * f64::pow(10., 2)) / 2.25; + + // Small-x double-logs with x^0 + let x0L6cff = 5.2235940 * 10. - 7.3744856 * nf_; + let x0L5cff = -2.9221399 * f64::pow(10., 2) + 1.8436214 * nf_; + let x0L4cff = + 7.3106077 * f64::pow(10., 3) - 3.7887135 * f64::pow(10., 2) * nf_ - 3.2438957 * 10. * nf2; + + // The resulting part of the function + let P3GQ01 = bfkl0 * (-(6. / (-1. + n).powu(4))) + + bfkl1 * 2. / (-1. + n).powu(3) + + x0L6cff * 720. / n.powu(7) + + x0L5cff * -120. / n.powu(6) + + x0L4cff * 24. / n.powu(5) + + x1L4cff * lm14(c) + + x1L5cff * lm15(c) + + y1L4cff * lm14m1(c) + + y1L5cff * lm15m1(c); + + // The selected approximations for nf = 3, 4, 5 + let P3gqApp1: Complex; + let P3gqApp2: Complex; + if nf == 3 { + P3gqApp1 = P3GQ01 + 6.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + - 744384.0 * 1. / ((-1. + n) * n) + + 2453640.0 * 1. / n + - 1540404.0 * (2. / (1. + n) + 1. / (2. + n)) + + 1933026.0 * -1. / n.powu(2) + + 1142069.0 * 2. / n.powu(3) + + 162196.0 * -6. / n.powu(4) + - 2172.1 * lm13(c) + - 93264.1 * lm12(c) + - 786973.0 * lm11(c) + + 875383.0 * lm12m1(c); + P3gqApp2 = + P3GQ01 + 3.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + 142414.0 * 1. / ((-1. + n) * n) + - 326525.0 * 1. / n + + 2159787.0 * ((3. + n) / (2. + 3. * n + n.powu(2))) + - 289064.0 * -1. / n.powu(2) + - 176358.0 * 2. / n.powu(3) + + 156541.0 * -6. / n.powu(4) + + 9016.5 * lm13(c) + + 136063.0 * lm12(c) + + 829482.0 * lm11(c) + - 2359050.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else if nf == 4 { + P3gqApp1 = P3GQ01 + 6.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + - 743535.0 * 1. / ((-1. + n) * n) + + 2125286.0 * 1. / n + - 1332472.0 * (2. / (1. + n) + 1. / (2. + n)) + + 1631173.0 * -1. / n.powu(2) + + 1015255.0 * 2. / n.powu(3) + + 142612.0 * -6. / n.powu(4) + - 1910.4 * lm13(c) + - 80851.0 * lm12(c) + - 680219.0 * lm11(c) + + 752733.0 * lm12m1(c); + P3gqApp2 = + P3GQ01 + 3.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + 160568.0 * 1. / ((-1. + n) * n) + - 361207.0 * 1. / n + + 2048948.0 * ((3. + n) / (2. + 3. * n + n.powu(2))) + - 245963.0 * -1. / n.powu(2) + - 171312.0 * 2. / n.powu(3) + + 163099.0 * -6. / n.powu(4) + + 8132.2 * lm13(c) + + 124425.0 * lm12(c) + + 762435.0 * lm11(c) + - 2193335.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else if nf == 5 { + P3gqApp1 = P3GQ01 + 6.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + - 785864.0 * 1. / ((-1. + n) * n) + + 285034.0 * 1. / n + - 131648.0 * (2. / (1. + n) + 1. / (2. + n)) + - 162840.0 * -1. / n.powu(2) + + 321220.0 * 2. / n.powu(3) + + 12688.0 * -6. / n.powu(4) + + 1423.4 * lm13(c) + + 1278.9 * lm12(c) + - 30919.9 * lm11(c) + + 47588.0 * lm12m1(c); + P3gqApp2 = + P3GQ01 + 3.0 * bfkl1 * (-(1. / (-1. + n).powu(2))) + 177094.0 * 1. / ((-1. + n) * n) + - 470694.0 * 1. / n + + 1348823.0 * ((3. + n) / (2. + 3. * n + n.powu(2))) + - 52985.0 * -1. / n.powu(2) + - 87354.0 * 2. / n.powu(3) + + 176885.0 * -6. / n.powu(4) + + 4748.8 * lm13(c) + + 65811.9 * lm12(c) + + 396390.0 * lm11(c) + - 1190212.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else { + panic!("nf=6 is not available at N3LO"); + } + + // We return (for now) one of the two error-band boundaries + // or the present best estimate, their average + let P3GQA = match variation { + 1 => P3gqApp1, + 2 => P3gqApp2, + _ => 0.5 * (P3gqApp1 + P3gqApp2), + }; + -P3GQA +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_reference_moments() { + fn gq3_moment(N: usize, nf: f64) -> f64 { + let nf2 = nf * nf; + let nf3 = nf2 * nf; + // From Eq. 9 of [\[Falcioni:2024xyt\]][crate::bib:Falcioni:2024xyt]. + let mom_list = [ + -16663.225488 + 4439.143749608238 * nf + - 202.55547919891168 * nf2 + - 6.375390720235586 * nf3, + -6565.7531450230645 + 1291.0067460871576 * nf + - 16.146190170051486 * nf2 + - 0.8397634037808341 * nf3, + -3937.479370556893 + 679.7185057363981 * nf + - 1.3720775271604673 * nf2 + - 0.13979432728276966 * nf3, + -2803.644107251366 + + 436.39305738710254 * nf + + 1.8149462465491055 * nf2 + + 0.07358858022119033 * nf3, + -2179.48761 + 310.063163 * nf + 2.65636842 * nf2 + 0.15719522 * nf3, + -1786.31231 + 234.383019 * nf + 2.82817592 * nf2 + 0.19211953 * nf3, + -1516.59810 + 184.745296 * nf + 2.78076831 * nf2 + 0.20536518 * nf3, + -1320.36106 + 150.076970 * nf + 2.66194730 * nf2 + 0.20798493 * nf3, + -1171.29329 + 124.717778 * nf + 2.52563073 * nf2 + 0.20512226 * nf3, + -1054.26140 + 105.497994 * nf + 2.39223358 * nf2 + 0.19938504 * nf3, + ]; + mom_list[(N - 2) / 2] + } + for variation in [0, 1, 2] { + for NF in [3, 4, 5] { + for N in [2.0, 4.0, 6.0, 8.0, 10., 12., 14., 16., 18., 20.] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_gq(&mut c, NF, variation); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(gq3_moment(N as usize, NF as f64), 0.), + rel = 4e-4 + ); + } + } + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsm.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsm.rs new file mode 100644 index 000000000..507a42af7 --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsm.rs @@ -0,0 +1,238 @@ +/// Compute the non-singlet-like non-singlet anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::cmplx; +use crate::constants::{CF, ZETA3}; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{lm11, lm11m1, lm12m1, lm13m1}; + +/// Compute the valence-like non-singlet anomalous dimension. +/// +/// The routine is taken from [\[Moch:2017uml\]][crate::bib::Moch2017uml]. +/// +/// The $n_f^{0,1}$ leading large-$N_c$ contributions and the $n_f^2$ part +/// are high-accuracy (0.1% or better) parametrizations of the exact +/// results. The $n_f^3$ expression is exact up to numerical truncations. +/// +/// The remaining $n_f^{0,1}$ terms are approximations based on the first +/// eight even moments together with small-x and large-x constraints. +/// The two sets spanning the error estimate are called via `variation = 1` +/// and ``variation = 2``. Any other value of `variation` invokes their average. +pub fn gamma_nsm(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + + // Leading large-n_c, nf^0 and nf^1, parametrized + #[rustfmt::skip] + let P3NSA0 = 360.0 / n.powu(7) + - 1920.0 / n.powu(6) + + 7147.812 / n.powu(5) + - 17179.356 / n.powu(4) + + 34241.9 / n.powu(3) + - 51671.329999999994 / n.powu(2) + + 19069.8 * lm11(c) + - (491664.8019540468 / n) + - 4533.0 / (1. + n).powu(3) + - 11825.0 / (1. + n).powu(2) + + 129203.0 / (1. + n) + - 254965.0 / (2. + n) + + 83377.5 / (3. + n) + - 45750.0 / (4. + n) + + (49150.0 * (6.803662258392675 + n) * S1) / (n.powu(2) * (1.0 + n)) + + (334400.0 * S2) / n; + #[rustfmt::skip] + let P3NSA1 = 160.0 / n.powu(6) + - 864.0 / n.powu(5) + + 2583.1848 / n.powu(4) + - 5834.624 / n.powu(3) + + 9239.374 / n.powu(2) + - 3079.76 * lm11(c) + - (114047.0 / n) + - 465.0 / (1. + n).powu(4) + - 1230.0 / (1. + n).powu(3) + + 7522.5 / (1. + n).powu(2) + + 55669.3 / (1. + n) + - 43057.8 / (2. + n) + + 13803.8 / (3. + n) + - 7896.0 / (4. + n) + - (120.0 * (-525.063 + n) * S1) / (n.powu(2) * (1.0 + n)) + + (63007.5 * S2) / n; + + // Nonleading large-n_c, nf^0 and nf^1: two approximations + #[rustfmt::skip] + let P3NMA01 = 0.4964335 * (720. / n.powu(7) - 720.0 / n.powu(6)) + - 13.5288 / n.powu(4) + + 1618.07 / n.powu(2) + - 2118.8669999999997 * lm11(c) + + 31897.8 * lm11m1(c) + + 4653.76 * lm12m1(c) + + 3902.3590000000004 / n + + 5992.88 / (1. + n) + + 19335.7 / (2. + n) + - 31321.4 / (3. + n); + #[rustfmt::skip] + let P3NMA02 = 0.4964335 * (720. / n.powu(7) - 2160.0 / n.powu(6)) + - 189.6138 / n.powu(4) + + 3065.92 / n.powu(3) + - 2118.8669999999997 * lm11(c) + - 3997.39 * lm11m1(c) + + 511.567 * lm13m1(c) + - (2099.268 / n) + + 4043.59 / (1. + n) + - 19430.190000000002 / (2. + n) + + 15386.6 / (3. + n); + + #[rustfmt::skip] + let P3NMA11 = 64.7083 / n.powu(5) + - 254.024 / n.powu(3) + + 337.931 * lm11(c) + + 1856.63 * lm11m1(c) + + 440.17 * lm12m1(c) + + 419.53485 / n + + 114.457 / (1. + n) + + 2341.816 / (2. + n) + - 2570.73 / (3. + n); + + #[rustfmt::skip] + let P3NMA12 = -17.0616 / n.powu(6) + - 19.53254 / n.powu(3) + + 337.931 * lm11(c) + - 1360.04 * lm11m1(c) + + 38.7337 * lm13m1(c) + - (367.64646999999997 / n) + + 335.995 / (1. + n) + - 1269.915 / (2. + n) + + 1605.91 / (3. + n); + + // nf^2 (parametrized) and nf^3 (exact) + #[rustfmt::skip] + let P3NSMA2 = -( + -193.84583328013258 + - 23.7037032 / n.powu(5) + + 117.5967 / n.powu(4) + - 256.5896 / n.powu(3) + + 437.881 / n.powu(2) + + 720.385709813466 / n + - 48.720000000000006 / (1. + n).powu(4) + + 189.51000000000002 / (1. + n).powu(3) + + 391.02500000000003 / (1. + n).powu(2) + + 367.4750000000001 / (1. + n) + + 404.47249999999997 / (2. + n) + - 2063.325 / ((1. + n).powu(2) * (2. + n)) + - (1375.55 * n) / ((1. + n).powu(2) * (2. + n)) + + 687.775 / ((1. + n) * (2. + n)) + - 81.71999999999998 / (3. + n) + + 114.9225 / (4. + n) + + 195.5772 * S1 + - (817.725 * S1) / n.powu(2) + + (714.46361 * S1) / n + - (687.775 * S1) / (1. + n) + - (817.725 * S2) / n + ); + let eta = 1. / n * 1. / (n + 1.); + #[rustfmt::skip] + let P3NSA3 = -CF * ( + -32. / 27. * ZETA3 * eta + - 16. / 9. * ZETA3 + - 16. / 27. * eta.powu(4) + - 16. / 81. * eta.powu(3) + + 80. / 27. * eta.powu(2) + - 320. / 81. * eta + + 32. / 27. * 1. / (n + 1.).powu(4) + + 128. / 27. * 1. / (n + 1.).powu(2) + + 64. / 27. * S1 * ZETA3 + - 32. / 81. * S1 + - 32. / 81. * S2 + - 160. / 81. * S3 + + 32. / 27. * S4 + + 131. / 81. + ); + + let mut result = cmplx!(0., 0.); + + // Assembly regular piece. + let nf = nf as f64; + let P3NSMAI = P3NSA0 + nf * P3NSA1 + nf.pow(3) * P3NSA3 + nf.pow(2) * P3NSMA2; + result += match variation { + 1 => P3NSMAI + P3NMA01 + nf * P3NMA11, + 2 => P3NSMAI + P3NMA02 + nf * P3NMA12, + _ => P3NSMAI + 0.5 * ((P3NMA01 + P3NMA02) + nf * (P3NMA11 + P3NMA12)), + }; + + // The singular piece. + let A4qI = 2.120902 * f64::pow(10., 4) - 5.179372 * f64::pow(10., 3) * nf + // + 1.955772 * f64::pow(10.,2) * nf.pow(2) + // + 3.272344 * nf.pow(3) + ; + let A4ap1 = -511.228 + 7.08645 * nf; + let A4ap2 = -502.481 + 7.82077 * nf; + let D1 = 1. / n - S1; + result += match variation { + 1 => (A4qI + A4ap1) * D1, + 2 => (A4qI + A4ap2) * D1, + _ => (A4qI + 0.5 * (A4ap1 + A4ap2)) * D1, + }; + + // ..The local piece. + let B4qI = 2.579609 * f64::pow(10., 4) + 0.08 - (5.818637 * f64::pow(10., 3) + 0.97) * nf + // + (1.938554 * f64::pow(10.,2) + 0.0037) * nf.pow(2) + // + 3.014982 * nf.pow(3) + ; + let B4ap1 = -2426.05 + 266.674 * nf - 0.05 * nf; + let B4ap2 = -2380.255 + 270.518 * nf - 0.05 * nf; + result += match variation { + 1 => B4qI + B4ap1, + 2 => B4qI + B4ap2, + _ => B4qI + 0.5 * (B4ap1 + B4ap2), + }; + + -result +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_quark_number_conservation() { + let NF = 5; + let mut c = Cache::new(cmplx!(1., 0.)); + let refs: [f64; 3] = [0.06776363, 0.064837, 0.07069]; + for var in [0, 1, 2] { + let test_value = gamma_nsm(&mut c, NF, var); + assert_approx_eq_cmplx!(f64, test_value, cmplx!(refs[var as usize], 0.), rel = 6e-5); + } + } + + #[test] + fn test_reference_moments() { + let NF = 4; + let nsm_nf4_refs: [f64; 7] = [ + 4322.890485339998, + 5491.581109692005, + 6221.256799360004, + 6774.606221595994, + 7229.056043916002, + 7618.358743427995, + 7960.658678124, + ]; + for N in [3.0, 5.0, 7.0, 9.0, 11.0, 13.0, 15.0] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_nsm(&mut c, NF, 0); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(nsm_nf4_refs[((N - 2.) / 2.) as usize], 0.), + rel = 8e-5 + ); + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsp.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsp.rs new file mode 100644 index 000000000..fafe956ef --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsp.rs @@ -0,0 +1,225 @@ +/// Compute the singlet-like non-singlet anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::cmplx; +use crate::constants::{CF, ZETA3}; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{lm11, lm11m1, lm12m1, lm13m1}; + +/// Compute the singlet-like non-singlet anomalous dimension. +/// +/// See [gamma_nsm][super::gamma_nsm] for implementation details. +pub fn gamma_nsp(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + + // Leading large-n_c, nf^0 and nf^1, parametrized + #[rustfmt::skip] + let P3NSA0 = 360.0 / n.powu(7) + - 1920.0 / n.powu(6) + + 7147.812 / n.powu(5) + - 17179.356 / n.powu(4) + + 34241.9 / n.powu(3) + - 51671.329999999994 / n.powu(2) + + 19069.8 * lm11(c) + - (491664.8019540468 / n) + - 4533.0 / (1. + n).powu(3) + - 11825.0 / (1. + n).powu(2) + + 129203.0 / (1. + n) + - 254965.0 / (2. + n) + + 83377.5 / (3. + n) + - 45750.0 / (4. + n) + + (49150.0 * (6.803662258392675 + n) * S1) / (n.powu(2) * (1.0 + n)) + + (334400.0 * S2) / n; + #[rustfmt::skip] + let P3NSA1 = 160.0 / n.powu(6) + - 864.0 / n.powu(5) + + 2583.1848 / n.powu(4) + - 5834.624 / n.powu(3) + + 9239.374 / n.powu(2) + - 3079.76 * lm11(c) + - (114047.0 / n) + - 465.0 / (1. + n).powu(4) + - 1230.0 / (1. + n).powu(3) + + 7522.5 / (1. + n).powu(2) + + 55669.3 / (1. + n) + - 43057.8 / (2. + n) + + 13803.8 / (3. + n) + - 7896.0 / (4. + n) + - (120.0 * (-525.063 + n) * S1) / (n.powu(2) * (1.0 + n)) + + (63007.5 * S2) / n; + + // Nonleading large-n_c, nf^0 and nf^1: two approximations + #[rustfmt::skip] + let P3NPA01 = - 107.16 / n.powu(7) + + 339.753 / n.powu(6) + - 1341.01 / n.powu(5) + + 2412.94 / n.powu(4) + - 3678.88 / n.powu(3) + - 2118.87 * lm11(c) + - 1777.27 * lm12m1(c) + - 204.183 * lm13m1(c) + + 1853.56 / n + - 8877.38 / (1. + n) + + 7393.83 / (2. + n) + - 2464.61 / (3. + n); + #[rustfmt::skip] + let P3NPA02 = - 107.16 / n.powu(7) + + 339.753 / n.powu(6) + - 1341.01 / n.powu(5) + + 379.152 / n.powu(3) + - 1389.73 / n.powu(2) + - 2118.87 * lm11(c) + - 173.936 * lm12m1(c) + + 223.078 * lm13m1(c) + - (2096.54 / n) + + 8698.39 / (1. + n) + - 19188.9 / (2. + n) + + 10490.5 / (3. + n); + #[rustfmt::skip] + let P3NPA11 = -33.5802 / n.powu(6) + + 111.802 / n.powu(5) + + 50.772 / n.powu(4) + - 118.608 / n.powu(3) + + 337.931 * lm11(c) + - 143.813 * lm11m1(c) + - 18.8803 * lm13m1(c) + + 304.82503 / n + - 1116.34 / (1. + n) + + 2187.58 / (2. + n) + - 1071.24 / (3. + n); + #[rustfmt::skip] + let P3NPA12 = - 33.5802 / n.powu(6) + + 111.802 / n.powu(5) + - 204.341 / n.powu(4) + + 267.404 / n.powu(3) + + 337.931 * lm11(c) + - 745.573 * lm11m1(c) + + 8.61438 * lm13m1(c) + - 385.52331999999996 / n + + 690.151 / (1. + n) + - 656.386 / (2. + n) + + 656.386 / (3. + n); + + // nf^2 (parametrized) and nf^3 (exact) + #[rustfmt::skip] + let P3NSPA2 = -( + -193.85906555742952 + - 18.962964 / n.powu(5) + + 99.1605 / n.powu(4) + - 225.141 / n.powu(3) + + 393.0056000000001 / n.powu(2) + - 403.50217685814835 / n + - 34.425000000000004 / (1. + n).powu(4) + + 108.42 / (1. + n).powu(3) + - 93.8225 / (1. + n).powu(2) + + 534.725 / (1. + n) + + 246.50250000000003 / (2. + n) + - 25.455 / ((1. + n).powu(2) * (2. + n)) + - (16.97 * n) / ((1. + n).powu(2) * (2. + n)) + + 8.485 / ((1. + n) * (2. + n)) + - 110.015 / (3. + n) + + 78.9875 / (4. + n) + + 195.5772 * S1 + - (101.0775 * S1) / n.powu(2) + + (35.17361 * S1) / n + - (8.485 * S1) / (1. + n) + - (101.0775 * S2) / n + ); + let eta = 1. / n * 1. / (n + 1.); + #[rustfmt::skip] + let P3NSA3 = -CF * ( + -32. / 27. * ZETA3 * eta + - 16. / 9. * ZETA3 + - 16. / 27. * eta.powu(4) + - 16. / 81. * eta.powu(3) + + 80. / 27. * eta.powu(2) + - 320. / 81. * eta + + 32. / 27. * 1. / (n + 1.).powu(4) + + 128. / 27. * 1. / (n + 1.).powu(2) + + 64. / 27. * S1 * ZETA3 + - 32. / 81. * S1 + - 32. / 81. * S2 + - 160. / 81. * S3 + + 32. / 27. * S4 + + 131. / 81. + ); + + let mut result = cmplx!(0., 0.); + + // Assembly regular piece. + let nf = nf as f64; + let P3NSPAI = P3NSA0 + nf * P3NSA1 + nf.pow(2) * P3NSPA2 + nf.pow(3) * P3NSA3; + result += match variation { + 1 => P3NSPAI + P3NPA01 + nf * P3NPA11, + 2 => P3NSPAI + P3NPA02 + nf * P3NPA12, + _ => P3NSPAI + 0.5 * ((P3NPA01 + P3NPA02) + nf * (P3NPA11 + P3NPA12)), + }; + + // The singular piece. + let A4qI = 2.120902 * f64::pow(10., 4) - 5.179372 * f64::pow(10., 3) * nf + // + 1.955772 * f64::pow(10.,2) * nf.pow(2) + // + 3.272344 * nf.pow(3) + ; + let A4ap1 = -507.152 + 7.33927 * nf; + let A4ap2 = -505.209 + 7.53662 * nf; + let D1 = 1. / n - S1; + result += match variation { + 1 => (A4qI + A4ap1) * D1, + 2 => (A4qI + A4ap2) * D1, + _ => (A4qI + 0.5 * (A4ap1 + A4ap2)) * D1, + }; + + // ..The local piece. + let B4qI = 2.579609 * f64::pow(10., 4) + 0.08 - (5.818637 * f64::pow(10., 3) + 0.97) * nf + // + (1.938554 * f64::pow(10.,2) + 0.0037) * nf.pow(2) + // + 3.014982 * nf.pow(3) + ; + let B4ap1 = -2405.03 + 267.965 * nf; + let B4ap2 = -2394.47 + 269.028 * nf; + result += match variation { + 1 => B4qI + B4ap1, + 2 => B4qI + B4ap2, + _ => B4qI + 0.5 * (B4ap1 + B4ap2), + }; + + -result +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_reference_moments() { + let NF = 4; + let nsp_nf4_refs: [f64; 8] = [ + 3679.6690577439995, + 5066.339235808004, + 5908.005605364002, + 6522.700744595994, + 7016.383458928004, + 7433.340927783997, + 7796.397038483998, + 8119.044600816003, + ]; + for N in [2.0, 4.0, 6.0, 8.0, 10.0, 12.0, 14.0, 16.0] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_nsp(&mut c, NF, 0); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(nsp_nf4_refs[((N - 2.) / 2.) as usize], 0.), + rel = 4e-5 + ); + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsv.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsv.rs new file mode 100644 index 000000000..07bb2bce5 --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gnsv.rs @@ -0,0 +1,120 @@ +/// Compute the valence-like non-singlet anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use super::gnsm; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{lm11m1, lm12m1, lm13m1}; + +/// Compute the sea-like non-singlet anomalous dimension. +fn gamma_nss(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + + // nf^1: two approximations + #[rustfmt::skip] + let P3NSA11 = 2880. / n.powu(7) + - 11672.4 / n.powu(6) + + 12802.560000000001 / n.powu(5) + - 7626.66 / n.powu(4) + + 6593.2 / n.powu(3) + - 3687.6 / n.powu(2) + + 4989.2 / (1. + n) + - 6596.93 / (2. + n) + + 1607.73 / (3. + n) + + 60.4 * lm12m1(c) + + 4.685 * lm13m1(c); + #[rustfmt::skip] + let P3NSA12 = -2880. / n.powu(7) + + 4066.32 / n.powu(6) + - 5682.24 / n.powu(5) + + 5540.88 / n.powu(4) + + 546.1 / n.powu(3) + - 2987.83 / n.powu(2) + + 2533.54 / n + - 1502.75 / (1. + n) + - 2297.56 / (2. + n) + + 1266.77 / (3. + n) + - 254.63 * lm11m1(c) + - 0.28953 * lm13m1(c); + + // nf^2 (parametrized) + #[rustfmt::skip] + let P3NSSA2 = 47.4074 / n.powu(6) + - 142.222 / n.powu(5) + + 32.1201 / n.powu(4) + - 132.824 / n.powu(3) + + 647.397 / n.powu(2) + + 19.7 * lm11m1(c) + - 3.43547 * lm12m1(c) + - 1262.0951538579698 / n + - 187.17000000000002 / (1. + n).powu(4) + + 453.885 / (1. + n).powu(3) + + 147.01749999999998 / (1. + n).powu(2) + + 1614.1000000000001 / (1. + n) + - 380.12500000000006 / (2. + n) + - 42.575 / (3. + n) + + (42.977500000000006 * S2) / n + + (0.0900000000000047 * (477.52777777775293 + n) * S1) / (n.powu(2) * (1. + n)); + + let nf = nf as f64; + let P3NSSA: Complex = match variation { + 1 => nf * P3NSA11 + nf.pow(2) * P3NSSA2, + 2 => nf * P3NSA12 + nf.pow(2) * P3NSSA2, + _ => 0.5 * nf * (P3NSA11 + P3NSA12) + nf.pow(2) * P3NSSA2, + }; + -P3NSSA +} + +/// Compute the valence non-singlet anomalous dimension. +/// +/// See [gamma_nsm][super::gamma_nsm] for implementation details. +pub fn gamma_nsv(c: &mut Cache, nf: u8, variation: u8) -> Complex { + gnsm::gamma_nsm(c, nf, variation) + gamma_nss(c, nf, variation) +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_quark_number_conservation() { + let NF = 5; + let mut c = Cache::new(cmplx!(1., 0.)); + let refs: [f64; 3] = [-0.01100459, -0. - 0.00779938, -0.0142098]; + for var in [0, 1, 2] { + let test_value = gamma_nss(&mut c, NF, var); + assert_approx_eq_cmplx!(f64, test_value, cmplx!(refs[var as usize], 0.), rel = 5e-6); + } + } + + #[test] + fn test_reference_moments() { + let NF = 4; + let nss_nf4_refs: [f64; 8] = [ + 50.10532524, + 39.001939964, + 21.141505811200002, + 12.4834195012, + 8.0006134908, + 5.4610639744, + 3.9114290952, + 2.90857799, + ]; + for N in [3.0, 5.0, 7.0, 9.0, 11.0, 13.0, 15.0, 17.0] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_nss(&mut c, NF, 0); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(nss_nf4_refs[((N - 2.) / 2.) as usize], 0.), + rel = 5e-5 + ); + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gps.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gps.rs new file mode 100644 index 000000000..70ae4d42b --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gps.rs @@ -0,0 +1,179 @@ +/// Compute the pure-singlet quark-to-quark anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::harmonics::cache::Cache; +use crate::harmonics::log_functions::{lm11m1, lm12m1, lm12m2, lm13m1, lm13m2, lm14m1, lm14m2}; + +/// Compute the pure-singlet quark-to-quark anomalous dimension. +/// +/// The routine is taken from [\[Falcioni:2023luc\]][crate::bib::Falcioni2023luc]. +/// +/// These are approximations for fixed `nf` = 3, 4 and 5 based on the +/// first ten even moments together with small-x/large-x constraints. +/// The two sets spanning the error estimate are called via `variation = 1` +/// and `variation = 2`. Any other value of `variation` invokes their average. +pub fn gamma_ps(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let nf_ = nf as f64; + let nf2 = nf_.pow(2); + let nf3 = nf_.pow(3); + let xm1lm1 = -(1. / (-1. + n).powu(2)) + 1. / n.powu(2); + + // Known large-x coefficients + let x1L4cff = -5.6460905 * 10. * nf_ + 3.6213992 * nf2; + let x1L3cff = -2.4755054 * f64::pow(10., 2) * nf_ + 4.0559671 * 10. * nf2 - 1.5802469 * nf3; + let y1L4cff = -1.3168724 * 10. * nf_; + let y1L3cff = -1.9911111 * f64::pow(10., 2) * nf_ + 1.3695473 * 10. * nf2; + + // Known small-x coefficients + let bfkl1 = 1.7492273 * f64::pow(10., 3) * nf_; + let x0L6cff = -7.5061728 * nf_ + 7.9012346 * 0.1 * nf2; + let x0L5cff = 2.8549794 * 10. * nf_ + 3.7925926 * nf2; + let x0L4cff = + -8.5480010 * f64::pow(10., 2) * nf_ + 7.7366255 * 10. * nf2 - 1.9753086 * 0.1 * nf3; + + // The resulting part of the function + let P3ps01 = bfkl1 * 2. / (-1. + n).powu(3) + + x0L6cff * 720. / n.powu(7) + + x0L5cff * -120. / n.powu(6) + + x0L4cff * 24. / n.powu(5) + + x1L3cff * lm13m1(c) + + x1L4cff * lm14m1(c) + + y1L3cff * lm13m2(c) + + y1L4cff * lm14m2(c); + + // The selected approximations for nf = 3, 4, 5 + let P3psApp1: Complex; + let P3psApp2: Complex; + if nf == 3 { + P3psApp1 = P3ps01 + 67731.0 * xm1lm1 + 274100.0 * 1. / ((-1. + n) * n) + - 104493.0 * (1. / n - n / (2. + 3. * n + n.powu(2))) + + 34403.0 * 1. / (6. + 5. * n + n.powu(2)) + + 353656.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 10620.0 * 2. / n.powu(3) + + 40006.0 * -6. / n.powu(4) + - 7412.1 * lm11m1(c) + - 2365.1 * lm12m1(c) + + 1533.0 * lm12m2(c); + P3psApp2 = P3ps01 + 54593.0 * xm1lm1 + 179748.0 * 1. / ((-1. + n) * n) + - 195263.0 * 1. / (n + n.powu(2)) + + 12789.0 * 2. / (3. + 4. * n + n.powu(2)) + + 4700.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + - 103604.0 * 2. / n.powu(3) + - 2758.3 * -6. / n.powu(4) + - 2801.2 * lm11m1(c) + - 1986.9 * lm12m1(c) + - 6005.9 * lm12m2(c); + } else if nf == 4 { + P3psApp1 = P3ps01 + 90154.0 * xm1lm1 + 359084.0 * 1. / ((-1. + n) * n) + - 136319.0 * (1. / n - n / (2. + 3. * n + n.powu(2))) + + 45379.0 * 1. / (6. + 5. * n + n.powu(2)) + + 461167.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 13869.0 * 2. / n.powu(3) + + 52525.0 * -6. / n.powu(4) + - 7498.2 * lm11m1(c) + - 2491.5 * lm12m1(c) + + 1727.2 * lm12m2(c); + P3psApp2 = P3ps01 + 72987.0 * xm1lm1 + 235802.0 * 1. / ((-1. + n) * n) + - 254921.0 * 1. / (n + n.powu(2)) + + 17138.0 * 2. / (3. + 4. * n + n.powu(2)) + + 5212.9 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + - 135378.0 * 2. / n.powu(3) + - 3350.9 * -6. / n.powu(4) + - 1472.7 * lm11m1(c) + - 1997.2 * lm12m1(c) + - 8123.3 * lm12m2(c); + } else if nf == 5 { + P3psApp1 = P3ps01 + 112481.0 * xm1lm1 + 440555.0 * 1. / ((-1. + n) * n) + - 166581.0 * (1. / n - n / (2. + 3. * n + n.powu(2))) + + 56087.0 * 1. / (6. + 5. * n + n.powu(2)) + + 562992.0 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + + 16882.0 * 2. / n.powu(3) + + 64577.0 * -6. / n.powu(4) + - 6570.1 * lm11m1(c) + - 2365.7 * lm12m1(c) + + 1761.7 * lm12m2(c); + P3psApp2 = P3ps01 + 91468.0 * xm1lm1 + 289658.0 * 1. / ((-1. + n) * n) + - 311749.0 * 1. / (n + n.powu(2)) + + 21521.0 * 2. / (3. + 4. * n + n.powu(2)) + + 4908.9 * (-(1. / n.powu(2)) + 1. / (1. + n).powu(2)) + - 165795.0 * 2. / n.powu(3) + - 3814.9 * -6. / n.powu(4) + + 804.5 * lm11m1(c) + - 1760.8 * lm12m1(c) + - 10295.0 * lm12m2(c); + } else { + panic!("nf=6 is not available at N3LO"); + } + // We return (for now) one of the two error-band boundaries + // or the present best estimate, their average + + let P3psA = match variation { + 1 => P3psApp1, + 2 => P3psApp2, + _ => 0.5 * (P3psApp1 + P3psApp2), + }; + -P3psA +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + // Test the prediction of N=22 wrt to [Falcioni:2023luc]. + fn test_gamma_ps_extrapolation() { + let n22_ref: [f64; 3] = [6.2478570, 10.5202730, 15.6913948]; + let mut c = Cache::new(cmplx!(22., 0.)); + for NF in [3, 4, 5] { + let test_value = gamma_ps(&mut c, NF, 0); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(n22_ref[(NF - 3) as usize], 0.), + rel = 3e-5 + ); + } + } + + #[test] + fn test_reference_moments() { + fn qq3ps_moment(N: usize, nf: f64) -> f64 { + let nf2 = nf * nf; + let nf3 = nf2 * nf; + // From Eq. 14 of [Falcioni:2023luc]. + let mom_list = [ + -691.5937093082381 * nf + 84.77398149891167 * nf2 + 4.4669568492355864 * nf3, + -109.33023358432462 * nf + 8.77688525974872 * nf2 + 0.3060771365698822 * nf3, + -46.030613749542226 * nf + 4.744075766957513 * nf2 + 0.042548957282380874 * nf3, + -24.01455020567638 * nf + 3.235193483272451 * nf2 - 0.007889256298951614 * nf3, + -13.730393879922417 * nf + 2.3750187592472374 * nf2 - 0.02102924056123573 * nf3, + -8.152592251923657 * nf + 1.8199581788320662 * nf2 - 0.024330231290833188 * nf3, + -4.8404471801109565 * nf + 1.4383273806219803 * nf2 - 0.024479943136069916 * nf3, + -2.7511363301137024 * nf + 1.164299642517469 * nf2 - 0.023546009234463816 * nf3, + -1.375969240387974 * nf + 0.9608733183576097 * nf2 - 0.022264393374041958 * nf3, + -0.4426815682220422 * nf + 0.8057453328332964 * nf2 - 0.02091826436475512 * nf3, + ]; + mom_list[(N - 2) / 2] + } + for variation in [0, 1, 2] { + for NF in [3, 4, 5] { + for N in [2.0, 4.0, 6.0, 8.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_ps(&mut c, NF, variation); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(qq3ps_moment(N as usize, NF as f64), 0.), + rel = 4e-4 + ); + } + } + } + } +} diff --git a/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gqg.rs b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gqg.rs new file mode 100644 index 000000000..19cb4dda9 --- /dev/null +++ b/crates/ekore/src/anomalous_dimensions/unpolarized/spacelike/as4/gqg.rs @@ -0,0 +1,171 @@ +/// Compute the singlet gluon-to-quark anomalous dimension. +use num::complex::Complex; +use num::traits::Pow; + +use crate::constants::ZETA2; +use crate::harmonics::cache::{Cache, K}; +use crate::harmonics::log_functions::{lm11, lm12, lm13, lm14, lm14m1, lm15, lm15m1}; + +/// Compute the singlet gluon-to-quark anomalous dimension. +/// +/// The routine is taken from [\[Falcioni:2023vqq\]][crate::bib::Falcioni2023vqq]. +/// +/// These are approximations for fixed `nf` = 3, 4 and 5 based on the +/// first 10 even moments together with small-x/large-x constraints. +/// The two sets indicating the error estimate are called via `variation = 1` +/// and `variation = 2`. Any other value of `variation` invokes their average. +pub fn gamma_qg(c: &mut Cache, nf: u8, variation: u8) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let nf_ = nf as f64; + let nf2 = nf_.pow(2); + let nf3 = nf_.pow(3); + + // Known large-x coefficients + let x1L5cff = 1.8518519 * nf_ - 4.1152263 * 0.1 * nf2; + let x1L4cff = 3.5687794 * 10. * nf_ - 3.5116598 * nf2 - 8.2304527 * 0.01 * nf3; + let y1L5cff = 2.8806584 * nf_ + 8.2304527 * 0.1 * nf2; + let y1L4cff = -4.0511391 * 10. * nf_ + 5.5418381 * nf2 + 1.6460905 * 0.1 * nf3; + + // Known small-x coefficients + let bfkl1 = 3.9357613 * f64::pow(10., 3) * nf_; + let x0L6cff = -1.9588477 * 10. * nf_ + 2.7654321 * nf2; + let x0L5cff = 2.1573663 * 10. * nf_ + 1.7244444 * 10. * nf2; + let x0L4cff = + -2.8667643 * f64::pow(10., 3) * nf_ + 3.0122403 * f64::pow(10., 2) * nf2 + 4.1316872 * nf3; + + // The resulting part of the function + let P3QG01 = bfkl1 * 2. / (-1. + n).powu(3) + + x0L6cff * 720. / n.powu(7) + + x0L5cff * -120. / n.powu(6) + + x0L4cff * 24. / n.powu(5) + + x1L4cff * lm14(c) + + x1L5cff * lm15(c) + + y1L4cff * lm14m1(c) + + y1L5cff * lm15m1(c); + + // The selected approximations for nf = 3, 4, 5 + let P3qgApp1: Complex; + let P3qgApp2: Complex; + if nf == 3 { + P3qgApp1 = P3QG01 + 187500.0 * -(1. / (-1. + n).powu(2)) + 826060.0 * 1. / ((-1. + n) * n) + - 150474.0 * 1. / n + + 226254.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 577733.0 * -1. / n.powu(2) + - 180747.0 * 2. / n.powu(3) + + 95411.0 * -6. / n.powu(4) + + 119.8 * lm13(c) + + 7156.3 * lm12(c) + + 45790.0 * lm11(c) + - 95682.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + P3qgApp2 = P3QG01 + 135000.0 * -(1. / (-1. + n).powu(2)) + 484742.0 * 1. / ((-1. + n) * n) + - 11627.0 * 1. / n + - 187478.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 413512.0 * -1. / n.powu(2) + - 82500.0 * 2. / n.powu(3) + + 29987.0 * -6. / n.powu(4) + - 850.1 * lm13(c) + - 11425.0 * lm12(c) + - 75323.0 * lm11(c) + + 282836.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else if nf == 4 { + P3qgApp1 = P3QG01 + 250000.0 * -(1. / (-1. + n).powu(2)) + 1089180.0 * 1. / ((-1. + n) * n) + - 241088.0 * 1. / n + + 342902.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 720081.0 * -1. / n.powu(2) + - 247071.0 * 2. / n.powu(3) + + 126405.0 * -6. / n.powu(4) + + 272.4 * lm13(c) + + 10911.0 * lm12(c) + + 60563.0 * lm11(c) + - 161448.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + P3qgApp2 = P3QG01 + 180000.0 * -(1. / (-1. + n).powu(2)) + 634090.0 * 1. / ((-1. + n) * n) + - 55958.0 * 1. / n + - 208744.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 501120.0 * -1. / n.powu(2) + - 116073.0 * 2. / n.powu(3) + + 39173.0 * -6. / n.powu(4) + - 1020.8 * lm13(c) + - 13864.0 * lm12(c) + - 100922.0 * lm11(c) + + 343243.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else if nf == 5 { + P3qgApp1 = P3QG01 + 312500.0 * -(1. / (-1. + n).powu(2)) + 1345700.0 * 1. / ((-1. + n) * n) + - 350466.0 * 1. / n + + 480028.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 837903.0 * -1. / n.powu(2) + - 315928.0 * 2. / n.powu(3) + + 157086.0 * -6. / n.powu(4) + + 472.7 * lm13(c) + + 15415.0 * lm12(c) + + 75644.0 * lm11(c) + - 244869.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + P3qgApp2 = P3QG01 + 225000.0 * -(1. / (-1. + n).powu(2)) + 776837.0 * 1. / ((-1. + n) * n) + - 119054.0 * 1. / n + - 209530.0 * (3. + n) / (2. + 3. * n + n.powu(2)) + + 564202.0 * -1. / n.powu(2) + - 152181.0 * 2. / n.powu(3) + + 48046.0 * -6. / n.powu(4) + - 1143.8 * lm13(c) + - 15553.0 * lm12(c) + - 126212.0 * lm11(c) + + 385995.0 * (S1 - n * (ZETA2 - S2)) / n.powu(2); + } else { + panic!("nf=6 is not available at N3LO"); + } + + // We return (for now) one of the two error-band boundaries + // or the present best estimate, their average + let P3QGA = match variation { + 1 => P3qgApp1, + 2 => P3qgApp2, + _ => 0.5 * (P3qgApp1 + P3qgApp2), + }; + -P3QGA +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + + #[test] + fn test_reference_moments() { + fn qg3_moment(N: usize, nf: f64) -> f64 { + let nf2 = nf * nf; + let nf3 = nf2 * nf; + // From Eq. 4 of [Falcioni:2023vqq]. + let mom_list = [ + -654.4627782205557 * nf + 245.61061978871788 * nf2 - 0.9249909688301847 * nf3, + 290.31106867034487 * nf - 76.51672403736478 * nf2 - 4.911625629947491 * nf3, + 335.80080466045274 * nf - 124.57102255718002 * nf2 - 4.193871425027802 * nf3, + 294.58768309440677 * nf - 135.3767647714609 * nf2 - 3.609775642729055 * nf3, + 241.6153399044715 * nf - 135.18742470907011 * nf2 - 3.189394834180898 * nf3, + 191.97124640777176 * nf - 131.16316638326697 * nf2 - 2.8771044305171913 * nf3, + 148.5682948286098 * nf - 125.82310814280595 * nf2 - 2.635918561148907 * nf3, + 111.34042526856348 * nf - 120.16819876888667 * nf2 - 2.4433790398202664 * nf3, + 79.51561588665083 * nf - 114.61713540075442 * nf2 - 2.28548686108789 * nf3, + 52.24329555231736 * nf - 109.34248910828198 * nf2 - 2.1531537251387527 * nf3, + ]; + mom_list[(N - 2) / 2] + } + for variation in [0, 1, 2] { + for NF in [3, 4, 5] { + for N in [2.0, 4.0, 6.0, 8.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0] { + let mut c = Cache::new(cmplx!(N, 0.)); + let test_value = gamma_qg(&mut c, NF, variation); + assert_approx_eq_cmplx!( + f64, + test_value, + cmplx!(qg3_moment(N as usize, NF as f64), 0.), + rel = 8e-4 + ); + } + } + } + } +} diff --git a/crates/ekore/src/bib.rs b/crates/ekore/src/bib.rs index 89ca13ccf..8d901981f 100644 --- a/crates/ekore/src/bib.rs +++ b/crates/ekore/src/bib.rs @@ -1,4 +1,4 @@ -//! List of References (autogenerated on 2025-01-13T15:14:10.528903). +//! List of References (autogenerated on 2025-02-24T12:35:24.669515). #[allow(non_snake_case)] /// The Three loop splitting functions in QCD: The Nonsinglet case @@ -155,3 +155,63 @@ pub fn deFlorian2015ujt() {} /// /// DOI: [10.1007/JHEP10(2016)056](https:dx.doi.org/10.1007/JHEP10(2016)056) pub fn deFlorian2016gvk() {} + +#[allow(non_snake_case)] +/// Four-Loop Non-Singlet Splitting Functions in the Planar Limit and Beyond +/// +/// Moch, S. and Ruijl, B. and Ueda, T. and Vermaseren, J. A. M. and Vogt, A. +/// +/// Published in: JHEP 10 (2017), 041 +/// +/// e-Print: [1707.08315](https://arxiv.org/abs/1707.08315) +/// +/// DOI: [10.1007/JHEP10(2017)041](https:dx.doi.org/10.1007/JHEP10(2017)041) +pub fn Moch2017uml() {} + +#[allow(non_snake_case)] +/// Four-loop splitting functions in QCD - the gluon-gluon case - +/// +/// Falcioni, G. and Herzog, F. and Moch, S. and Pelloni, A. and Vogt, A. +/// +/// Published in: Phys. Lett. B 860 (2025), 139194 +/// +/// e-Print: [2410.08089](https://arxiv.org/abs/2410.08089) +/// +/// DOI: [10.1016/j.physletb.2024.139194](https:dx.doi.org/10.1016/j.physletb.2024.139194) +pub fn Falcioni2024qpd() {} + +#[allow(non_snake_case)] +/// Four-loop splitting functions in QCD - The quark-to-gluon case +/// +/// Falcioni, G. and Herzog, F. and Moch, S. and Pelloni, A. and Vogt, A. +/// +/// Published in: Phys. Lett. B 856 (2024), 138906 +/// +/// e-Print: [2404.09701](https://arxiv.org/abs/2404.09701) +/// +/// DOI: [10.1016/j.physletb.2024.138906](https:dx.doi.org/10.1016/j.physletb.2024.138906) +pub fn Falcioni2024xyt() {} + +#[allow(non_snake_case)] +/// Four-loop splitting functions in QCD - The quark-quark case +/// +/// Falcioni, G. and Herzog, F. and Moch, S. and Vogt, A. +/// +/// Published in: Phys. Lett. B 842 (2023), 137944 +/// +/// e-Print: [2302.07593](https://arxiv.org/abs/2302.07593) +/// +/// DOI: [10.1016/j.physletb.2023.137944](https:dx.doi.org/10.1016/j.physletb.2023.137944) +pub fn Falcioni2023luc() {} + +#[allow(non_snake_case)] +/// Four-loop splitting functions in QCD - The gluon-to-quark case +/// +/// Falcioni, G. and Herzog, F. and Moch, S. and Vogt, A. +/// +/// Published in: Phys. Lett. B 846 (2023), 138215 +/// +/// e-Print: [2307.04158](https://arxiv.org/abs/2307.04158) +/// +/// DOI: [10.1016/j.physletb.2023.138215](https:dx.doi.org/10.1016/j.physletb.2023.138215) +pub fn Falcioni2023vqq() {} diff --git a/crates/ekore/src/constants.rs b/crates/ekore/src/constants.rs index 98023e016..a66cf71be 100644 --- a/crates/ekore/src/constants.rs +++ b/crates/ekore/src/constants.rs @@ -43,6 +43,9 @@ pub const ZETA3: f64 = 1.2020569031595942; /// $\zeta(4) = \pi^4 / 90$. pub const ZETA4: f64 = 1.082323233711138; +/// Riemann zeta function at z = 5. +pub const ZETA5: f64 = 1.03692775514337; + /// singlet-like non-singlet |PID|. pub const PID_NSP: u16 = 10101; diff --git a/crates/ekore/src/harmonics.rs b/crates/ekore/src/harmonics.rs index a8a3b3241..75dcfd41b 100644 --- a/crates/ekore/src/harmonics.rs +++ b/crates/ekore/src/harmonics.rs @@ -2,8 +2,10 @@ pub mod cache; pub mod g_functions; +pub mod log_functions; pub mod polygamma; pub mod w1; pub mod w2; pub mod w3; pub mod w4; +pub mod w5; diff --git a/crates/ekore/src/harmonics/cache.rs b/crates/ekore/src/harmonics/cache.rs index da0f4800b..595834017 100644 --- a/crates/ekore/src/harmonics/cache.rs +++ b/crates/ekore/src/harmonics/cache.rs @@ -3,7 +3,7 @@ use num::{complex::Complex, Zero}; use std::collections::HashMap; -use crate::harmonics::{g_functions, w1, w2, w3, w4}; +use crate::harmonics::{g_functions, w1, w2, w3, w4, w5}; /// List of available elements. #[derive(Debug, PartialEq, Eq, Hash)] @@ -16,6 +16,8 @@ pub enum K { S3, /// $S_4(N)$ S4, + /// $S_5(N)$ + S5, /// $S_1(N/2)$ S1h, /// $S_2(N/2)$ @@ -83,21 +85,22 @@ impl Cache { K::S2 => w2::S2(self.n), K::S3 => w3::S3(self.n), K::S4 => w4::S4(self.n), + K::S5 => w5::S5(self.n), K::S1h => w1::S1(self.n / 2.), K::S2h => w2::S2(self.n / 2.), K::S3h => w3::S3(self.n / 2.), K::S1mh => w1::S1((self.n - 1.) / 2.), K::S2mh => w2::S2((self.n - 1.) / 2.), K::S3mh => w3::S3((self.n - 1.) / 2.), - K::G3 => g_functions::g3(self.n, self.get(K::S1)), - K::Sm1e => w1::Sm1e(self.get(K::S1), self.get(K::S1h)), - K::Sm1o => w1::Sm1o(self.get(K::S1), self.get(K::S1mh)), - K::Sm2e => w2::Sm2e(self.get(K::S2), self.get(K::S2h)), - K::Sm2o => w2::Sm2o(self.get(K::S2), self.get(K::S2mh)), - K::Sm3e => w3::Sm3e(self.get(K::S3), self.get(K::S3h)), - K::Sm3o => w3::Sm3o(self.get(K::S3), self.get(K::S3mh)), - K::Sm21e => w3::Sm21e(self.n, self.get(K::S1), self.get(K::Sm1e)), - K::Sm21o => w3::Sm21o(self.n, self.get(K::S1), self.get(K::Sm1o)), + K::G3 => g_functions::g3(self), + K::Sm1e => w1::Sm1e(self), + K::Sm1o => w1::Sm1o(self), + K::Sm2e => w2::Sm2e(self), + K::Sm2o => w2::Sm2o(self), + K::Sm3e => w3::Sm3e(self), + K::Sm3o => w3::Sm3o(self), + K::Sm21e => w3::Sm21e(self), + K::Sm21o => w3::Sm21o(self), }; // insert self.m.insert(k, val); @@ -140,7 +143,7 @@ mod tests { #[test] fn test_recursive_harmonic_sum() { - const SX: [fn(Complex) -> Complex; 4] = [w1::S1, w2::S2, w3::S3, w4::S4]; + const SX: [fn(Complex) -> Complex; 5] = [w1::S1, w2::S2, w3::S3, w4::S4, w5::S5]; const NS: [Complex; 2] = [cmplx!(1.0, 0.0), cmplx!(2.34, 3.45)]; const ITERS: [usize; 2] = [1, 2]; for sit in SX.iter().enumerate() { diff --git a/crates/ekore/src/harmonics/g_functions.rs b/crates/ekore/src/harmonics/g_functions.rs index 9b5348476..71264a977 100644 --- a/crates/ekore/src/harmonics/g_functions.rs +++ b/crates/ekore/src/harmonics/g_functions.rs @@ -2,6 +2,7 @@ use crate::constants::ZETA2; use crate::harmonics::cache::recursive_harmonic_sum as s; +use crate::harmonics::cache::{Cache, K}; use num::{complex::Complex, Zero}; /// Compute the Mellin transform of $\text{Li}_2(x)/(1+x)$. @@ -11,7 +12,9 @@ use num::{complex::Complex, Zero}; /// /// We use the name from [\[MuselliPhD\]](crate::bib::MuselliPhD), but not his implementation - rather we use the /// Pegasus [\[Vogt:2004ns\]](crate::bib::Vogt2004ns) implementation. -pub fn g3(N: Complex, S1: Complex) -> Complex { +pub fn g3(c: &mut Cache) -> Complex { + let N = c.n(); + let S1 = c.get(K::S1); const CS: [f64; 7] = [ 1.0000e0, -0.9992e0, 0.9851e0, -0.9005e0, 0.6621e0, -0.3174e0, 0.0699e0, ]; @@ -25,8 +28,8 @@ pub fn g3(N: Complex, S1: Complex) -> Complex { #[cfg(test)] mod tests { + use crate::harmonics::cache::Cache; use crate::harmonics::g_functions::g3; - use crate::harmonics::w1; use crate::{assert_approx_eq_cmplx, cmplx}; use num::complex::Complex; @@ -41,9 +44,9 @@ mod tests { ]; for it in NS.iter().enumerate() { let n = *it.1; - let s1 = w1::S1(n); + let mut c = Cache::new(n); let refval = REFVALS[it.0]; - let g3 = g3(n, s1); + let g3 = g3(&mut c); assert_approx_eq_cmplx!(f64, g3, refval, epsilon = 1e-6); } } diff --git a/crates/ekore/src/harmonics/log_functions.rs b/crates/ekore/src/harmonics/log_functions.rs new file mode 100644 index 000000000..1db2bcf9f --- /dev/null +++ b/crates/ekore/src/harmonics/log_functions.rs @@ -0,0 +1,285 @@ +//! Mellin transformation of logarithms. +//! +//! We provide transforms of: +//! +//! * $\ln^k(1-x), \quad k = 1,2,3,4$ +//! * $(1-x)\ln^k(1-x), \quad k = 1,2,3$ +//! * $(1-x)^2\ln^k(1-x), \quad k = 1,2,3$ +use super::cache::{Cache, K}; +use num::complex::Complex; +use num::pow; + +/// Mellin transform of $(1-x)\ln(1-x)$. +pub fn lm11m1(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + 1. / (1. + n).powu(2) - S1 / (1. + n).powu(2) - S1 / (n * (1. + n).powu(2)) +} + +/// Mellin transform of $(1-x)\ln^2(1-x)$. +pub fn lm12m1(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + -2. / (1. + n).powu(3) - (2. * S1) / (1. + n).powu(2) + S1.powu(2) / n - S1.powu(2) / (1. + n) + + S2 / n + - S2 / (1. + n) +} + +/// Mellin transform of $(1-x)\ln^3(1-x)$. +pub fn lm13m1(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + (3. * n * (1. + n).powu(2) * S1.powu(2) - (1. + n).powu(3) * S1.powu(3) + + 3. * n * (1. + n).powu(2) * S2 + - 3. * (1. + n) * S1 * (-2. * n + (1. + n).powu(2) * S2) + - 2. * (-3. * n + (1. + n).powu(3) * S3)) + / (n * (1. + n).powu(4)) +} + +/// Mellin transform of $(1-x)\ln^4(1-x)$. +pub fn lm14m1(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + (-24. * n - 4. * n * (1. + n).powu(3) * S1.powu(3) + (1. + n).powu(4) * S1.powu(4) + - 12. * n * (1. + n).powu(2) * S2 + + 3. * (1. + n).powu(4) * S2.powu(2) + + 6. * (1. + n).powu(2) * S1.powu(2) * (-2. * n + (1. + n).powu(2) * S2) + - 8. * n * S3 + - 24. * n.powu(2) * S3 + - 24. * n.powu(3) * S3 + - 8. * n.powu(4) * S3 + - 4. * (1. + n) + * S1 + * (3. * n * (1. + n).powu(2) * S2 - 2. * (-3. * n + (1. + n).powu(3) * S3)) + + 6. * S4 + + 24. * n * S4 + + 36. * n.powu(2) * S4 + + 24. * n.powu(3) * S4 + + 6. * n.powu(4) * S4) + / (n * (1. + n).powu(5)) +} + +/// Mellin transform of $(1-x)\ln^5(1-x)$. +pub fn lm15m1(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + let S5 = c.get(K::S5); + (1. / (n * (1. + n).powu(6))) + * (5. * n * (1. + n).powu(4) * S1.powu(4) - (1. + n).powu(5) * S1.powu(5) + + 15. * n * (1. + n).powu(4) * S2.powu(2) + - 10. * (1. + n).powu(3) * S1.powu(3) * (-2. * n + (1. + n).powu(2) * S2) + + 40. * n * (1. + n).powu(3) * S3 + - 20. * (1. + n).powu(2) * S2 * (-3. * n + (1. + n).powu(3) * S3) + + 10. + * S1.powu(2) + * (3. * n * (1. + n).powu(4) * S2 + - 2. * (1. + n).powu(2) * (-3. * n + (1. + n).powu(3) * S3)) + + 30. * n * (1. + n).powu(4) * S4 + - 5. * (1. + n) + * S1 + * (-12. * n * (1. + n).powu(2) * S2 + 3. * (1. + n).powu(4) * S2.powu(2) + - 8. * n * (1. + n).powu(3) * S3 + + 6. * (-4. * n + (1. + n).powu(4) * S4)) + - 24. * (-5. * n + (1. + n).powu(5) * S5)) +} + +/// Mellin transform of $\ln(1-x)$. +pub fn lm11(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + -S1 / n +} + +/// Mellin transform of $\ln^2(1-x)$. +pub fn lm12(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + (S1.powu(2) + S2) / n +} + +/// Mellin transform of $\ln^3(1-x)$. +pub fn lm13(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + -((S1.powu(3) + (3. * S1) * S2 + 2. * S3) / n) +} + +/// Mellin transform of $\ln^4(1-x)$. +pub fn lm14(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + (S1.powu(4) + 6. * S1.powu(2) * S2 + 3. * S2.powu(2) + 8. * S1 * S3 + 6. * S4) / n +} + +/// Mellin transform of $\ln^5(1-x)$. +pub fn lm15(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + let S5 = c.get(K::S5); + -(S1.powu(5) + + 10. * S1.powu(3) * S2 + + 20. * S1.powu(2) * S3 + + 15. * S1 * (S2.powu(2) + 2. * S4) + + 4. * (5. * S2 * S3 + 6. * S5)) + / n +} + +/// Mellin transform of $(1-x)^2\ln(1-x)$. +pub fn lm11m2(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + (5. + 3. * n - (2. * (1. + n) * (2. + n) * S1) / n) / ((1. + n).powu(2) * (2. + n).powu(2)) +} + +/// Mellin transform of $(1-x)^2\ln^2(1-x)$. +pub fn lm12m2(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + (2. * (n * (-9. - 8. * n + n.powu(3)) + - n * (10. + 21. * n + 14. * n.powu(2) + 3. * n.powu(3)) * S1 + + pow(2. + 3. * n + n.powu(2), 2) * S1.powu(2) + + pow(2. + 3. * n + n.powu(2), 2) * S2)) + / (n * (1. + n).powu(3) * (2. + n).powu(3)) +} + +/// Mellin transform of $(1-x)^2\ln^3(1-x)$. +pub fn lm13m2(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + (-6. * n * (-17. - 21. * n - 2. * n.powu(2) + 6. * n.powu(3) + 2. * n.powu(4)) + + 3. * n * (5. + 3. * n) * pow(2. + 3. * n + n.powu(2), 2) * S1.powu(2) + - 2. * pow(2. + 3. * n + n.powu(2), 3) * S1.powu(3) + + 3. * n * (5. + 3. * n) * pow(2. + 3. * n + n.powu(2), 2) * S2 + - 6. * (2. + 3. * n + n.powu(2)) + * S1 + * (n * (-9. - 8. * n + n.powu(3)) + pow(2. + 3. * n + n.powu(2), 2) * S2) + - 4. * pow(2. + 3. * n + n.powu(2), 3) * S3) + / (n * (1. + n).powu(4) * (2. + n).powu(4)) +} + +/// Mellin transform of $(1-x)^2\ln^4(1-x)$. +pub fn lm14m2(c: &mut Cache) -> Complex { + let n = c.n(); + let S1 = c.get(K::S1); + let S2 = c.get(K::S2); + let S3 = c.get(K::S3); + let S4 = c.get(K::S4); + 2. / (n * (1. + n).powu(5) * pow(2. + n, 5)) + * (12. * n * (-33. + n * (-54. + n * (-15. + n * (20. + 3. * n * (5. + n))))) + - 2. * n * (1. + n).powu(3) * pow(2. + n, 3) * (5. + 3. * n) * S1.powu(3) + + (1. + n).powu(4) * pow(2. + n, 4) * S1.powu(4) + + 6. * n * (1. + n).powu(2) * pow(2. + n, 2) * (-9. - 8. * n + n.powu(3)) * S2 + + 3. * (1. + n).powu(4) * pow(2. + n, 4) * S2.powu(2) + + 6. * (1. + n).powu(2) + * pow(2. + n, 2) + * S1.powu(2) + * (n * (-9. - 8. * n + n.powu(3)) + (1. + n).powu(2) * pow(2. + n, 2) * S2) + - 4. * n * (1. + n).powu(3) * pow(2. + n, 3) * (5. + 3. * n) * S3 + + 2. * (1. + n) + * (2. + n) + * S1 + * (6. * n * (-17. + n * (-21. + 2. * n * (-1. + n * (3. + n)))) + - 3. * n * (1. + n).powu(2) * pow(2. + n, 2) * (5. + 3. * n) * S2 + + 4. * (1. + n).powu(3) * pow(2. + n, 3) * S3) + + 6. * (1. + n).powu(4) * pow(2. + n, 4) * S4) +} + +#[cfg(test)] +mod tests { + + use super::*; + use crate::harmonics::cache::Cache; + use crate::{assert_approx_eq_cmplx, cmplx}; + use num::complex::Complex; + use num::traits::Pow; + use rgsl::integration::qng; + + #[test] + fn test_lm1Xm1() { + fn mellin_lm1Xm1(x: f64, k: f64, n: f64) -> f64 { + f64::pow(x, n - 1.0) * (1.0 - x) * f64::pow((1.0 - x).ln(), k) + } + + const NS: [f64; 5] = [1.0, 1.5, 2.0, 2.34, 56.789]; + for N in NS { + let mut c = Cache::new(cmplx!(N, 0.)); + + let ref_values = [ + lm11m1(&mut c), + lm12m1(&mut c), + lm13m1(&mut c), + lm14m1(&mut c), + ]; + + for (k, &ref_value) in ref_values.iter().enumerate() { + let test_value = qng( + |x| mellin_lm1Xm1(x, k as f64 + 1.0, N), + 0.0, + 1.0, + 1e-3, + 1e-4, + ) + .unwrap(); + assert_approx_eq_cmplx!(f64, cmplx!(test_value.0, 0.0), ref_value, epsilon = 1e-6); + } + + // here the integration is poorly convergent + let ref_values5 = lm15m1(&mut c); + let test_value = + qng(|x| mellin_lm1Xm1(x, 4.0 + 1.0, N), 0.0, 0.9999, 1e-2, 1e-4).unwrap(); + assert_approx_eq_cmplx!(f64, cmplx!(test_value.0, 0.0), ref_values5, epsilon = 1e-3); + } + } + + #[test] + fn test_lmXm2() { + fn mellin_lm1Xm2(x: f64, k: f64, n: f64) -> f64 { + f64::pow(x, n - 1.0) * f64::pow(1.0 - x, 2) * f64::pow((1.0 - x).ln(), k) + } + + const NS: [f64; 5] = [1.0, 1.5, 2.0, 2.34, 56.789]; + for N in NS { + let mut c = Cache::new(cmplx!(N, 0.)); + + let ref_values = [ + lm11m2(&mut c), + lm12m2(&mut c), + lm13m2(&mut c), + lm14m2(&mut c), + ]; + + for (k, &ref_value) in ref_values.iter().enumerate() { + let test_value = qng( + |x| mellin_lm1Xm2(x, k as f64 + 1.0, N), + 0.0, + 1.0, + 1e-4, + 1e-4, + ) + .unwrap(); + assert_approx_eq_cmplx!(f64, cmplx!(test_value.0, 0.0), ref_value, epsilon = 1e-6); + } + } + } +} diff --git a/crates/ekore/src/harmonics/w1.rs b/crates/ekore/src/harmonics/w1.rs index 64431e6d0..2421d12d8 100644 --- a/crates/ekore/src/harmonics/w1.rs +++ b/crates/ekore/src/harmonics/w1.rs @@ -1,6 +1,7 @@ //! Harmonic sums of weight 1. use num::complex::Complex; +use crate::harmonics::cache::{Cache, K}; use crate::harmonics::polygamma::cern_polygamma; /// Compute the harmonic sum $S_1(N)$. @@ -12,11 +13,15 @@ pub fn S1(N: Complex) -> Complex { } /// Analytic continuation of harmonic sum $S_{-1}(N)$ for even moments. -pub fn Sm1e(hS1: Complex, hS1h: Complex) -> Complex { +pub fn Sm1e(c: &mut Cache) -> Complex { + let hS1 = c.get(K::S1); + let hS1h = c.get(K::S1h); hS1h - hS1 } /// Analytic continuation of harmonic sum $S_{-1}(N)$ for odd moments. -pub fn Sm1o(hS1: Complex, hS1mh: Complex) -> Complex { +pub fn Sm1o(c: &mut Cache) -> Complex { + let hS1 = c.get(K::S1); + let hS1mh = c.get(K::S1mh); hS1mh - hS1 } diff --git a/crates/ekore/src/harmonics/w2.rs b/crates/ekore/src/harmonics/w2.rs index a6d7cc89f..643ed19ed 100644 --- a/crates/ekore/src/harmonics/w2.rs +++ b/crates/ekore/src/harmonics/w2.rs @@ -2,6 +2,7 @@ use num::complex::Complex; use crate::constants::ZETA2; +use crate::harmonics::cache::{Cache, K}; use crate::harmonics::polygamma::cern_polygamma; /// Compute the harmonic sum $S_2(N)$. @@ -13,11 +14,15 @@ pub fn S2(N: Complex) -> Complex { } /// Analytic continuation of harmonic sum $S_{-2}(N)$ for even moments. -pub fn Sm2e(hS2: Complex, hS2h: Complex) -> Complex { +pub fn Sm2e(c: &mut Cache) -> Complex { + let hS2 = c.get(K::S2); + let hS2h = c.get(K::S2h); 1. / 2. * hS2h - hS2 } /// Analytic continuation of harmonic sum $S_{-2}(N)$ for odd moments. -pub fn Sm2o(hS2: Complex, hS2mh: Complex) -> Complex { +pub fn Sm2o(c: &mut Cache) -> Complex { + let hS2 = c.get(K::S2); + let hS2mh = c.get(K::S2mh); 1. / 2. * hS2mh - hS2 } diff --git a/crates/ekore/src/harmonics/w3.rs b/crates/ekore/src/harmonics/w3.rs index bc7286dcb..81bc0317c 100644 --- a/crates/ekore/src/harmonics/w3.rs +++ b/crates/ekore/src/harmonics/w3.rs @@ -4,6 +4,7 @@ use num::traits::Pow; use std::f64::consts::LN_2; use crate::constants::{ZETA2, ZETA3}; +use crate::harmonics::cache::{Cache, K}; use crate::harmonics::g_functions::g3; use crate::harmonics::polygamma::cern_polygamma; @@ -16,26 +17,29 @@ pub fn S3(N: Complex) -> Complex { } /// Analytic continuation of harmonic sum $S_{-3}(N)$ for even moments. -pub fn Sm3e(hS3: Complex, hS3h: Complex) -> Complex { +pub fn Sm3e(c: &mut Cache) -> Complex { + let hS3 = c.get(K::S3); + let hS3h = c.get(K::S3h); 1. / (2.).pow(2) * hS3h - hS3 } /// Analytic continuation of harmonic sum $S_{-3}(N)$ for odd moments. -pub fn Sm3o(hS3: Complex, hS3mh: Complex) -> Complex { +pub fn Sm3o(c: &mut Cache) -> Complex { + let hS3 = c.get(K::S3); + let hS3mh = c.get(K::S3mh); 1. / (2.).pow(2) * hS3mh - hS3 } /// Analytic continuation of harmonic sum $S_{-2,1}(N)$ for even moments. -pub fn Sm21e(N: Complex, hS1: Complex, hSm1: Complex) -> Complex { - Sm21(N, 1., hS1, hSm1) +pub fn Sm21e(c: &mut Cache) -> Complex { + let hSm1 = c.get(K::Sm1e); + let mut cp1 = Cache::new(c.n() + 1.); + ZETA2 * hSm1 - 5. / 8. * ZETA3 + ZETA2 * LN_2 - g3(&mut cp1) } /// Analytic continuation of harmonic sum $S_{-2,1}(N)$ for odd moments. -pub fn Sm21o(N: Complex, hS1: Complex, hSm1: Complex) -> Complex { - Sm21(N, -1., hS1, hSm1) -} - -/// Analytic continuation of harmonic sum $S_{-2,1}(N)$ for odd moments. -fn Sm21(N: Complex, eta: f64, hS1: Complex, hSm1: Complex) -> Complex { - -eta * g3(N + 1., hS1 + 1. / (N + 1.)) + ZETA2 * hSm1 - 5. / 8. * ZETA3 + ZETA2 * LN_2 +pub fn Sm21o(c: &mut Cache) -> Complex { + let hSm1 = c.get(K::Sm1o); + let mut cp1 = Cache::new(c.n() + 1.); + ZETA2 * hSm1 - 5. / 8. * ZETA3 + ZETA2 * LN_2 + g3(&mut cp1) } diff --git a/crates/ekore/src/harmonics/w5.rs b/crates/ekore/src/harmonics/w5.rs new file mode 100644 index 000000000..e77730a5b --- /dev/null +++ b/crates/ekore/src/harmonics/w5.rs @@ -0,0 +1,13 @@ +//! Harmonic sums of weight 5. +use num::complex::Complex; + +use crate::constants::ZETA5; +use crate::harmonics::polygamma::cern_polygamma; + +/// Compute the harmonic sum $S_5(N)$. +/// +/// $$S_5(N) = \sum\limits_{j=1}^N \frac 1 {j^5} = - \frac 1 24 \psi_4(N+1)+\zeta(5)$$ +/// with $\psi_4(N)$ the 4rd polygamma function and $\zeta$ the Riemann zeta function. +pub fn S5(N: Complex) -> Complex { + ZETA5 + 1.0 / 24.0 * cern_polygamma(N + 1.0, 4) +} diff --git a/poetry.lock b/poetry.lock index fbb303a53..a46d17559 100644 --- a/poetry.lock +++ b/poetry.lock @@ -159,6 +159,24 @@ charset-normalizer = ["charset-normalizer"] html5lib = ["html5lib"] lxml = ["lxml"] +[[package]] +name = "bibtexparser" +version = "2.0.0b8" +description = "Bibtex parser for python 3" +optional = false +python-versions = "*" +files = [ + {file = "bibtexparser-2.0.0b8-py3-none-any.whl", hash = "sha256:eda1e6dbc7bbc040d351c926ead97a36f4734858656f3514a933259878bd37f1"}, + {file = "bibtexparser-2.0.0b8.tar.gz", hash = "sha256:11eaaa1fe80154f70092e91fa7e3fed5e5c81021f6c6c09151e36b8533f36994"}, +] + +[package.dependencies] +pylatexenc = ">=2.10" + +[package.extras] +docs = ["sphinx"] +test = ["jupyter", "pytest", "pytest-cov", "pytest-xdist"] + [[package]] name = "bleach" version = "6.2.0" @@ -2304,6 +2322,16 @@ files = [ [package.extras] windows-terminal = ["colorama (>=0.4.6)"] +[[package]] +name = "pylatexenc" +version = "2.10" +description = "Simple LaTeX parser providing latex-to-unicode and unicode-to-latex conversion" +optional = false +python-versions = "*" +files = [ + {file = "pylatexenc-2.10.tar.gz", hash = "sha256:3dd8fd84eb46dc30bee1e23eaab8d8fb5a7f507347b23e5f38ad9675c84f40d3"}, +] + [[package]] name = "pylint" version = "3.3.4" diff --git a/pyproject.toml b/pyproject.toml index a08eba08e..21976847b 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -61,6 +61,7 @@ sphinx-rtd-theme = "^1.0.0" sphinxcontrib-bibtex = "^2.4.1" nbsphinx = "^0.8.8" ipykernel = "^6.13.0" +bibtexparser = ">=2.0.0b8" [tool.poetry.group.test] optional = true diff --git a/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_as4_fhmv.py b/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_as4_fhmv.py index 4f778c0dd..30891649e 100644 --- a/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_as4_fhmv.py +++ b/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_as4_fhmv.py @@ -269,4 +269,4 @@ def test_gamma_ps_extrapolation(): my_res = [] for nf in [3, 4, 5]: my_res.append(gps.gamma_ps(N, nf, sx_cache, 0)) - np.testing.assert_allclose(n22_ref, n22_ref) + np.testing.assert_allclose(my_res, n22_ref, rtol=3e-5)